Actual source code: ex63f.F
petsc-3.12.0 2019-09-29
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: #include <petsc/finclude/petscmat.h>
10: use petscmat
11: implicit none
13: PetscErrorCode ierr
14: PetscInt row,col,ten
15: PetscMPIInt rank
16: PetscScalar v
17: Mat A
18: PetscViewer view
20: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
21: if (ierr .ne. 0) then
22: print*,'Unable to initialize PETSc'
23: stop
24: endif
26: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
27: !
28: ! Proc-0 Create a seq-dense matrix and write it to a file
29: !
30: if (rank .eq. 0) then
31: ten = 10
32: call MatCreateSeqDense(PETSC_COMM_SELF,ten,ten, &
33: & PETSC_NULL_SCALAR,A,ierr)
34: v = 1.0
35: do row=0,9
36: do col=0,9
37: call MatSetValue(A,row,col,v,INSERT_VALUES,ierr)
38: v = v + 1.0
39: enddo
40: enddo
42: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
43: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
45: call PetscObjectSetName(A,'Original Matrix',ierr)
46: call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
47: !
48: ! Now Write this matrix to a binary file
49: !
50: call PetscViewerBinaryOpen(PETSC_COMM_SELF,'dense.mat', &
51: & FILE_MODE_WRITE,view,ierr)
52: call MatView(A,view,ierr)
53: call PetscViewerDestroy(view,ierr)
54: call MatDestroy(A,ierr)
55: !
56: ! Read this matrix into a SeqDense matrix
58: call PetscViewerBinaryOpen(PETSC_COMM_SELF,'dense.mat', &
59: & FILE_MODE_READ,view,ierr)
60: call MatCreate(PETSC_COMM_SELF,A,ierr)
61: call MatSetType(A, MATSEQDENSE,ierr)
62: call MatLoad(A,view,ierr)
64: call PetscObjectSetName(A,'SeqDense Matrix read in from file', &
65: & ierr)
66: call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
67: call MatDestroy(A,ierr)
68: call PetscViewerDestroy(view,ierr)
69: endif
71: !
72: ! Use a barrier, so that the procs do not try opening the file before
73: ! it is created.
74: !
75: call MPI_Barrier(PETSC_COMM_WORLD,ierr)
77: call PetscViewerBinaryOpen(PETSC_COMM_WORLD,'dense.mat', &
78: & FILE_MODE_READ,view,ierr)
79: call MatCreate(PETSC_COMM_WORLD,A,ierr)
80: call MatSetType(A, MATMPIDENSE,ierr)
81: call MatLoad(A,view,ierr)
83: call PetscObjectSetName(A, 'MPIDense Matrix read in from file', &
84: & ierr)
85: call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
86: call MatDestroy(A,ierr)
87: call PetscViewerDestroy(view,ierr)
88: call PetscFinalize(ierr)
89: end
91: !/*TEST
92: !
93: ! test:
94: ! nsize: 2
95: ! output_file: output/ex63_1.out
96: !
97: !TEST*/