Actual source code: ex16f90.F90
petsc-3.9.2 2018-05-20
1: !
2: !
3: ! Tests MatDenseGetArray()
4: !
6: program main
7: #include <petsc/finclude/petscmat.h>
8: use petscmat
9: implicit none
11: Mat A
12: PetscErrorCode ierr
13: PetscInt i,j,m,n,iar(1),jar(1)
14: PetscInt rstart,rend
15: PetscInt one
16: PetscScalar v(1)
17: PetscScalar, pointer :: array(:,:)
20: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
21: if (ierr .ne. 0) then
22: print*,'Unable to initialize PETSc'
23: stop
24: endif
26: m = 3
27: n = 2
28: one = 1
29: !
30: ! Create a parallel dense matrix shared by all processors
31: !
32: call MatCreateDense(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,m,n,PETSC_NULL_SCALAR,A,ierr);CHKERRA(ierr)
34: !
35: ! Set values into the matrix. All processors set all values.
36: !
37: do 10, i=0,m-1
38: iar(1) = i
39: do 20, j=0,n-1
40: jar(1) = j
41: v(1) = 9.0/real(i+j+1)
42: call MatSetValues(A,one,iar,one,jar,v,INSERT_VALUES,ierr);CHKERRA(ierr)
43: 20 continue
44: 10 continue
46: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr);CHKERRA(ierr)
47: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr);CHKERRA(ierr)
49: !
50: ! Print the matrix to the screen
51: !
52: call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRA(ierr)
55: !
56: ! Print the local portion of the matrix to the screen
57: !
58: call MatDenseGetArrayF90(A,array,ierr);CHKERRA(ierr)
59: call MatGetOwnershipRange(A,rstart,rend,ierr);CHKERRA(ierr)
60: call PetscSequentialPhaseBegin(PETSC_COMM_WORLD,1,ierr);CHKERRA(ierr)
61: !
62: ! Fortran IO may not come out in the correct order since each process
63: ! is individually doing IO
64: ! do 30 i=1,rend-rstart
65: ! write(6,100) (PetscRealPart(array(i,j)),j=1,n)
66: ! 30 continue
67: ! 100 format(2F6.2)
69: call PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1,ierr);CHKERRA(ierr)
71: call MatDenseRestoreArrayF90(A,array,ierr);CHKERRA(ierr)
72: !
73: ! Free the space used by the matrix
74: !
75: call MatDestroy(A,ierr);CHKERRA(ierr)
76: call PetscFinalize(ierr)
77: end
80: !/*TEST
81: !
82: ! test:
83: ! nsize: 2
84: ! filter: sort -b
85: ! filter_output: sort -b
86: !
87: !TEST*/