Actual source code: ex36f.F

  1: !
  2: !
  3: !   This program demonstrates use of PETSc dense matrices.
  4: !
  5:       program main

 7:  #include include/finclude/petsc.h

  9:       PetscErrorCode ierr

 11:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)

 13: !  Demo of PETSc-allocated dense matrix storage
 14:       call Demo1()

 16: !  Demo of user-allocated dense matrix storage
 17:       call Demo2()

 19:       call PetscFinalize(ierr)
 20:       end

 22: ! -----------------------------------------------------------------
 23: !
 24: !  Demo1 -  This subroutine demonstrates the use of PETSc-allocated dense
 25: !  matrix storage.  Here MatGetArray() is used for direct access to the
 26: !  array that stores the dense matrix.  The user declares an array (aa(1))
 27: !  and index variable (ia), which are then used together to manipulate
 28: !  the array contents.
 29: !
 30: !  Note the use of PETSC_NULL_SCALAR in MatCreateSeqDense() to indicate that no
 31: !  storage is being provided by the user. (Do NOT pass a zero in that
 32: !  location.)
 33: !
 34:       subroutine Demo1()

 36:  #include include/finclude/petsc.h
 37:  #include include/finclude/petscmat.h
 38:  #include include/finclude/petscviewer.h

 40:       Mat         A
 41:       PetscInt   n,m
 42:       PetscErrorCode ierr
 43:       PetscScalar aa(1)
 44:       PetscOffset ia

 46:       n = 4
 47:       m = 5

 49: !  Create matrix
 50:       call MatCreateSeqDense(PETSC_COMM_SELF,m,n,PETSC_NULL_SCALAR,     &
 51:      &     A,ierr)
 52: !      call MatCreateMPIDense(PETSC_COMM_WORLD,m,n,m,n,PETSC_NULL_SCALAR,A,ierr)

 54: !  Access array storage
 55:       call MatGetArray(A,aa,ia,ierr)

 57: !  Set matrix values directly
 58:       call FillUpMatrix(m,n,aa(ia+1))

 60:       call MatRestoreArray(A,aa,ia,ierr)

 62: !  Finalize matrix assembly
 63:       call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
 64:       call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)

 66: !  View matrix
 67:       call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)

 69: !  Clean up
 70:       call MatDestroy(A,ierr)
 71:       return
 72:       end

 74: ! -----------------------------------------------------------------
 75: !
 76: !  Demo2 -  This subroutine demonstrates the use of user-allocated dense
 77: !  matrix storage.
 78: !
 79:       subroutine Demo2()

 81:  #include include/finclude/petsc.h
 82:  #include include/finclude/petscmat.h
 83:  #include include/finclude/petscviewer.h

 85:       PetscInt   n,m
 86:       PetscErrorCode ierr
 87:       parameter (m=5,n=4)
 88:       Mat       A
 89:       PetscScalar    aa(m,n)

 91: !  Create matrix
 92:       call MatCreateSeqDense(PETSC_COMM_SELF,m,n,aa,A,ierr)

 94: !  Set matrix values directly
 95:       call FillUpMatrix(m,n,aa)

 97: !  Finalize matrix assembly
 98:       call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
 99:       call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)

101: !  View matrix
102:       call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)

104: !  Clean up
105:       call MatDestroy(A,ierr)
106:       return
107:       end

109: ! -----------------------------------------------------------------

111:       subroutine FillUpMatrix(m,n,X)
112:       PetscInt          m,n,i,j
113:       PetscScalar      X(m,n)

115:       do 10, j=1,n
116:         do 20, i=1,m
117:           X(i,j) = 1.0/(i+j-1)
118:  20     continue
119:  10   continue
120:       return
121:       end