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,bs,PETSC_DECIDE,PETSC_DECIDE,M,N,PETSC_DECIDE,NULL,PETSC_DECIDE,NULL,&A);CHKERRQ(ierr); 50 51 /* 52 Set values into the matrix 53 */ 54 ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr); 55 ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr); 56 for (i=rstart; i<rend; i++) { 57 for (j=0; j<N; j++) { 58 PetscReal v = MakeValue(i,j,M); 59 if (PetscAbsReal(v) > 0) { 60 ierr = MatSetValue(A,i,j,v,INSERT_VALUES);CHKERRQ(ierr); 61 } 62 } 63 } 64 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 65 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 66 ierr = MatViewFromOptions(A,NULL,"-mat_base_view");CHKERRQ(ierr); 67 68 /* 69 Store the binary matrix to a file 70 */ 71 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);CHKERRQ(ierr); 72 for (i=0; i<3; i++) { 73 ierr = MatView(A,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,MATSBAIJ);CHKERRQ(ierr); 84 for (i=0; i<3; i++) { 85 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 86 ierr = MatLoad(A,view);CHKERRQ(ierr); 87 ierr = CheckValuesAIJ(A);CHKERRQ(ierr); 88 } 89 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 90 ierr = MatViewFromOptions(A,NULL,"-mat_load_view");CHKERRQ(ierr); 91 ierr = MatDestroy(&A);CHKERRQ(ierr); 92 93 /* 94 Reload in SEQSBAIJ matrix and check its values 95 */ 96 ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 97 ierr = MatCreate(PETSC_COMM_SELF,&A);CHKERRQ(ierr); 98 ierr = MatSetType(A,MATSEQSBAIJ);CHKERRQ(ierr); 99 for (i=0; i<3; i++) { 100 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 101 ierr = MatLoad(A,view);CHKERRQ(ierr); 102 ierr = CheckValuesAIJ(A);CHKERRQ(ierr); 103 } 104 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 105 ierr = MatDestroy(&A);CHKERRQ(ierr); 106 107 /* 108 Reload in MPISBAIJ matrix and check its values 109 */ 110 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 111 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 112 ierr = MatSetType(A,MATMPISBAIJ);CHKERRQ(ierr); 113 for (i=0; i<3; i++) { 114 if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 115 ierr = MatLoad(A,view);CHKERRQ(ierr); 116 ierr = CheckValuesAIJ(A);CHKERRQ(ierr); 117 } 118 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 119 ierr = MatDestroy(&A);CHKERRQ(ierr); 120 121 ierr = PetscFinalize(); 122 return ierr; 123 } 124 125 /*TEST 126 127 testset: 128 args: -viewer_binary_mpiio 0 129 output_file: output/ex50.out 130 test: 131 suffix: stdio_1 132 nsize: 1 133 test: 134 suffix: stdio_2 135 nsize: 2 136 test: 137 suffix: stdio_3 138 nsize: 3 139 test: 140 suffix: stdio_4 141 nsize: 4 142 test: 143 suffix: stdio_5 144 nsize: 4 145 146 testset: 147 requires: mpiio 148 args: -viewer_binary_mpiio 1 149 output_file: output/ex50.out 150 test: 151 suffix: mpiio_1 152 nsize: 1 153 test: 154 suffix: mpiio_2 155 nsize: 2 156 test: 157 suffix: mpiio_3 158 nsize: 3 159 test: 160 suffix: mpiio_4 161 nsize: 4 162 test: 163 suffix: mpiio_5 164 nsize: 5 165 166 TEST*/ 167