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: ! petscsys.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 finclude/petscsys.h
22: #include finclude/petscvec.h
23: #include finclude/petscmat.h
24: #include "finclude/petscmat.h90"
26: Mat A
27: PetscErrorCode ierr
28: PetscInt i,j,m,n,iar(1),jar(1)
29: PetscInt rstart,rend
30: PetscInt one
31: PetscScalar v(1)
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: iar(1) = i
51: do 20, j=0,n-1
52: jar(1) = j
53: v(1) = 9.d0/(i+j+1)
54: call MatSetValues(A,one,iar,one,jar,v,INSERT_VALUES,ierr)
55: 20 continue
56: 10 continue
58: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
59: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
61: !
62: ! Print the matrix to the screen
63: !
64: call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
67: !
68: ! Print the local portion of the matrix to the screen
69: !
70: call MatGetArrayF90(A,array,ierr)
71: call MatGetOwnershipRange(A,rstart,rend,ierr)
72: call PetscSequentialPhaseBegin(PETSC_COMM_WORLD,1,ierr)
73: do 30 i=1,rend-rstart
74: write(6,100) (PetscRealPart(array(i,j)),j=1,n)
75: 30 continue
76: 100 format(2F6.2)
78: call PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1,ierr)
80: call MatRestoreArrayF90(A,array,ierr)
81: !
82: ! Free the space used by the matrix
83: !
84: call MatDestroy(A,ierr)
85: call PetscFinalize(ierr)
86: end