1 static char help[] = "Tests MatView()/MatLoad() with binary viewers for BAIJ matrices.\n\n";
2
3 #include <petscmat.h>
4 #include <petscviewer.h>
5
6 #include <petsc/private/hashtable.h>
MakeValue(PetscInt i,PetscInt j,PetscInt M)7 static PetscReal MakeValue(PetscInt i,PetscInt j,PetscInt M)
8 {
9 PetscHash_t h = PetscHashCombine(PetscHashInt(i),PetscHashInt(j));
10 return (PetscReal) ((h % 5 == 0) ? (1 + i + j*M) : 0);
11 }
12
CheckValuesAIJ(Mat A)13 static PetscErrorCode CheckValuesAIJ(Mat A)
14 {
15 PetscInt M,N,rstart,rend,i,j;
16 PetscReal v,w;
17 PetscScalar val;
18 PetscErrorCode ierr;
19
20 PetscFunctionBegin;
21 ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr);
22 ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr);
23 for (i=rstart; i<rend; i++) {
24 for (j=0; j<N; j++) {
25 ierr = MatGetValue(A,i,j,&val);CHKERRQ(ierr);
26 v = MakeValue(i,j,M); w = PetscRealPart(val);
27 if (PetscAbsReal(v-w) > 0) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Matrix entry (%D,%D) should be %g, got %g",i,j,(double)v,(double)w);
28 }
29 }
30 PetscFunctionReturn(0);
31 }
32
main(int argc,char ** args)33 int main(int argc,char **args)
34 {
35 Mat A;
36 PetscInt M = 24,N = 48,bs = 2;
37 PetscInt rstart,rend,i,j;
38 PetscErrorCode ierr;
39 PetscViewer view;
40
41 ierr = PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr;
42 /*
43 Create a parallel BAIJ matrix shared by all processors
44 */
45 ierr = MatCreateBAIJ(PETSC_COMM_WORLD,
46 bs,
47 PETSC_DECIDE,PETSC_DECIDE,
48 M,N,
49 PETSC_DECIDE,NULL,
50 PETSC_DECIDE,NULL,
51 &A);CHKERRQ(ierr);
52
53 /*
54 Set values into the matrix
55 */
56 ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr);
57 ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr);
58 for (i=rstart; i<rend; i++) {
59 for (j=0; j<N; j++) {
60 PetscReal v = MakeValue(i,j,M);
61 if (PetscAbsReal(v) > 0) {
62 ierr = MatSetValue(A,i,j,v,INSERT_VALUES);CHKERRQ(ierr);
63 }
64 }
65 }
66 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
67 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
68 ierr = MatViewFromOptions(A,NULL,"-mat_base_view");CHKERRQ(ierr);
69
70 /*
71 Store the binary matrix to a file
72 */
73 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);CHKERRQ(ierr);
74 for (i=0; i<3; i++) {
75 ierr = MatView(A,view);CHKERRQ(ierr);
76 }
77 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
78 ierr = MatDestroy(&A);CHKERRQ(ierr);
79
80 /*
81 Now reload the matrix and check its values
82 */
83 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr);
84 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
85 ierr = MatSetType(A,MATBAIJ);CHKERRQ(ierr);
86 for (i=0; i<3; i++) {
87 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);}
88 ierr = MatLoad(A,view);CHKERRQ(ierr);
89 ierr = CheckValuesAIJ(A);CHKERRQ(ierr);
90 }
91 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
92 ierr = MatViewFromOptions(A,NULL,"-mat_load_view");CHKERRQ(ierr);
93 ierr = MatDestroy(&A);CHKERRQ(ierr);
94
95 /*
96 Reload in SEQBAIJ matrix and check its values
97 */
98 ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr);
99 ierr = MatCreate(PETSC_COMM_SELF,&A);CHKERRQ(ierr);
100 ierr = MatSetType(A,MATSEQBAIJ);CHKERRQ(ierr);
101 for (i=0; i<3; i++) {
102 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);}
103 ierr = MatLoad(A,view);CHKERRQ(ierr);
104 ierr = CheckValuesAIJ(A);CHKERRQ(ierr);
105 }
106 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
107 ierr = MatDestroy(&A);CHKERRQ(ierr);
108
109 /*
110 Reload in MPIBAIJ matrix and check its values
111 */
112 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr);
113 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
114 ierr = MatSetType(A,MATMPIBAIJ);CHKERRQ(ierr);
115 for (i=0; i<3; i++) {
116 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);}
117 ierr = MatLoad(A,view);CHKERRQ(ierr);
118 ierr = CheckValuesAIJ(A);CHKERRQ(ierr);
119 }
120 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
121 ierr = MatDestroy(&A);CHKERRQ(ierr);
122
123 ierr = PetscFinalize();
124 return ierr;
125 }
126
127 /*TEST
128
129 testset:
130 args: -viewer_binary_mpiio 0
131 output_file: output/ex45.out
132 test:
133 suffix: stdio_1
134 nsize: 1
135 test:
136 suffix: stdio_2
137 nsize: 2
138 test:
139 suffix: stdio_3
140 nsize: 3
141 test:
142 suffix: stdio_4
143 nsize: 4
144 test:
145 suffix: stdio_5
146 nsize: 4
147
148 testset:
149 requires: mpiio
150 args: -viewer_binary_mpiio 1
151 output_file: output/ex45.out
152 test:
153 suffix: mpiio_1
154 nsize: 1
155 test:
156 suffix: mpiio_2
157 nsize: 2
158 test:
159 suffix: mpiio_3
160 nsize: 3
161 test:
162 suffix: mpiio_4
163 nsize: 4
164 test:
165 suffix: mpiio_5
166 nsize: 5
167
168 TEST*/
169