Actual source code: ex36f.F

  1: !
  2: !      "$Id: ex36.F,v 1.30 2001/08/07 03:03:07 balay Exp $";
  3: !
  4: !   This program demonstrates use of PETSc dense matrices.
  5: !
  6:       program main

 8:  #include include/finclude/petsc.h

 10:       integer err

 12:       call PetscInitialize(PETSC_NULL_CHARACTER,err)

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

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

 20:       call PetscFinalize(err)
 21:       end

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

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

 41:       Mat         A
 42:       integer     n,m,err
 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,err)
 52: !      call MatCreateMPIDense(PETSC_COMM_WORLD,m,n,m,n,PETSC_NULL_SCALAR,A,err)

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

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

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

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

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

 69: !  Clean up
 70:       call MatDestroy(A,err)
 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:       integer   n,m,err
 86:       parameter (m=5,n=4)
 87:       Mat       A
 88:       PetscScalar    aa(m,n)

 90: !  Create matrix
 91:       call MatCreateSeqDense(PETSC_COMM_SELF,m,n,aa,A,err)

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

 96: !  Finalize matrix assembly
 97:       call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,err)
 98:       call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,err)

100: !  View matrix
101:       call MatView(A,PETSC_VIEWER_STDOUT_SELF,err)

103: !  Clean up
104:       call MatDestroy(A,err)
105:       return
106:       end

108: ! -----------------------------------------------------------------

110:       subroutine FillUpMatrix(m,n,X)
111:       integer          m,n,i,j
112:       PetscScalar      X(m,n)

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