Actual source code: ex16f90.F
1: !
2: !
3: ! Tests MatGetArray() on a dense matrix
4: !
6: program main
7: implicit none
10: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
11: ! Include files
12: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
13: !
14: ! The following include statements are required for Fortran programs
15: ! that use PETSc vectors:
16: ! petsc.h - base PETSc routines
17: ! petscvec.h - defines (INSERT_VALUES)
18: ! petscmat.h - matrices
19: ! petscmat.h90 - to allow access to Fortran 90 features of matrices
21: #include include/finclude/petsc.h
22: #include include/finclude/petscvec.h
23: #include include/finclude/petscmat.h
24: #include "include/finclude/petscmat.h90"
26: Mat A
27: PetscErrorCode ierr
28: PetscInt i,j,m,n
29: PetscInt rstart,rend
30: PetscInt one
31: PetscScalar v
32: PetscScalar, pointer :: array(:,:)
35: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
37: m = 3
38: n = 2
39: one = 1
40: !
41: ! Create a parallel dense matrix shared by all processors
42: !
43: call MatCreateMPIDense(PETSC_COMM_WORLD,PETSC_DECIDE, &
44: & PETSC_DECIDE,m,n,PETSC_NULL_SCALAR,A,ierr)
46: !
47: ! Set values into the matrix. All processors set all values.
48: !
49: do 10, i=0,m-1
50: do 20, j=0,n-1
51: v = 9.d0/(i+j+1)
52: call MatSetValues(A,one,i,one,j,v,INSERT_VALUES,ierr)
53: 20 continue
54: 10 continue
56: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
57: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
59: !
60: ! Print the matrix to the screen
61: !
62: call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
65: !
66: ! Print the local portion of the matrix to the screen
67: !
68: call MatGetArrayF90(A,array,ierr)
69: call MatGetOwnershipRange(A,rstart,rend,ierr)
70: call PetscSequentialPhaseBegin(PETSC_COMM_WORLD,1,ierr)
71: do 30 i=1,rend-rstart
72: write(6,100) (PetscRealPart(array(i,j)),j=1,n)
73: 30 continue
74: 100 format(2F6.2)
76: call PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1,ierr)
78: call MatRestoreArrayF90(A,array,ierr)
79: !
80: ! Free the space used by the matrix
81: !
82: call MatDestroy(A,ierr)
83: call PetscFinalize(ierr)
84: end