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