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 PetscErrorCode ierr; 19 20 PetscFunctionBegin; 21 ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr); 22 ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr); 23 for (i=rstart; i<rend; i++) { 24 for (j=0; j<N; j++) { 25 ierr = MatGetValue(A,i,j,&val);CHKERRQ(ierr); 26 v = MakeValue(i,j,M); w = PetscRealPart(val); 27 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); 28 } 29 } 30 PetscFunctionReturn(0); 31 } 32 33 int main(int argc,char **args) 34 { 35 Mat A; 36 PetscInt M = 11,N = 13; 37 PetscInt rstart,rend,i,j; 38 PetscErrorCode ierr; 39 PetscViewer view; 40 41 ierr = PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr; 42 /* 43 Create a parallel AIJ matrix shared by all processors 44 */ 45 ierr = MatCreateAIJ(PETSC_COMM_WORLD, 46 PETSC_DECIDE,PETSC_DECIDE, 47 M,N, 48 PETSC_DECIDE,NULL, 49 PETSC_DECIDE,NULL, 50 &A);CHKERRQ(ierr); 51 52 /* 53 Set values into the matrix 54 */ 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,MATAIJ);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 SEQAIJ 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,MATSEQAIJ);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 MPIAIJ 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,MATMPIAIJ);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/ex44.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_15 144 nsize: 15 145 146 testset: 147 requires: mpiio 148 args: -viewer_binary_mpiio 1 149 output_file: output/ex44.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_15 164 nsize: 15 165 166 TEST*/ 167