1 static char help[] = "Tests MatDenseGetArray() and MatView()/MatLoad() with binary viewers.\n\n"; 2 3 #include <petscmat.h> 4 #include <petscviewer.h> 5 6 static PetscErrorCode CheckValues(Mat A, PetscBool one) 7 { 8 const PetscScalar *array; 9 PetscInt M, N, rstart, rend, lda, i, j; 10 11 PetscFunctionBegin; 12 PetscCall(MatDenseGetArrayRead(A, &array)); 13 PetscCall(MatDenseGetLDA(A, &lda)); 14 PetscCall(MatGetSize(A, &M, &N)); 15 PetscCall(MatGetOwnershipRange(A, &rstart, &rend)); 16 for (i = rstart; i < rend; i++) { 17 for (j = 0; j < N; j++) { 18 PetscInt ii = i - rstart, jj = j; 19 PetscReal v = (PetscReal)(one ? 1 : (1 + i + j * M)); 20 PetscReal w = PetscRealPart(array[ii + jj * lda]); 21 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); 22 } 23 } 24 PetscCall(MatDenseRestoreArrayRead(A, &array)); 25 PetscFunctionReturn(PETSC_SUCCESS); 26 } 27 28 #define CheckValuesIJ(A) CheckValues(A, PETSC_FALSE) 29 #define CheckValuesOne(A) CheckValues(A, PETSC_TRUE) 30 31 int main(int argc, char **args) 32 { 33 Mat A; 34 PetscInt i, j, M = 4, N = 3, rstart, rend; 35 PetscScalar *array; 36 char mattype[256]; 37 PetscViewer view; 38 39 PetscFunctionBeginUser; 40 PetscCall(PetscInitialize(&argc, &args, NULL, help)); 41 PetscCall(PetscStrncpy(mattype, MATMPIDENSE, sizeof(mattype))); 42 PetscCall(PetscOptionsGetString(NULL, NULL, "-mat_type", mattype, sizeof(mattype), NULL)); 43 /* 44 Create a parallel dense matrix shared by all processors 45 */ 46 PetscCall(MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, M, N, NULL, &A)); 47 PetscCall(MatConvert(A, mattype, MAT_INPLACE_MATRIX, &A)); 48 /* 49 Set values into the matrix 50 */ 51 for (i = 0; i < M; i++) { 52 for (j = 0; j < N; j++) { 53 PetscScalar v = (PetscReal)(1 + i + j * M); 54 PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES)); 55 } 56 } 57 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 58 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 59 PetscCall(MatScale(A, 2.0)); 60 PetscCall(MatScale(A, 1.0 / 2.0)); 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 < 2; i++) { 67 PetscCall(MatView(A, view)); 68 PetscCall(PetscViewerPushFormat(view, PETSC_VIEWER_NATIVE)); 69 PetscCall(MatView(A, view)); 70 PetscCall(PetscViewerPopFormat(view)); 71 } 72 PetscCall(PetscViewerDestroy(&view)); 73 PetscCall(MatDestroy(&A)); 74 75 /* 76 Now reload the matrix and check its values 77 */ 78 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view)); 79 PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 80 PetscCall(MatSetType(A, mattype)); 81 for (i = 0; i < 4; i++) { 82 if (i > 0) PetscCall(MatZeroEntries(A)); 83 PetscCall(MatLoad(A, view)); 84 PetscCall(CheckValuesIJ(A)); 85 } 86 PetscCall(PetscViewerDestroy(&view)); 87 88 PetscCall(MatGetOwnershipRange(A, &rstart, &rend)); 89 PetscCall(PetscMalloc1((rend - rstart) * N, &array)); 90 for (i = 0; i < (rend - rstart) * N; i++) array[i] = (PetscReal)1; 91 PetscCall(MatDensePlaceArray(A, array)); 92 PetscCall(MatScale(A, 2.0)); 93 PetscCall(MatScale(A, 1.0 / 2.0)); 94 PetscCall(CheckValuesOne(A)); 95 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view)); 96 PetscCall(MatView(A, view)); 97 PetscCall(MatDenseResetArray(A)); 98 PetscCall(PetscFree(array)); 99 PetscCall(CheckValuesIJ(A)); 100 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_TRUE)); 101 PetscCall(MatView(A, view)); 102 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_FALSE)); 103 PetscCall(PetscViewerDestroy(&view)); 104 PetscCall(MatDestroy(&A)); 105 106 PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 107 PetscCall(MatSetType(A, mattype)); 108 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view)); 109 PetscCall(MatLoad(A, view)); 110 PetscCall(CheckValuesOne(A)); 111 PetscCall(MatZeroEntries(A)); 112 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_TRUE)); 113 PetscCall(MatLoad(A, view)); 114 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_FALSE)); 115 PetscCall(CheckValuesIJ(A)); 116 PetscCall(PetscViewerDestroy(&view)); 117 PetscCall(MatDestroy(&A)); 118 119 { 120 PetscInt m = PETSC_DECIDE, n = PETSC_DECIDE; 121 PetscCall(PetscSplitOwnership(PETSC_COMM_WORLD, &m, &M)); 122 PetscCall(PetscSplitOwnership(PETSC_COMM_WORLD, &n, &N)); 123 /* TODO: MatCreateDense requires data!=NULL at all processes! */ 124 PetscCall(PetscMalloc1(m * N + 1, &array)); 125 126 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_READ, &view)); 127 PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, n, M, N, array, &A)); 128 PetscCall(MatLoad(A, view)); 129 PetscCall(CheckValuesOne(A)); 130 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_TRUE)); 131 PetscCall(MatLoad(A, view)); 132 PetscCall(PetscViewerBinarySetSkipHeader(view, PETSC_FALSE)); 133 PetscCall(CheckValuesIJ(A)); 134 PetscCall(MatDestroy(&A)); 135 PetscCall(PetscViewerDestroy(&view)); 136 137 PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, n, M, N, array, &A)); 138 PetscCall(CheckValuesIJ(A)); 139 PetscCall(MatDestroy(&A)); 140 141 PetscCall(PetscFree(array)); 142 } 143 144 PetscCall(PetscFinalize()); 145 return 0; 146 } 147 148 /*TEST 149 150 testset: 151 args: -viewer_binary_mpiio 0 152 output_file: output/ex16.out 153 test: 154 suffix: stdio_1 155 nsize: 1 156 args: -mat_type seqdense 157 test: 158 suffix: stdio_2 159 nsize: 2 160 test: 161 suffix: stdio_3 162 nsize: 3 163 test: 164 suffix: stdio_4 165 nsize: 4 166 test: 167 suffix: stdio_5 168 nsize: 5 169 test: 170 requires: cuda 171 args: -mat_type seqdensecuda 172 suffix: stdio_cuda_1 173 nsize: 1 174 test: 175 requires: cuda 176 args: -mat_type mpidensecuda 177 suffix: stdio_cuda_2 178 nsize: 2 179 test: 180 requires: cuda 181 args: -mat_type mpidensecuda 182 suffix: stdio_cuda_3 183 nsize: 3 184 test: 185 requires: cuda 186 args: -mat_type mpidensecuda 187 suffix: stdio_cuda_4 188 nsize: 4 189 test: 190 requires: cuda 191 args: -mat_type mpidensecuda 192 suffix: stdio_cuda_5 193 nsize: 5 194 195 testset: 196 requires: mpiio 197 args: -viewer_binary_mpiio 1 198 output_file: output/ex16.out 199 test: 200 suffix: mpiio_1 201 nsize: 1 202 test: 203 suffix: mpiio_2 204 nsize: 2 205 test: 206 suffix: mpiio_3 207 nsize: 3 208 test: 209 suffix: mpiio_4 210 nsize: 4 211 test: 212 suffix: mpiio_5 213 nsize: 5 214 test: 215 requires: cuda 216 args: -mat_type mpidensecuda 217 suffix: mpiio_cuda_1 218 nsize: 1 219 test: 220 requires: cuda 221 args: -mat_type mpidensecuda 222 suffix: mpiio_cuda_2 223 nsize: 2 224 test: 225 requires: cuda 226 args: -mat_type mpidensecuda 227 suffix: mpiio_cuda_3 228 nsize: 3 229 test: 230 requires: cuda 231 args: -mat_type mpidensecuda 232 suffix: mpiio_cuda_4 233 nsize: 4 234 test: 235 requires: cuda 236 args: -mat_type mpidensecuda 237 suffix: mpiio_cuda_5 238 nsize: 5 239 240 TEST*/ 241