Actual source code: ex45.c
1: static char help[] = "Tests MatView()/MatLoad() with binary viewers for BAIJ 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 = 24, N = 48, bs = 2;
37: PetscInt rstart, rend, i, j;
38: PetscViewer view;
40: PetscFunctionBeginUser;
41: PetscCall(PetscInitialize(&argc, &args, NULL, help));
42: /*
43: Create a parallel BAIJ matrix shared by all processors
44: */
45: PetscCall(MatCreateBAIJ(PETSC_COMM_WORLD, bs, PETSC_DECIDE, PETSC_DECIDE, M, N, PETSC_DECIDE, NULL, PETSC_DECIDE, NULL, &A));
47: /*
48: Set values into the matrix
49: */
50: PetscCall(MatGetSize(A, &M, &N));
51: PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
52: for (i = rstart; i < rend; i++) {
53: for (j = 0; j < N; j++) {
54: PetscReal v = MakeValue(i, j, M);
55: if (PetscAbsReal(v) > 0) PetscCall(MatSetValue(A, i, j, v, INSERT_VALUES));
56: }
57: }
58: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
59: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
60: PetscCall(MatViewFromOptions(A, NULL, "-mat_base_view"));
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 < 3; i++) PetscCall(MatView(A, view));
67: PetscCall(PetscViewerDestroy(&view));
68: PetscCall(MatDestroy(&A));
70: /*
71: Now reload the matrix and check its values
72: */
73: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view));
74: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
75: PetscCall(MatSetType(A, MATBAIJ));
76: for (i = 0; i < 3; i++) {
77: if (i > 0) PetscCall(MatZeroEntries(A));
78: PetscCall(MatLoad(A, view));
79: PetscCall(CheckValuesAIJ(A));
80: }
81: PetscCall(PetscViewerDestroy(&view));
82: PetscCall(MatViewFromOptions(A, NULL, "-mat_load_view"));
83: PetscCall(MatDestroy(&A));
85: /*
86: Reload in SEQBAIJ matrix and check its values
87: */
88: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF, "matrix.dat", FILE_MODE_READ, &view));
89: PetscCall(MatCreate(PETSC_COMM_SELF, &A));
90: PetscCall(MatSetType(A, MATSEQBAIJ));
91: for (i = 0; i < 3; i++) {
92: if (i > 0) PetscCall(MatZeroEntries(A));
93: PetscCall(MatLoad(A, view));
94: PetscCall(CheckValuesAIJ(A));
95: }
96: PetscCall(PetscViewerDestroy(&view));
97: PetscCall(MatDestroy(&A));
99: /*
100: Reload in MPIBAIJ matrix and check its values
101: */
102: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view));
103: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
104: PetscCall(MatSetType(A, MATMPIBAIJ));
105: for (i = 0; i < 3; i++) {
106: if (i > 0) PetscCall(MatZeroEntries(A));
107: PetscCall(MatLoad(A, view));
108: PetscCall(CheckValuesAIJ(A));
109: }
110: PetscCall(PetscViewerDestroy(&view));
111: PetscCall(MatDestroy(&A));
113: PetscCall(PetscFinalize());
114: return 0;
115: }
117: /*TEST
119: testset:
120: args: -viewer_binary_mpiio 0
121: output_file: output/ex45.out
122: test:
123: suffix: stdio_1
124: nsize: 1
125: test:
126: suffix: stdio_2
127: nsize: 2
128: test:
129: suffix: stdio_3
130: nsize: 3
131: test:
132: suffix: stdio_4
133: nsize: 4
134: test:
135: suffix: stdio_5
136: nsize: 4
138: testset:
139: requires: mpiio
140: args: -viewer_binary_mpiio 1
141: output_file: output/ex45.out
142: test:
143: suffix: mpiio_1
144: nsize: 1
145: test:
146: suffix: mpiio_2
147: nsize: 2
148: test:
149: suffix: mpiio_3
150: nsize: 3
151: test:
152: suffix: mpiio_4
153: nsize: 4
154: test:
155: suffix: mpiio_5
156: nsize: 5
158: TEST*/