1: static char help[] = "Tests MatView()/MatLoad() with binary viewers for AIJ matrices.\n\n";

  3: #include <petscmat.h>
  4: #include <petscviewer.h>

  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: }

 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;

 21:   MatGetSize(A,&M,&N);
 22:   MatGetOwnershipRange(A,&rstart,&rend);
 23:   for (i=rstart; i<rend; i++) {
 24:     for (j=0; j<N; j++) {
 25:       MatGetValue(A,i,j,&val);
 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:   return(0);
 31: }

 33: int main(int argc,char **args)
 34: {
 35:   Mat            A;
 36:   PetscInt       M = 11,N = 13;
 37:   PetscInt       rstart,rend,i,j;
 39:   PetscViewer    view;

 41:   PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr;
 42:   /*
 43:       Create a parallel AIJ matrix shared by all processors
 44:   */
 45:   MatCreateAIJ(PETSC_COMM_WORLD,
 46:                       PETSC_DECIDE,PETSC_DECIDE,
 47:                       M,N,
 48:                       PETSC_DECIDE,NULL,
 49:                       PETSC_DECIDE,NULL,
 50:                       &A);

 52:   /*
 53:       Set values into the matrix
 54:   */
 55:   MatGetOwnershipRange(A,&rstart,&rend);
 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:         MatSetValue(A,i,j,v,INSERT_VALUES);
 61:       }
 62:     }
 63:   }
 64:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 65:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 66:   MatViewFromOptions(A,NULL,"-mat_base_view");

 68:   /*
 69:       Store the binary matrix to a file
 70:   */
 71:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);
 72:   for (i=0; i<3; i++) {
 73:     MatView(A,view);
 74:   }
 75:   PetscViewerDestroy(&view);
 76:   MatDestroy(&A);

 78:   /*
 79:       Now reload the matrix and check its values
 80:   */
 81:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);
 82:   MatCreate(PETSC_COMM_WORLD,&A);
 83:   MatSetType(A,MATAIJ);
 84:   for (i=0; i<3; i++) {
 85:     if (i > 0) {MatZeroEntries(A);}
 86:     MatLoad(A,view);
 87:     CheckValuesAIJ(A);
 88:   }
 89:   PetscViewerDestroy(&view);
 90:   MatViewFromOptions(A,NULL,"-mat_load_view");
 91:   MatDestroy(&A);

 93:   /*
 94:       Reload in SEQAIJ matrix and check its values
 95:   */
 96:   PetscViewerBinaryOpen(PETSC_COMM_SELF,"matrix.dat",FILE_MODE_READ,&view);
 97:   MatCreate(PETSC_COMM_SELF,&A);
 98:   MatSetType(A,MATSEQAIJ);
 99:   for (i=0; i<3; i++) {
100:     if (i > 0) {MatZeroEntries(A);}
101:     MatLoad(A,view);
102:     CheckValuesAIJ(A);
103:   }
104:   PetscViewerDestroy(&view);
105:   MatDestroy(&A);

107:   /*
108:      Reload in MPIAIJ matrix and check its values
109:   */
110:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);
111:   MatCreate(PETSC_COMM_WORLD,&A);
112:   MatSetType(A,MATMPIAIJ);
113:   for (i=0; i<3; i++) {
114:     if (i > 0) {MatZeroEntries(A);}
115:     MatLoad(A,view);
116:     CheckValuesAIJ(A);
117:   }
118:   PetscViewerDestroy(&view);
119:   MatDestroy(&A);

121:   PetscFinalize();
122:   return ierr;
123: }


126: /*TEST

128:    testset:
129:       args: -viewer_binary_mpiio 0
130:       output_file: output/ex44.out
131:       test:
132:         suffix: stdio_1
133:         nsize: 1
134:       test:
135:         suffix: stdio_2
136:         nsize: 2
137:       test:
138:         suffix: stdio_3
139:         nsize: 3
140:       test:
141:         suffix: stdio_4
142:         nsize: 4
143:       test:
144:         suffix: stdio_15
145:         nsize: 15

147:    testset:
148:       requires: mpiio
149:       args: -viewer_binary_mpiio 1
150:       output_file: output/ex44.out
151:       test:
152:         suffix: mpiio_1
153:         nsize: 1
154:       test:
155:         suffix: mpiio_2
156:         nsize: 2
157:       test:
158:         suffix: mpiio_3
159:         nsize: 3
160:       test:
161:         suffix: mpiio_4
162:         nsize: 4
163:       test:
164:         suffix: mpiio_15
165:         nsize: 15

167: TEST*/