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