1 static char help[] = "Tests MatView()/MatLoad() with binary viewers for BAIJ matrices.\n\n"; 2 3 #include <petscmat.h> 4 #include <petscviewer.h> 5 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 } 12 13 static PetscErrorCode CheckValuesAIJ(Mat A) 14 { 15 PetscInt M, N, rstart, rend, i, j; 16 PetscReal v, w; 17 PetscScalar val; 18 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 } 32 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; 39 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)); 46 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")); 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 < 3; i++) PetscCall(MatView(A, view)); 67 PetscCall(PetscViewerDestroy(&view)); 68 PetscCall(MatDestroy(&A)); 69 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)); 84 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)); 98 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)); 112 113 PetscCall(PetscFinalize()); 114 return 0; 115 } 116 117 /*TEST 118 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 137 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 157 158 TEST*/ 159