Actual source code: ex63.F
1: !
2: ! "$Id: ex63.F,v 1.11 2001/08/07 03:03:07 balay Exp $";
3: !
4: ! This program tests storage of PETSc Dense matrix.
5: ! It Creates a Dense matrix, and stores it into a file,
6: ! and then reads it back in as a SeqDense and MPIDense
7: ! matrix, and prints out the contents.
8: !
9: program main
10: implicit none
11: #include include/finclude/petsc.h
12: #include include/finclude/petscmat.h
13: #include include/finclude/petscvec.h
14: #include include/finclude/petscviewer.h
15: integer ierr,row,col,rank
16: PetscScalar v
17: Mat A
18: PetscViewer view
20: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
21:
22: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
23: !
24: ! Proc-0 Create a seq-dense matrix and write it to a file
25: !
26: if (rank .eq. 0) then
27: call MatCreateSeqDense(PETSC_COMM_SELF,10,10, &
28: & PETSC_NULL_SCALAR,A,ierr)
29: v = 1.0d0
30: do row=0,9
31: do col=0,9
32: call MatSetValue(A,row,col,v,INSERT_VALUES,ierr)
33: v = v + 1.0d0
34: enddo
35: enddo
36:
37: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
38: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
39: write (6,100)
40: 100 format('The Contents of the Original Matrix')
41: call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
42: !
43: ! Now Write this matrix to a binary file
44: !
45: call PetscViewerBinaryOpen(PETSC_COMM_SELF,'dense.mat', &
46: & PETSC_BINARY_CREATE,view,ierr)
47: call PetscViewerSetFormat(view,PETSC_VIEWER_BINARY_NATIVE,ierr)
48: call MatView(A,view,ierr)
49: call PetscViewerDestroy(view,ierr)
50: call MatDestroy(A,ierr)
51: !
52: ! Read this matrix into a SeqDense matrix
54: call PetscViewerBinaryOpen(PETSC_COMM_SELF,'dense.mat', &
55: & PETSC_BINARY_RDONLY,view,ierr)
56: call MatLoad(view,MATSEQDENSE,A,ierr)
57: write (6,200)
58: 200 format('The Contents of SeqDense Matrix read in from file')
59: call MatView(A,PETSC_VIEWER_STDOUT_SELF,ierr)
60: call MatDestroy(A,ierr)
61: call PetscViewerDestroy(view,ierr)
62: endif
64: !
65: ! Use a barrier, so that the procs do not try opening the file before
66: ! it is created.
67: !
68: call MPI_Barrier(PETSC_COMM_WORLD,ierr)
70: call PetscViewerBinaryOpen(PETSC_COMM_WORLD,'dense.mat', &
71: & PETSC_BINARY_RDONLY,view,ierr)
72: call MatLoad(view,MATMPIDENSE,A,ierr)
73: if (rank .eq. 0) then
74: write (6,300)
75: 300 format('The Contents of MPIDense Matrix read in from file')
76: endif
77: call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
78: call MatDestroy(A,ierr)
79: call PetscViewerDestroy(view,ierr)
80: call PetscFinalize(ierr)
81: end