Actual source code: ex16f90.F

  1: !
  2: !  "$Id: ex16f90.F,v 1.17 2001/08/07 03:03:07 balay Exp $";
  3: !
  4: !
  5: !  Tests MatGetArray() on a dense matrix
  6: !

  8:       program main
  9:       implicit none


 12: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 13: !                    Include files
 14: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 15: !
 16: !  The following include statements are required for Fortran programs
 17: !  that use PETSc vectors:
 18: !     petsc.h  - base PETSc routines
 19: !     petscvec.h - defines (INSERT_VALUES)
 20: !     petscmat.h    - matrices
 21: !     petscmat.h90  - to allow access to Fortran 90 features of matrices

 23:  #include include/finclude/petsc.h
 24:  #include include/finclude/petscvec.h
 25:  #include include/finclude/petscmat.h
 26: #include "include/finclude/petscmat.h90"

 28:       Mat A
 29:       integer i,j,m,n,ierr
 30:       integer rstart,rend
 31:       PetscScalar  v
 32:       PetscScalar, pointer :: array(:,:)


 35:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)

 37:       m = 3
 38:       n = 2
 39: !
 40: !      Create a parallel dense matrix shared by all processors
 41: !
 42:       call MatCreateMPIDense(PETSC_COMM_WORLD,PETSC_DECIDE,                &
 43:      &                       PETSC_DECIDE,m,n,PETSC_NULL_SCALAR,A,ierr)

 45: !
 46: !     Set values into the matrix. All processors set all values.
 47: !
 48:       do 10, i=0,m-1
 49:         do 20, j=0,n-1
 50:           v = 9.d0/(i+j+1)
 51:           call MatSetValues(A,1,i,1,j,v,INSERT_VALUES,ierr)
 52:  20     continue
 53:  10   continue

 55:       call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
 56:       call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)

 58: !
 59: !       Print the matrix to the screen
 60: !
 61:       call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)


 64: !
 65: !      Print the local portion of the matrix to the screen
 66: !
 67:       call MatGetArrayF90(A,array,ierr)
 68:       call MatGetOwnershipRange(A,rstart,rend,ierr)
 69:       call PetscSequentialPhaseBegin(PETSC_COMM_WORLD,1,ierr)
 70:       do 30 i=1,rend-rstart
 71:          write(6,100) (PetscRealPart(array(i,j)),j=1,n)
 72:  30   continue
 73:  100  format(2F6.2)

 75:       call PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1,ierr)

 77:       call MatRestoreArrayF90(A,array,ierr)
 78: !
 79: !      Free the space used by the matrix
 80: !
 81:       call MatDestroy(A,ierr)
 82:       call PetscFinalize(ierr)
 83:       end