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

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

  6: static PetscErrorCode CheckValues(Mat A,PetscBool one)
  7: {
  8:   const PetscScalar *array;
  9:   PetscInt          M,N,rstart,rend,lda,i,j;
 10:   PetscErrorCode    ierr;

 13:   MatDenseGetArrayRead(A,&array);
 14:   MatDenseGetLDA(A,&lda);
 15:   MatGetSize(A,&M,&N);
 16:   MatGetOwnershipRange(A,&rstart,&rend);
 17:   for (i=rstart; i<rend; i++) {
 18:     for (j=0; j<N; j++) {
 19:       PetscInt ii = i - rstart, jj = j;
 20:       PetscReal v = (PetscReal)(one ? 1 : (1 + i + j*M));
 21:       PetscReal w = PetscRealPart(array[ii + jj*lda]);
 22:       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);
 23:     }
 24:   }
 25:   MatDenseRestoreArrayRead(A,&array);
 26:   return(0);
 27: }

 29: #define CheckValuesIJ(A)  CheckValues(A,PETSC_FALSE)
 30: #define CheckValuesOne(A) CheckValues(A,PETSC_TRUE)

 32: int main(int argc,char **args)
 33: {
 34:   Mat            A;
 35:   PetscInt       i,j,M = 4,N = 3,rstart,rend;
 37:   PetscScalar    *array;
 38:   PetscViewer    view;

 40:   PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr;
 41:   /*
 42:       Create a parallel dense matrix shared by all processors
 43:   */
 44:   MatCreateDense(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,M,N,NULL,&A);

 46:   /*
 47:      Set values into the matrix
 48:   */
 49:   for (i=0; i<M; i++) {
 50:     for (j=0; j<N; j++) {
 51:       PetscScalar v = (PetscReal)(1 + i + j*M);
 52:       MatSetValues(A,1,&i,1,&j,&v,INSERT_VALUES);
 53:     }
 54:   }
 55:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 56:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);

 58:   /*
 59:       Store the binary matrix to a file
 60:   */
 61:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);
 62:   for (i=0; i<2; i++) {
 63:     MatView(A,view);
 64:     PetscViewerPushFormat(view,PETSC_VIEWER_NATIVE);
 65:     MatView(A,view);
 66:     PetscViewerPopFormat(view);
 67:   }
 68:   PetscViewerDestroy(&view);
 69:   MatDestroy(&A);

 71:   /*
 72:       Now reload the matrix and check its values
 73:   */
 74:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);
 75:   MatCreate(PETSC_COMM_WORLD,&A);
 76:   MatSetType(A,MATMPIDENSE);
 77:   for (i=0; i<4; i++) {
 78:     if (i > 0) {MatZeroEntries(A);}
 79:     MatLoad(A,view);
 80:     CheckValuesIJ(A);
 81:   }
 82:   PetscViewerDestroy(&view);

 84:   MatGetOwnershipRange(A,&rstart,&rend);
 85:   PetscMalloc1((rend-rstart)*N,&array);
 86:   for (i=0; i<(rend-rstart)*N; i++) array[i] = (PetscReal)1;
 87:   MatDensePlaceArray(A,array);
 88:   CheckValuesOne(A);
 89:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_WRITE,&view);
 90:   MatView(A,view);
 91:   MatDenseResetArray(A);
 92:   PetscFree(array);
 93:   CheckValuesIJ(A);
 94:   PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);
 95:   MatView(A,view);
 96:   PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);
 97:   PetscViewerDestroy(&view);
 98:   MatDestroy(&A);

100:   MatCreate(PETSC_COMM_WORLD,&A);
101:   MatSetType(A,MATMPIDENSE);
102:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);
103:   MatLoad(A,view);
104:   CheckValuesOne(A);
105:   MatZeroEntries(A);
106:   PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);
107:   MatLoad(A,view);
108:   PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);
109:   CheckValuesIJ(A);
110:   PetscViewerDestroy(&view);
111:   MatDestroy(&A);

113:   {
114:     PetscInt m = PETSC_DECIDE, n = PETSC_DECIDE;
115:     PetscSplitOwnership(PETSC_COMM_WORLD,&m,&M);
116:     PetscSplitOwnership(PETSC_COMM_WORLD,&n,&N);
117:     /* TODO: MatCreateDense requires data!=NULL at all processes! */
118:     PetscMalloc1(m*N+1,&array);

120:     PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);
121:     MatCreateDense(PETSC_COMM_WORLD,m,n,M,N,array,&A);
122:     MatLoad(A,view);
123:     CheckValuesOne(A);
124:     PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);
125:     MatLoad(A,view);
126:     PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);
127:     CheckValuesIJ(A);
128:     MatDestroy(&A);
129:     PetscViewerDestroy(&view);

131:     MatCreateDense(PETSC_COMM_WORLD,m,n,M,N,array,&A);
132:     CheckValuesIJ(A);
133:     MatDestroy(&A);

135:     PetscFree(array);
136:   }

138:   PetscFinalize();
139:   return ierr;
140: }


143: /*TEST

145:    testset:
146:       args: -viewer_binary_mpiio 0
147:       output_file: output/ex16.out
148:       test:
149:         suffix: stdio_1
150:         nsize: 1
151:       test:
152:         suffix: stdio_2
153:         nsize: 2
154:       test:
155:         suffix: stdio_3
156:         nsize: 3
157:       test:
158:         suffix: stdio_4
159:         nsize: 4
160:       test:
161:         suffix: stdio_5
162:         nsize: 5

164:    testset:
165:       requires: mpiio
166:       args: -viewer_binary_mpiio 1
167:       output_file: output/ex16.out
168:       test:
169:         suffix: mpiio_1
170:         nsize: 1
171:       test:
172:         suffix: mpiio_2
173:         nsize: 2
174:       test:
175:         suffix: mpiio_3
176:         nsize: 3
177:       test:
178:         suffix: mpiio_4
179:         nsize: 4
180:       test:
181:         suffix: mpiio_5
182:         nsize: 5


185: TEST*/