1 static char help[] = "Tests MatDenseGetArray() and MatView()/MatLoad() with binary viewers.\n\n";
2
3 #include <petscmat.h>
4 #include <petscviewer.h>
5
CheckValues(Mat A,PetscBool one)6 static PetscErrorCode CheckValues(Mat A, PetscBool one)
7 {
8 const PetscScalar *array;
9 PetscInt M, N, rstart, rend, lda, i, j;
10
11 PetscFunctionBegin;
12 PetscCall(MatDenseGetArrayRead(A, &array));
13 PetscCall(MatDenseGetLDA(A, &lda));
14 PetscCall(MatGetSize(A, &M, &N));
15 PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
16 for (i = rstart; i < rend; i++) {
17 for (j = 0; j < N; j++) {
18 PetscInt ii = i - rstart, jj = j;
19 PetscReal v = (PetscReal)(one ? 1 : (1 + i + j * M));
20 PetscReal w = PetscRealPart(array[ii + jj * lda]);
21 PetscCheck(PetscAbsReal(v - w) <= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Matrix entry (%" PetscInt_FMT ",%" PetscInt_FMT ") should be %g, got %g", i, j, (double)v, (double)w);
22 }
23 }
24 PetscCall(MatDenseRestoreArrayRead(A, &array));
25 PetscFunctionReturn(PETSC_SUCCESS);
26 }
27
28 #define CheckValuesIJ(A) CheckValues(A, PETSC_FALSE)
29 #define CheckValuesOne(A) CheckValues(A, PETSC_TRUE)
30
main(int argc,char ** args)31 int main(int argc, char **args)
32 {
33 Mat A;
34 PetscInt i, j, M = 4, N = 3, rstart, rend;
35 PetscScalar *array;
36 char mattype[256];
37 PetscViewer view;
38
39 PetscFunctionBeginUser;
40 PetscCall(PetscInitialize(&argc, &args, NULL, help));
41 PetscCall(PetscStrncpy(mattype, MATMPIDENSE, sizeof(mattype)));
42 PetscCall(PetscOptionsGetString(NULL, NULL, "-mat_type", mattype, sizeof(mattype), NULL));
43 /*
44 Create a parallel dense matrix shared by all processors
45 */
46 PetscCall(MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, M, N, NULL, &A));
47 PetscCall(MatConvert(A, mattype, MAT_INPLACE_MATRIX, &A));
48 /*
49 Set values into the matrix
50 */
51 for (i = 0; i < M; i++) {
52 for (j = 0; j < N; j++) {
53 PetscScalar v = (PetscReal)(1 + i + j * M);
54 PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES));
55 }
56 }
57 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
58 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
59 PetscCall(MatScale(A, 2.0));
60 PetscCall(MatScale(A, 1.0 / 2.0));
61
62 /*
63 Store the binary matrix to a file
64 */
65 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view));
66 for (i = 0; i < 2; i++) {
67 PetscCall(MatView(A, view));
68 PetscCall(PetscViewerPushFormat(view, PETSC_VIEWER_NATIVE));
69 PetscCall(MatView(A, view));
70 PetscCall(PetscViewerPopFormat(view));
71 }
72 PetscCall(PetscViewerDestroy(&view));
73 PetscCall(MatDestroy(&A));
74
75 /*
76 Now reload the matrix and check its values
77 */
78 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view));
79 PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
80 PetscCall(MatSetType(A, mattype));
81 for (i = 0; i < 4; i++) {
82 if (i > 0) PetscCall(MatZeroEntries(A));
83 PetscCall(MatLoad(A, view));
84 PetscCall(CheckValuesIJ(A));
85 }
86 PetscCall(PetscViewerDestroy(&view));
87
88 PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
89 PetscCall(PetscMalloc1((rend - rstart) * N, &array));
90 for (i = 0; i < (rend - rstart) * N; i++) array[i] = (PetscReal)1;
91 PetscCall(MatDensePlaceArray(A, array));
92 PetscCall(MatScale(A, 2.0));
93 PetscCall(MatScale(A, 1.0 / 2.0));
94 PetscCall(CheckValuesOne(A));
95 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view));
96 PetscCall(MatView(A, view));
97 PetscCall(MatDenseResetArray(A));
98 PetscCall(PetscFree(array));
99 PetscCall(CheckValuesIJ(A));
100 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_TRUE));
101 PetscCall(MatView(A, view));
102 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_FALSE));
103 PetscCall(PetscViewerDestroy(&view));
104 PetscCall(MatDestroy(&A));
105
106 PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
107 PetscCall(MatSetType(A, mattype));
108 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view));
109 PetscCall(MatLoad(A, view));
110 PetscCall(CheckValuesOne(A));
111 PetscCall(MatZeroEntries(A));
112 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_TRUE));
113 PetscCall(MatLoad(A, view));
114 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_FALSE));
115 PetscCall(CheckValuesIJ(A));
116 PetscCall(PetscViewerDestroy(&view));
117 PetscCall(MatDestroy(&A));
118
119 {
120 PetscInt m = PETSC_DECIDE, n = PETSC_DECIDE;
121 PetscCall(PetscSplitOwnership(PETSC_COMM_WORLD, &m, &M));
122 PetscCall(PetscSplitOwnership(PETSC_COMM_WORLD, &n, &N));
123 /* TODO: MatCreateDense requires data!=NULL at all processes! */
124 PetscCall(PetscMalloc1(m * N + 1, &array));
125
126 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view));
127 PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, n, M, N, array, &A));
128 PetscCall(MatLoad(A, view));
129 PetscCall(CheckValuesOne(A));
130 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_TRUE));
131 PetscCall(MatLoad(A, view));
132 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_FALSE));
133 PetscCall(CheckValuesIJ(A));
134 PetscCall(MatDestroy(&A));
135 PetscCall(PetscViewerDestroy(&view));
136
137 PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, n, M, N, array, &A));
138 PetscCall(CheckValuesIJ(A));
139 PetscCall(MatDestroy(&A));
140
141 PetscCall(PetscFree(array));
142 }
143
144 PetscCall(PetscFinalize());
145 return 0;
146 }
147
148 /*TEST
149
150 testset:
151 args: -viewer_binary_mpiio 0
152 output_file: output/empty.out
153 test:
154 suffix: stdio_1
155 nsize: 1
156 args: -mat_type seqdense
157 test:
158 suffix: stdio_2
159 nsize: 2
160 test:
161 suffix: stdio_3
162 nsize: 3
163 test:
164 suffix: stdio_4
165 nsize: 4
166 test:
167 suffix: stdio_5
168 nsize: 5
169 test:
170 requires: cuda
171 args: -mat_type seqdensecuda
172 suffix: stdio_cuda_1
173 nsize: 1
174 test:
175 requires: cuda
176 args: -mat_type mpidensecuda
177 suffix: stdio_cuda_2
178 nsize: 2
179 test:
180 requires: cuda
181 args: -mat_type mpidensecuda
182 suffix: stdio_cuda_3
183 nsize: 3
184 test:
185 requires: cuda
186 args: -mat_type mpidensecuda
187 suffix: stdio_cuda_4
188 nsize: 4
189 test:
190 requires: cuda
191 args: -mat_type mpidensecuda
192 suffix: stdio_cuda_5
193 nsize: 5
194
195 testset:
196 requires: mpiio
197 args: -viewer_binary_mpiio 1
198 output_file: output/empty.out
199 test:
200 suffix: mpiio_1
201 nsize: 1
202 test:
203 suffix: mpiio_2
204 nsize: 2
205 test:
206 suffix: mpiio_3
207 nsize: 3
208 test:
209 suffix: mpiio_4
210 nsize: 4
211 test:
212 suffix: mpiio_5
213 nsize: 5
214 test:
215 requires: cuda
216 args: -mat_type mpidensecuda
217 suffix: mpiio_cuda_1
218 nsize: 1
219 test:
220 requires: cuda
221 args: -mat_type mpidensecuda
222 suffix: mpiio_cuda_2
223 nsize: 2
224 test:
225 requires: cuda
226 args: -mat_type mpidensecuda
227 suffix: mpiio_cuda_3
228 nsize: 3
229 test:
230 requires: cuda
231 args: -mat_type mpidensecuda
232 suffix: mpiio_cuda_4
233 nsize: 4
234 test:
235 requires: cuda
236 args: -mat_type mpidensecuda
237 suffix: mpiio_cuda_5
238 nsize: 5
239
240 TEST*/
241