Actual source code: ex36f.F
1: !
2: !
3: ! This program demonstrates use of PETSc dense matrices.
4: !
5: program main
7: #include finclude/petscsys.h
9: PetscErrorCode ierr
11: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
13: ! Demo of PETSc-allocated dense matrix storage
14: call Demo1()
16: ! Demo of user-allocated dense matrix storage
17: call Demo2()
19: call PetscFinalize(ierr)
20: end
22: ! -----------------------------------------------------------------
23: !
24: ! Demo1 - This subroutine demonstrates the use of PETSc-allocated dense
25: ! matrix storage. Here MatGetArray() is used for direct access to the
26: ! array that stores the dense matrix. The user declares an array (aa(1))
27: ! and index variable (ia), which are then used together to manipulate
28: ! the array contents.
29: !
30: ! Note the use of PETSC_NULL_SCALAR in MatCreateSeqDense() to indicate that no
31: ! storage is being provided by the user. (Do NOT pass a zero in that
32: ! location.)
33: !
34: subroutine Demo1()
36: #include finclude/petscsys.h
37: #include finclude/petscmat.h
38: #include finclude/petscviewer.h
40: Mat A
41: PetscInt n,m
42: PetscErrorCode ierr
43: PetscScalar aa(1)
44: PetscOffset ia
46: n = 4
47: m = 5
49: ! Create matrix
51: ! call MatCreateSeqDense(PETSC_COMM_SELF,m,n,PETSC_NULL_SCALAR, &
52: ! & A,ierr)
54: ! Using MatCreate() instead of MatCreateSeqDense() as above to avoid Nag F90 errors
55: ! However both cases are equivalent
57: call MatCreate(PETSC_COMM_SELF,A,ierr)
58: call MatSetSizes(A,m,n,m,n,ierr)
59: call MatSetType(A,MATSEQDENSE,ierr)
61: ! Access array storage
62: call MatGetArray(A,aa,ia,ierr)
64: ! Set matrix values directly
65: call FillUpMatrix(m,n,aa(ia+1))
67: call MatRestoreArray(A,aa,ia,ierr)
69: ! Finalize matrix assembly
70: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
71: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
73: ! View matrix
74: call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
76: ! Clean up
77: call MatDestroy(A,ierr)
78: return
79: end
81: ! -----------------------------------------------------------------
82: !
83: ! Demo2 - This subroutine demonstrates the use of user-allocated dense
84: ! matrix storage.
85: !
86: subroutine Demo2()
88: #include finclude/petscsys.h
89: #include finclude/petscmat.h
90: #include finclude/petscviewer.h
92: PetscInt n,m
93: PetscErrorCode ierr
94: parameter (m=5,n=4)
95: Mat A
96: PetscScalar aa(m,n)
98: ! Create matrix
99: call MatCreateSeqDense(PETSC_COMM_SELF,m,n,aa,A,ierr)
101: ! Set matrix values directly
102: call FillUpMatrix(m,n,aa)
104: ! Finalize matrix assembly
105: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
106: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
108: ! View matrix
109: call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
111: ! Clean up
112: call MatDestroy(A,ierr)
113: return
114: end
116: ! -----------------------------------------------------------------
118: subroutine FillUpMatrix(m,n,X)
119: PetscInt m,n,i,j
120: PetscScalar X(m,n)
122: do 10, j=1,n
123: do 20, i=1,m
124: X(i,j) = 1.d0/(i+j-1)
125: 20 continue
126: 10 continue
127: return
128: end