Actual source code: ex36f.F
petsc-3.7.5 2017-01-01
1: !
2: !
3: ! This program demonstrates use of PETSc dense matrices.
4: !
5: program main
7: #include <petsc/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 MatDenseGetArray() 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 <petsc/finclude/petscsys.h>
37: #include <petsc/finclude/petscmat.h>
38: #include <petsc/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)
60: call MatSetUp(A,ierr)
62: ! Access array storage
63: call MatDenseGetArray(A,aa,ia,ierr)
65: ! Set matrix values directly
66: call FillUpMatrix(m,n,aa(ia+1))
68: call MatDenseRestoreArray(A,aa,ia,ierr)
70: ! Finalize matrix assembly
71: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
72: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
74: ! View matrix
75: call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
77: ! Clean up
78: call MatDestroy(A,ierr)
79: return
80: end
82: ! -----------------------------------------------------------------
83: !
84: ! Demo2 - This subroutine demonstrates the use of user-allocated dense
85: ! matrix storage.
86: !
87: subroutine Demo2()
89: #include <petsc/finclude/petscsys.h>
90: #include <petsc/finclude/petscmat.h>
91: #include <petsc/finclude/petscviewer.h>
93: PetscInt n,m
94: PetscErrorCode ierr
95: parameter (m=5,n=4)
96: Mat A
97: PetscScalar aa(m,n)
99: ! Create matrix
100: call MatCreateSeqDense(PETSC_COMM_SELF,m,n,aa,A,ierr)
102: ! Set matrix values directly
103: call FillUpMatrix(m,n,aa)
105: ! Finalize matrix assembly
106: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
107: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
109: ! View matrix
110: call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
112: ! Clean up
113: call MatDestroy(A,ierr)
114: return
115: end
117: ! -----------------------------------------------------------------
119: subroutine FillUpMatrix(m,n,X)
120: PetscInt m,n,i,j
121: PetscScalar X(m,n)
123: do 10, j=1,n
124: do 20, i=1,m
125: X(i,j) = 1.0/real(i+j-1)
126: 20 continue
127: 10 continue
128: return
129: end