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(0); 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(PetscStrcpy(mattype,MATMPIDENSE)); 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