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:   MatGetSize(A, &M, &N);
 20:   MatGetOwnershipRange(A, &rstart, &rend);
 21:   for (i = rstart; i < rend; i++) {
 22:     for (j = 0; j < N; j++) {
 23:       MatGetValue(A, i, j, &val);
 24:       v = MakeValue(i, j, M);
 25:       w = PetscRealPart(val);
 27:     }
 28:   }
 29:   return 0;
 30: }

 32: int main(int argc, char **args)
 33: {
 34:   Mat         A;
 35:   PetscInt    M = 11, N = 13;
 36:   PetscInt    rstart, rend, i, j;
 37:   PetscViewer view;

 40:   PetscInitialize(&argc, &args, NULL, help);
 41:   /*
 42:       Create a parallel AIJ matrix shared by all processors
 43:   */
 44:   MatCreateAIJ(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, M, N, PETSC_DECIDE, NULL, PETSC_DECIDE, NULL, &A);

 46:   /*
 47:       Set values into the matrix
 48:   */
 49:   MatGetOwnershipRange(A, &rstart, &rend);
 50:   for (i = rstart; i < rend; i++) {
 51:     for (j = 0; j < N; j++) {
 52:       PetscReal v = MakeValue(i, j, M);
 53:       if (PetscAbsReal(v) > 0) MatSetValue(A, i, j, v, INSERT_VALUES);
 54:     }
 55:   }
 56:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 57:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
 58:   MatViewFromOptions(A, NULL, "-mat_base_view");

 60:   /*
 61:       Store the binary matrix to a file
 62:   */
 63:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);
 64:   for (i = 0; i < 3; i++) MatView(A, view);
 65:   PetscViewerDestroy(&view);
 66:   MatDestroy(&A);

 68:   /*
 69:       Now reload the matrix and check its values
 70:   */
 71:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view);
 72:   MatCreate(PETSC_COMM_WORLD, &A);
 73:   MatSetType(A, MATAIJ);
 74:   for (i = 0; i < 3; i++) {
 75:     if (i > 0) MatZeroEntries(A);
 76:     MatLoad(A, view);
 77:     CheckValuesAIJ(A);
 78:   }
 79:   PetscViewerDestroy(&view);
 80:   MatViewFromOptions(A, NULL, "-mat_load_view");
 81:   MatDestroy(&A);

 83:   /*
 84:       Reload in SEQAIJ matrix and check its values
 85:   */
 86:   PetscViewerBinaryOpen(PETSC_COMM_SELF, "matrix.dat", FILE_MODE_READ, &view);
 87:   MatCreate(PETSC_COMM_SELF, &A);
 88:   MatSetType(A, MATSEQAIJ);
 89:   for (i = 0; i < 3; i++) {
 90:     if (i > 0) MatZeroEntries(A);
 91:     MatLoad(A, view);
 92:     CheckValuesAIJ(A);
 93:   }
 94:   PetscViewerDestroy(&view);
 95:   MatDestroy(&A);

 97:   /*
 98:      Reload in MPIAIJ matrix and check its values
 99:   */
100:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view);
101:   MatCreate(PETSC_COMM_WORLD, &A);
102:   MatSetType(A, MATMPIAIJ);
103:   for (i = 0; i < 3; i++) {
104:     if (i > 0) MatZeroEntries(A);
105:     MatLoad(A, view);
106:     CheckValuesAIJ(A);
107:   }
108:   PetscViewerDestroy(&view);
109:   MatDestroy(&A);

111:   PetscFinalize();
112:   return 0;
113: }

115: /*TEST

117:    testset:
118:       args: -viewer_binary_mpiio 0
119:       output_file: output/ex44.out
120:       test:
121:         suffix: stdio_1
122:         nsize: 1
123:       test:
124:         suffix: stdio_2
125:         nsize: 2
126:       test:
127:         suffix: stdio_3
128:         nsize: 3
129:       test:
130:         suffix: stdio_4
131:         nsize: 4
132:       test:
133:         suffix: stdio_15
134:         nsize: 15

136:    testset:
137:       requires: mpiio
138:       args: -viewer_binary_mpiio 1
139:       output_file: output/ex44.out
140:       test:
141:         suffix: mpiio_1
142:         nsize: 1
143:       test:
144:         suffix: mpiio_2
145:         nsize: 2
146:       test:
147:         suffix: mpiio_3
148:         nsize: 3
149:       test:
150:         suffix: mpiio_4
151:         nsize: 4
152:       test:
153:         suffix: mpiio_15
154:         nsize: 15

156: TEST*/