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: integer i,j,m,n,ierr
28: integer rstart,rend
29: PetscScalar v
30: PetscScalar, pointer :: array(:,:)
33: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
35: m = 3
36: n = 2
37: !
38: ! Create a parallel dense matrix shared by all processors
39: !
40: call MatCreateMPIDense(PETSC_COMM_WORLD,PETSC_DECIDE, &
41: & PETSC_DECIDE,m,n,PETSC_NULL_SCALAR,A,ierr)
43: !
44: ! Set values into the matrix. All processors set all values.
45: !
46: do 10, i=0,m-1
47: do 20, j=0,n-1
48: v = 9.d0/(i+j+1)
49: call MatSetValues(A,1,i,1,j,v,INSERT_VALUES,ierr)
50: 20 continue
51: 10 continue
53: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
54: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
56: !
57: ! Print the matrix to the screen
58: !
59: call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
62: !
63: ! Print the local portion of the matrix to the screen
64: !
65: call MatGetArrayF90(A,array,ierr)
66: call MatGetOwnershipRange(A,rstart,rend,ierr)
67: call PetscSequentialPhaseBegin(PETSC_COMM_WORLD,1,ierr)
68: do 30 i=1,rend-rstart
69: write(6,100) (PetscRealPart(array(i,j)),j=1,n)
70: 30 continue
71: 100 format(2F6.2)
73: call PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1,ierr)
75: call MatRestoreArrayF90(A,array,ierr)
76: !
77: ! Free the space used by the matrix
78: !
79: call MatDestroy(A,ierr)
80: call PetscFinalize(ierr)
81: end