Actual source code: ex50.c

  1: static char help[] = "Tests MatView()/MatLoad() with binary viewers for SBAIJ 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:   PetscBool       seqsbaij,mpisbaij,sbaij;

 20:   MatGetSize(A,&M,&N);
 21:   MatGetOwnershipRange(A,&rstart,&rend);
 22:   PetscObjectTypeCompare((PetscObject)A,MATSEQSBAIJ,&seqsbaij);
 23:   PetscObjectTypeCompare((PetscObject)A,MATMPISBAIJ,&mpisbaij);
 24:   sbaij = (seqsbaij||mpisbaij) ? PETSC_TRUE : PETSC_FALSE;
 25:   for (i=rstart; i<rend; i++) {
 26:     for (j=(sbaij?i:0); j<N; j++) {
 27:       MatGetValue(A,i,j,&val);
 28:       v = MakeValue(i,j,M); w = PetscRealPart(val);
 30:     }
 31:   }
 32:   return 0;
 33: }

 35: int main(int argc,char **args)
 36: {
 37:   Mat            A;
 38:   PetscInt       M = 24,N = 24,bs = 3;
 39:   PetscInt       rstart,rend,i,j;
 40:   PetscViewer    view;

 42:   PetscInitialize(&argc,&args,NULL,help);
 43:   /*
 44:       Create a parallel SBAIJ matrix shared by all processors
 45:   */
 46:   MatCreateSBAIJ(PETSC_COMM_WORLD,bs,PETSC_DECIDE,PETSC_DECIDE,M,N,PETSC_DECIDE,NULL,PETSC_DECIDE,NULL,&A);

 48:   /*
 49:       Set values into the matrix
 50:   */
 51:   MatGetSize(A,&M,&N);
 52:   MatGetOwnershipRange(A,&rstart,&rend);
 53:   for (i=rstart; i<rend; i++) {
 54:     for (j=0; j<N; j++) {
 55:       PetscReal v = MakeValue(i,j,M);
 56:       if (PetscAbsReal(v) > 0) {
 57:         MatSetValue(A,i,j,v,INSERT_VALUES);
 58:       }
 59:     }
 60:   }
 61:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 62:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 63:   MatViewFromOptions(A,NULL,"-mat_base_view");

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

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

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

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

118:   PetscFinalize();
119:   return 0;
120: }

122: /*TEST

124:    testset:
125:       args: -viewer_binary_mpiio 0
126:       output_file: output/ex50.out
127:       test:
128:         suffix: stdio_1
129:         nsize: 1
130:       test:
131:         suffix: stdio_2
132:         nsize: 2
133:       test:
134:         suffix: stdio_3
135:         nsize: 3
136:       test:
137:         suffix: stdio_4
138:         nsize: 4
139:       test:
140:         suffix: stdio_5
141:         nsize: 4

143:    testset:
144:       requires: mpiio
145:       args: -viewer_binary_mpiio 1
146:       output_file: output/ex50.out
147:       test:
148:         suffix: mpiio_1
149:         nsize: 1
150:       test:
151:         suffix: mpiio_2
152:         nsize: 2
153:       test:
154:         suffix: mpiio_3
155:         nsize: 3
156:       test:
157:         suffix: mpiio_4
158:         nsize: 4
159:       test:
160:         suffix: mpiio_5
161:         nsize: 5

163: TEST*/