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