Actual source code: ex44.c
1: static char help[] = "Tests MatView()/MatLoad() with binary viewers for AIJ matrices.\n\n";
3: #include <petscmat.h>
4: #include <petscviewer.h>
6: #include <petsc/private/hashtable.h>
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: }
13: static PetscErrorCode CheckValuesAIJ(Mat A)
14: {
15: PetscInt M, N, rstart, rend, i, j;
16: PetscReal v, w;
17: PetscScalar val;
19: PetscFunctionBegin;
20: PetscCall(MatGetSize(A, &M, &N));
21: PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
22: for (i = rstart; i < rend; i++) {
23: for (j = 0; j < N; j++) {
24: PetscCall(MatGetValue(A, i, j, &val));
25: v = MakeValue(i, j, M);
26: w = PetscRealPart(val);
27: 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);
28: }
29: }
30: PetscFunctionReturn(PETSC_SUCCESS);
31: }
33: int main(int argc, char **args)
34: {
35: Mat A;
36: PetscInt M = 11, N = 13;
37: PetscInt rstart, rend, i, j;
38: PetscViewer view;
40: PetscFunctionBeginUser;
41: PetscCall(PetscInitialize(&argc, &args, NULL, help));
42: /*
43: Create a parallel AIJ matrix shared by all processors
44: */
45: PetscCall(MatCreateAIJ(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, M, N, PETSC_DECIDE, NULL, PETSC_DECIDE, NULL, &A));
47: /*
48: Set values into the matrix
49: */
50: PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
51: for (i = rstart; i < rend; i++) {
52: for (j = 0; j < N; j++) {
53: PetscReal v = MakeValue(i, j, M);
54: if (PetscAbsReal(v) > 0) PetscCall(MatSetValue(A, i, j, v, INSERT_VALUES));
55: }
56: }
57: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
58: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
59: PetscCall(MatViewFromOptions(A, NULL, "-mat_base_view"));
61: /*
62: Store the binary matrix to a file
63: */
64: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view));
65: for (i = 0; i < 3; i++) PetscCall(MatView(A, view));
66: PetscCall(PetscViewerDestroy(&view));
67: PetscCall(MatDestroy(&A));
69: /*
70: Now reload the matrix and check its values
71: */
72: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view));
73: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
74: PetscCall(MatSetType(A, MATAIJ));
75: for (i = 0; i < 3; i++) {
76: if (i > 0) PetscCall(MatZeroEntries(A));
77: PetscCall(MatLoad(A, view));
78: PetscCall(CheckValuesAIJ(A));
79: }
80: PetscCall(PetscViewerDestroy(&view));
81: PetscCall(MatViewFromOptions(A, NULL, "-mat_load_view"));
82: PetscCall(MatDestroy(&A));
84: /*
85: Reload in SEQAIJ matrix and check its values
86: */
87: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF, "matrix.dat", FILE_MODE_READ, &view));
88: PetscCall(MatCreate(PETSC_COMM_SELF, &A));
89: PetscCall(MatSetType(A, MATSEQAIJ));
90: for (i = 0; i < 3; i++) {
91: if (i > 0) PetscCall(MatZeroEntries(A));
92: PetscCall(MatLoad(A, view));
93: PetscCall(CheckValuesAIJ(A));
94: }
95: PetscCall(PetscViewerDestroy(&view));
96: PetscCall(MatDestroy(&A));
98: /*
99: Reload in MPIAIJ matrix and check its values
100: */
101: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view));
102: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
103: PetscCall(MatSetType(A, MATMPIAIJ));
104: for (i = 0; i < 3; i++) {
105: if (i > 0) PetscCall(MatZeroEntries(A));
106: PetscCall(MatLoad(A, view));
107: PetscCall(CheckValuesAIJ(A));
108: }
109: PetscCall(PetscViewerDestroy(&view));
110: PetscCall(MatDestroy(&A));
112: PetscCall(PetscFinalize());
113: return 0;
114: }
116: /*TEST
118: testset:
119: args: -viewer_binary_mpiio 0
120: output_file: output/ex44.out
121: test:
122: suffix: stdio_1
123: nsize: 1
124: test:
125: suffix: stdio_2
126: nsize: 2
127: test:
128: suffix: stdio_3
129: nsize: 3
130: test:
131: suffix: stdio_4
132: nsize: 4
133: test:
134: suffix: stdio_15
135: nsize: 15
137: testset:
138: requires: mpiio
139: args: -viewer_binary_mpiio 1
140: output_file: output/ex44.out
141: test:
142: suffix: mpiio_1
143: nsize: 1
144: test:
145: suffix: mpiio_2
146: nsize: 2
147: test:
148: suffix: mpiio_3
149: nsize: 3
150: test:
151: suffix: mpiio_4
152: nsize: 4
153: test:
154: suffix: mpiio_15
155: nsize: 15
157: TEST*/