Actual source code: ex63f.F

  1: !
  2: !
  3: !   This program tests storage of PETSc Dense matrix.
  4: !   It Creates a Dense matrix, and stores it into a file,
  5: !   and then reads it back in as a SeqDense and MPIDense
  6: !   matrix, and prints out the contents.
  7: !
  8:       program main
  9:       implicit none
 10:  #include finclude/petscsys.h
 11:  #include finclude/petscmat.h
 12:  #include finclude/petscvec.h
 13:  #include finclude/petscviewer.h
 14:       PetscErrorCode ierr
 15:       PetscInt row,col,ten
 16:       PetscMPIInt rank
 17:       PetscScalar  v
 18:       Mat     A
 19:       PetscViewer  view

 21:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 22: 
 23:       call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
 24: !
 25: !     Proc-0 Create a seq-dense matrix and write it to a file
 26: !
 27:       if (rank .eq. 0) then
 28:          ten = 10
 29:          call MatCreateSeqDense(PETSC_COMM_SELF,ten,ten,                  &
 30:      &        PETSC_NULL_SCALAR,A,ierr)
 31:          v = 1.0d0
 32:          do row=0,9
 33:             do col=0,9
 34:                call MatSetValue(A,row,col,v,INSERT_VALUES,ierr)
 35:                v = v + 1.0d0
 36:             enddo
 37:          enddo
 38: 
 39:          call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
 40:          call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
 41:          write (6,100)
 42:  100     format('The Contents of the Original Matrix')
 43:          call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
 44: !
 45: !        Now Write this matrix to a binary file
 46: !
 47:          call PetscViewerBinaryOpen(PETSC_COMM_SELF,'dense.mat',           &
 48:      &        FILE_MODE_WRITE,view,ierr)
 49:          call PetscViewerSetFormat(view,PETSC_VIEWER_NATIVE,ierr)
 50:          call MatView(A,view,ierr)
 51:          call PetscViewerDestroy(view,ierr)
 52:          call MatDestroy(A,ierr)
 53: !
 54: !        Read this matrix into a SeqDense matrix

 56:          call PetscViewerBinaryOpen(PETSC_COMM_SELF,'dense.mat',            &
 57:      &        FILE_MODE_READ,view,ierr)
 58:          call MatLoad(view,MATSEQDENSE,A,ierr)
 59:          write (6,200)
 60:  200     format('The Contents of SeqDense Matrix read in from file')
 61:          call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
 62:          call MatDestroy(A,ierr)
 63:          call PetscViewerDestroy(view,ierr)
 64:       endif

 66: !
 67: !     Use a barrier, so that the procs do not try opening the file before
 68: !     it is created.
 69: !
 70:       call MPI_Barrier(PETSC_COMM_WORLD,ierr)

 72:       call PetscViewerBinaryOpen(PETSC_COMM_WORLD,'dense.mat',               &
 73:      &     FILE_MODE_READ,view,ierr)
 74:       call MatLoad(view,MATMPIDENSE,A,ierr)
 75:       if (rank .eq. 0) then
 76:          write (6,300)
 77:  300     format('The Contents of MPIDense Matrix read in from file')
 78:       endif
 79:       call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
 80:       call MatDestroy(A,ierr)
 81:       call PetscViewerDestroy(view,ierr)
 82:       call PetscFinalize(ierr)
 83:       end