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