Actual source code: ex50.c
1: static char help[] = "Tests MatView()/MatLoad() with binary viewers for SBAIJ 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;
18: PetscBool seqsbaij, mpisbaij, sbaij;
20: PetscFunctionBegin;
21: PetscCall(MatGetSize(A, &M, &N));
22: PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
23: PetscCall(PetscObjectTypeCompare((PetscObject)A, MATSEQSBAIJ, &seqsbaij));
24: PetscCall(PetscObjectTypeCompare((PetscObject)A, MATMPISBAIJ, &mpisbaij));
25: sbaij = (seqsbaij || mpisbaij) ? PETSC_TRUE : PETSC_FALSE;
26: for (i = rstart; i < rend; i++) {
27: for (j = (sbaij ? i : 0); j < N; j++) {
28: PetscCall(MatGetValue(A, i, j, &val));
29: v = MakeValue(i, j, M);
30: w = PetscRealPart(val);
31: 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);
32: }
33: }
34: PetscFunctionReturn(PETSC_SUCCESS);
35: }
37: int main(int argc, char **args)
38: {
39: Mat A;
40: PetscInt M = 24, N = 24, bs = 3;
41: PetscInt rstart, rend, i, j;
42: PetscViewer view;
44: PetscFunctionBeginUser;
45: PetscCall(PetscInitialize(&argc, &args, NULL, help));
46: /*
47: Create a parallel SBAIJ matrix shared by all processors
48: */
49: PetscCall(MatCreateSBAIJ(PETSC_COMM_WORLD, bs, PETSC_DECIDE, PETSC_DECIDE, M, N, PETSC_DECIDE, NULL, PETSC_DECIDE, NULL, &A));
51: /*
52: Set values into the matrix
53: */
54: PetscCall(MatGetSize(A, &M, &N));
55: PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
56: for (i = rstart; i < rend; i++) {
57: for (j = 0; j < N; j++) {
58: PetscReal v = MakeValue(i, j, M);
59: if (PetscAbsReal(v) > 0) PetscCall(MatSetValue(A, i, j, v, INSERT_VALUES));
60: }
61: }
62: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
63: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
64: PetscCall(MatViewFromOptions(A, NULL, "-mat_base_view"));
66: /*
67: Store the binary matrix to a file
68: */
69: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view));
70: for (i = 0; i < 3; i++) PetscCall(MatView(A, view));
71: PetscCall(PetscViewerDestroy(&view));
72: PetscCall(MatDestroy(&A));
74: /*
75: Now reload the matrix and check its values
76: */
77: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view));
78: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
79: PetscCall(MatSetType(A, MATSBAIJ));
80: for (i = 0; i < 3; i++) {
81: if (i > 0) PetscCall(MatZeroEntries(A));
82: PetscCall(MatLoad(A, view));
83: PetscCall(CheckValuesAIJ(A));
84: }
85: PetscCall(PetscViewerDestroy(&view));
86: PetscCall(MatViewFromOptions(A, NULL, "-mat_load_view"));
87: PetscCall(MatDestroy(&A));
89: /*
90: Reload in SEQSBAIJ matrix and check its values
91: */
92: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF, "matrix.dat", FILE_MODE_READ, &view));
93: PetscCall(MatCreate(PETSC_COMM_SELF, &A));
94: PetscCall(MatSetType(A, MATSEQSBAIJ));
95: for (i = 0; i < 3; i++) {
96: if (i > 0) PetscCall(MatZeroEntries(A));
97: PetscCall(MatLoad(A, view));
98: PetscCall(CheckValuesAIJ(A));
99: }
100: PetscCall(PetscViewerDestroy(&view));
101: PetscCall(MatDestroy(&A));
103: /*
104: Reload in MPISBAIJ matrix and check its values
105: */
106: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view));
107: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
108: PetscCall(MatSetType(A, MATMPISBAIJ));
109: for (i = 0; i < 3; i++) {
110: if (i > 0) PetscCall(MatZeroEntries(A));
111: PetscCall(MatLoad(A, view));
112: PetscCall(CheckValuesAIJ(A));
113: }
114: PetscCall(PetscViewerDestroy(&view));
115: PetscCall(MatDestroy(&A));
117: PetscCall(PetscFinalize());
118: return 0;
119: }
121: /*TEST
123: testset:
124: args: -viewer_binary_mpiio 0
125: output_file: output/ex50.out
126: test:
127: suffix: stdio_1
128: nsize: 1
129: test:
130: suffix: stdio_2
131: nsize: 2
132: test:
133: suffix: stdio_3
134: nsize: 3
135: test:
136: suffix: stdio_4
137: nsize: 4
138: test:
139: suffix: stdio_5
140: nsize: 4
142: testset:
143: requires: mpiio
144: args: -viewer_binary_mpiio 1
145: output_file: output/ex50.out
146: test:
147: suffix: mpiio_1
148: nsize: 1
149: test:
150: suffix: mpiio_2
151: nsize: 2
152: test:
153: suffix: mpiio_3
154: nsize: 3
155: test:
156: suffix: mpiio_4
157: nsize: 4
158: test:
159: suffix: mpiio_5
160: nsize: 5
162: TEST*/