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