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 PetscViewer view; 39 40 ierr = PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr; 41 /* 42 Create a parallel dense matrix shared by all processors 43 */ 44 ierr = MatCreateDense(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,M,N,NULL,&A);CHKERRQ(ierr); 45 46 /* 47 Set values into the matrix 48 */ 49 for (i=0; i<M; i++) { 50 for (j=0; j<N; j++) { 51 PetscScalar v = (PetscReal)(1 + i + j*M); 52 ierr = MatSetValues(A,1,&i,1,&j,&v,INSERT_VALUES);CHKERRQ(ierr); 53 } 54 } 55 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 56 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 57 58 /* 59 Store the binary matrix to a file 60 */ 61 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);CHKERRQ(ierr); 62 for (i=0; i<2; i++) { 63 ierr = MatView(A,view);CHKERRQ(ierr); 64 ierr = PetscViewerPushFormat(view,PETSC_VIEWER_NATIVE);CHKERRQ(ierr); 65 ierr = MatView(A,view);CHKERRQ(ierr); 66 ierr = PetscViewerPopFormat(view);CHKERRQ(ierr); 67 } 68 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 69 ierr = MatDestroy(&A);CHKERRQ(ierr); 70 71 /* 72 Now reload the matrix and check its values 73 */ 74 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 75 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 76 ierr = MatSetType(A,MATMPIDENSE);CHKERRQ(ierr); 77 for (i=0; i<4; i++) { 78 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 79 ierr = MatLoad(A,view);CHKERRQ(ierr); 80 ierr = CheckValuesIJ(A);CHKERRQ(ierr); 81 } 82 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 83 84 ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr); 85 ierr = PetscMalloc1((rend-rstart)*N,&array);CHKERRQ(ierr); 86 for (i=0; i<(rend-rstart)*N; i++) array[i] = (PetscReal)1; 87 ierr = MatDensePlaceArray(A,array);CHKERRQ(ierr); 88 ierr = CheckValuesOne(A);CHKERRQ(ierr); 89 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_WRITE,&view);CHKERRQ(ierr); 90 ierr = MatView(A,view);CHKERRQ(ierr); 91 ierr = MatDenseResetArray(A);CHKERRQ(ierr); 92 ierr = PetscFree(array);CHKERRQ(ierr); 93 ierr = CheckValuesIJ(A);CHKERRQ(ierr); 94 ierr = PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);CHKERRQ(ierr); 95 ierr = MatView(A,view);CHKERRQ(ierr); 96 ierr = PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);CHKERRQ(ierr); 97 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 98 ierr = MatDestroy(&A);CHKERRQ(ierr); 99 100 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 101 ierr = MatSetType(A,MATMPIDENSE);CHKERRQ(ierr); 102 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 103 ierr = MatLoad(A,view);CHKERRQ(ierr); 104 ierr = CheckValuesOne(A);CHKERRQ(ierr); 105 ierr = MatZeroEntries(A);CHKERRQ(ierr); 106 ierr = PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);CHKERRQ(ierr); 107 ierr = MatLoad(A,view);CHKERRQ(ierr); 108 ierr = PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);CHKERRQ(ierr); 109 ierr = CheckValuesIJ(A);CHKERRQ(ierr); 110 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 111 ierr = MatDestroy(&A);CHKERRQ(ierr); 112 113 { 114 PetscInt m = PETSC_DECIDE, n = PETSC_DECIDE; 115 ierr = PetscSplitOwnership(PETSC_COMM_WORLD,&m,&M);CHKERRQ(ierr); 116 ierr = PetscSplitOwnership(PETSC_COMM_WORLD,&n,&N);CHKERRQ(ierr); 117 /* TODO: MatCreateDense requires data!=NULL at all processes! */ 118 ierr = PetscMalloc1(m*N+1,&array);CHKERRQ(ierr); 119 120 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 121 ierr = MatCreateDense(PETSC_COMM_WORLD,m,n,M,N,array,&A);CHKERRQ(ierr); 122 ierr = MatLoad(A,view);CHKERRQ(ierr); 123 ierr = CheckValuesOne(A);CHKERRQ(ierr); 124 ierr = PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);CHKERRQ(ierr); 125 ierr = MatLoad(A,view);CHKERRQ(ierr); 126 ierr = PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);CHKERRQ(ierr); 127 ierr = CheckValuesIJ(A);CHKERRQ(ierr); 128 ierr = MatDestroy(&A);CHKERRQ(ierr); 129 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 130 131 ierr = MatCreateDense(PETSC_COMM_WORLD,m,n,M,N,array,&A);CHKERRQ(ierr); 132 ierr = CheckValuesIJ(A);CHKERRQ(ierr); 133 ierr = MatDestroy(&A);CHKERRQ(ierr); 134 135 ierr = PetscFree(array);CHKERRQ(ierr); 136 } 137 138 ierr = PetscFinalize(); 139 return ierr; 140 } 141 142 143 /*TEST 144 145 testset: 146 args: -viewer_binary_mpiio 0 147 output_file: output/ex16.out 148 test: 149 suffix: stdio_1 150 nsize: 1 151 test: 152 suffix: stdio_2 153 nsize: 2 154 test: 155 suffix: stdio_3 156 nsize: 3 157 test: 158 suffix: stdio_4 159 nsize: 4 160 test: 161 suffix: stdio_5 162 nsize: 5 163 164 testset: 165 requires: mpiio 166 args: -viewer_binary_mpiio 1 167 output_file: output/ex16.out 168 test: 169 suffix: mpiio_1 170 nsize: 1 171 test: 172 suffix: mpiio_2 173 nsize: 2 174 test: 175 suffix: mpiio_3 176 nsize: 3 177 test: 178 suffix: mpiio_4 179 nsize: 4 180 test: 181 suffix: mpiio_5 182 nsize: 5 183 184 185 TEST*/ 186