1 static char help[] = "Tests MatView()/MatLoad() with binary viewers for SBAIJ matrices.\n\n"; 2 3 #include <petscmat.h> 4 #include <petscviewer.h> 5 6 #include <petsc/private/hashtable.h> 7 static PetscReal MakeValue(PetscInt i,PetscInt j,PetscInt M) 8 { 9 PetscHash_t h = PetscHashCombine(PetscHashInt(i),PetscHashInt(j)); 10 return (PetscReal) ((h % 5 == 0) ? (1 + i + j*M) : 0); 11 } 12 13 static PetscErrorCode CheckValuesAIJ(Mat A) 14 { 15 PetscInt M,N,rstart,rend,i,j; 16 PetscReal v,w; 17 PetscScalar val; 18 PetscBool seqsbaij,mpisbaij,sbaij; 19 PetscErrorCode ierr; 20 21 PetscFunctionBegin; 22 ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr); 23 ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr); 24 ierr = PetscObjectTypeCompare((PetscObject)A,MATSEQSBAIJ,&seqsbaij);CHKERRQ(ierr); 25 ierr = PetscObjectTypeCompare((PetscObject)A,MATMPISBAIJ,&mpisbaij);CHKERRQ(ierr); 26 sbaij = (seqsbaij||mpisbaij) ? PETSC_TRUE : PETSC_FALSE; 27 for (i=rstart; i<rend; i++) { 28 for (j=(sbaij?i:0); j<N; j++) { 29 ierr = MatGetValue(A,i,j,&val);CHKERRQ(ierr); 30 v = MakeValue(i,j,M); w = PetscRealPart(val); 31 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); 32 } 33 } 34 PetscFunctionReturn(0); 35 } 36 37 int main(int argc,char **args) 38 { 39 Mat A; 40 PetscInt M = 24,N = 24,bs = 3; 41 PetscInt rstart,rend,i,j; 42 PetscErrorCode ierr; 43 PetscViewer view; 44 45 ierr = PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr; 46 /* 47 Create a parallel SBAIJ matrix shared by all processors 48 */ 49 ierr = MatCreateSBAIJ(PETSC_COMM_WORLD, 50 bs, 51 PETSC_DECIDE,PETSC_DECIDE, 52 M,N, 53 PETSC_DECIDE,NULL, 54 PETSC_DECIDE,NULL, 55 &A);CHKERRQ(ierr); 56 57 /* 58 Set values into the matrix 59 */ 60 ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr); 61 ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr); 62 for (i=rstart; i<rend; i++) { 63 for (j=0; j<N; j++) { 64 PetscReal v = MakeValue(i,j,M); 65 if (PetscAbsReal(v) > 0) { 66 ierr = MatSetValue(A,i,j,v,INSERT_VALUES);CHKERRQ(ierr); 67 } 68 } 69 } 70 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 71 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 72 ierr = MatViewFromOptions(A,NULL,"-mat_base_view");CHKERRQ(ierr); 73 74 /* 75 Store the binary matrix to a file 76 */ 77 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);CHKERRQ(ierr); 78 for (i=0; i<3; i++) { 79 ierr = MatView(A,view);CHKERRQ(ierr); 80 } 81 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 82 ierr = MatDestroy(&A);CHKERRQ(ierr); 83 84 /* 85 Now reload the matrix and check its values 86 */ 87 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 88 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 89 ierr = MatSetType(A,MATSBAIJ);CHKERRQ(ierr); 90 for (i=0; i<3; i++) { 91 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 92 ierr = MatLoad(A,view);CHKERRQ(ierr); 93 ierr = CheckValuesAIJ(A);CHKERRQ(ierr); 94 } 95 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 96 ierr = MatViewFromOptions(A,NULL,"-mat_load_view");CHKERRQ(ierr); 97 ierr = MatDestroy(&A);CHKERRQ(ierr); 98 99 /* 100 Reload in SEQSBAIJ matrix and check its values 101 */ 102 ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 103 ierr = MatCreate(PETSC_COMM_SELF,&A);CHKERRQ(ierr); 104 ierr = MatSetType(A,MATSEQSBAIJ);CHKERRQ(ierr); 105 for (i=0; i<3; i++) { 106 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 107 ierr = MatLoad(A,view);CHKERRQ(ierr); 108 ierr = CheckValuesAIJ(A);CHKERRQ(ierr); 109 } 110 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 111 ierr = MatDestroy(&A);CHKERRQ(ierr); 112 113 /* 114 Reload in MPISBAIJ matrix and check its values 115 */ 116 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 117 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 118 ierr = MatSetType(A,MATMPISBAIJ);CHKERRQ(ierr); 119 for (i=0; i<3; i++) { 120 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 121 ierr = MatLoad(A,view);CHKERRQ(ierr); 122 ierr = CheckValuesAIJ(A);CHKERRQ(ierr); 123 } 124 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 125 ierr = MatDestroy(&A);CHKERRQ(ierr); 126 127 ierr = PetscFinalize(); 128 return ierr; 129 } 130 131 /*TEST 132 133 testset: 134 args: -viewer_binary_mpiio 0 135 output_file: output/ex50.out 136 test: 137 suffix: stdio_1 138 nsize: 1 139 test: 140 suffix: stdio_2 141 nsize: 2 142 test: 143 suffix: stdio_3 144 nsize: 3 145 test: 146 suffix: stdio_4 147 nsize: 4 148 test: 149 suffix: stdio_5 150 nsize: 4 151 152 testset: 153 requires: mpiio 154 args: -viewer_binary_mpiio 1 155 output_file: output/ex50.out 156 test: 157 suffix: mpiio_1 158 nsize: 1 159 test: 160 suffix: mpiio_2 161 nsize: 2 162 test: 163 suffix: mpiio_3 164 nsize: 3 165 test: 166 suffix: mpiio_4 167 nsize: 4 168 test: 169 suffix: mpiio_5 170 nsize: 5 171 172 TEST*/ 173