Actual source code: ex63f.F
petsc-3.5.4 2015-05-23
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)
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
39: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
40: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
42: call PetscObjectSetName(A,"Original Matrix",ierr)
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 MatView(A,view,ierr)
50: call PetscViewerDestroy(view,ierr)
51: call MatDestroy(A,ierr)
52: !
53: ! Read this matrix into a SeqDense matrix
55: call PetscViewerBinaryOpen(PETSC_COMM_SELF,'dense.mat', &
56: & FILE_MODE_READ,view,ierr)
57: call MatCreate(PETSC_COMM_SELF,A,ierr)
58: call MatSetType(A, MATSEQDENSE,ierr)
59: call MatLoad(A,view,ierr)
61: call PetscObjectSetName(A,"SeqDense Matrix read in from file", &
62: & ierr)
63: call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
64: call MatDestroy(A,ierr)
65: call PetscViewerDestroy(view,ierr)
66: endif
68: !
69: ! Use a barrier, so that the procs do not try opening the file before
70: ! it is created.
71: !
72: call MPI_Barrier(PETSC_COMM_WORLD,ierr)
74: call PetscViewerBinaryOpen(PETSC_COMM_WORLD,'dense.mat', &
75: & FILE_MODE_READ,view,ierr)
76: call MatCreate(PETSC_COMM_WORLD,A,ierr)
77: call MatSetType(A, MATMPIDENSE,ierr)
78: call MatLoad(A,view,ierr)
80: call PetscObjectSetName(A, "MPIDense Matrix read in from file", &
81: & ierr)
82: call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
83: call MatDestroy(A,ierr)
84: call PetscViewerDestroy(view,ierr)
85: call PetscFinalize(ierr)
86: end