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