Actual source code: ex16.c

  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;

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

 27: #define CheckValuesIJ(A)  CheckValues(A,PETSC_FALSE)
 28: #define CheckValuesOne(A) CheckValues(A,PETSC_TRUE)

 30: int main(int argc,char **args)
 31: {
 32:   Mat            A;
 33:   PetscInt       i,j,M = 4,N = 3,rstart,rend;
 34:   PetscScalar    *array;
 35:   char           mattype[256];
 36:   PetscViewer    view;

 38:   PetscInitialize(&argc,&args,NULL,help);
 39:   PetscStrcpy(mattype,MATMPIDENSE);
 40:   PetscOptionsGetString(NULL,NULL,"-mat_type",mattype,sizeof(mattype),NULL);
 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);
 45:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 46:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 47:   MatConvert(A,mattype,MAT_INPLACE_MATRIX,&A);
 48:   /*
 49:      Set values into the matrix
 50:   */
 51:   for (i=0; i<M; i++) {
 52:     for (j=0; j<N; j++) {
 53:       PetscScalar v = (PetscReal)(1 + i + j*M);
 54:       MatSetValues(A,1,&i,1,&j,&v,INSERT_VALUES);
 55:     }
 56:   }
 57:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 58:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 59:   MatScale(A,2.0);
 60:   MatScale(A,1.0/2.0);

 62:   /*
 63:       Store the binary matrix to a file
 64:   */
 65:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);
 66:   for (i=0; i<2; i++) {
 67:     MatView(A,view);
 68:     PetscViewerPushFormat(view,PETSC_VIEWER_NATIVE);
 69:     MatView(A,view);
 70:     PetscViewerPopFormat(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,mattype);
 81:   for (i=0; i<4; i++) {
 82:     if (i > 0) MatZeroEntries(A);
 83:     MatLoad(A,view);
 84:     CheckValuesIJ(A);
 85:   }
 86:   PetscViewerDestroy(&view);

 88:   MatGetOwnershipRange(A,&rstart,&rend);
 89:   PetscMalloc1((rend-rstart)*N,&array);
 90:   for (i=0; i<(rend-rstart)*N; i++) array[i] = (PetscReal)1;
 91:   MatDensePlaceArray(A,array);
 92:   MatScale(A,2.0);
 93:   MatScale(A,1.0/2.0);
 94:   CheckValuesOne(A);
 95:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_WRITE,&view);
 96:   MatView(A,view);
 97:   MatDenseResetArray(A);
 98:   PetscFree(array);
 99:   CheckValuesIJ(A);
100:   PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);
101:   MatView(A,view);
102:   PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);
103:   PetscViewerDestroy(&view);
104:   MatDestroy(&A);

106:   MatCreate(PETSC_COMM_WORLD,&A);
107:   MatSetType(A,mattype);
108:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);
109:   MatLoad(A,view);
110:   CheckValuesOne(A);
111:   MatZeroEntries(A);
112:   PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);
113:   MatLoad(A,view);
114:   PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);
115:   CheckValuesIJ(A);
116:   PetscViewerDestroy(&view);
117:   MatDestroy(&A);

119:   {
120:     PetscInt m = PETSC_DECIDE, n = PETSC_DECIDE;
121:     PetscSplitOwnership(PETSC_COMM_WORLD,&m,&M);
122:     PetscSplitOwnership(PETSC_COMM_WORLD,&n,&N);
123:     /* TODO: MatCreateDense requires data!=NULL at all processes! */
124:     PetscMalloc1(m*N+1,&array);

126:     PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);
127:     MatCreateDense(PETSC_COMM_WORLD,m,n,M,N,array,&A);
128:     MatLoad(A,view);
129:     CheckValuesOne(A);
130:     PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);
131:     MatLoad(A,view);
132:     PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);
133:     CheckValuesIJ(A);
134:     MatDestroy(&A);
135:     PetscViewerDestroy(&view);

137:     MatCreateDense(PETSC_COMM_WORLD,m,n,M,N,array,&A);
138:     CheckValuesIJ(A);
139:     MatDestroy(&A);

141:     PetscFree(array);
142:   }

144:   PetscFinalize();
145:   return 0;
146: }

148: /*TEST

150:    testset:
151:       args: -viewer_binary_mpiio 0
152:       output_file: output/ex16.out
153:       test:
154:         suffix: stdio_1
155:         nsize: 1
156:         args: -mat_type seqdense
157:       test:
158:         suffix: stdio_2
159:         nsize: 2
160:       test:
161:         suffix: stdio_3
162:         nsize: 3
163:       test:
164:         suffix: stdio_4
165:         nsize: 4
166:       test:
167:         suffix: stdio_5
168:         nsize: 5
169:       test:
170:         requires: cuda
171:         args: -mat_type seqdensecuda
172:         suffix: stdio_cuda_1
173:         nsize: 1
174:       test:
175:         requires: cuda
176:         args: -mat_type mpidensecuda
177:         suffix: stdio_cuda_2
178:         nsize: 2
179:       test:
180:         requires: cuda
181:         args: -mat_type mpidensecuda
182:         suffix: stdio_cuda_3
183:         nsize: 3
184:       test:
185:         requires: cuda
186:         args: -mat_type mpidensecuda
187:         suffix: stdio_cuda_4
188:         nsize: 4
189:       test:
190:         requires: cuda
191:         args: -mat_type mpidensecuda
192:         suffix: stdio_cuda_5
193:         nsize: 5

195:    testset:
196:       requires: mpiio
197:       args: -viewer_binary_mpiio 1
198:       output_file: output/ex16.out
199:       test:
200:         suffix: mpiio_1
201:         nsize: 1
202:       test:
203:         suffix: mpiio_2
204:         nsize: 2
205:       test:
206:         suffix: mpiio_3
207:         nsize: 3
208:       test:
209:         suffix: mpiio_4
210:         nsize: 4
211:       test:
212:         suffix: mpiio_5
213:         nsize: 5
214:       test:
215:         requires: cuda
216:         args: -mat_type mpidensecuda
217:         suffix: mpiio_cuda_1
218:         nsize: 1
219:       test:
220:         requires: cuda
221:         args: -mat_type mpidensecuda
222:         suffix: mpiio_cuda_2
223:         nsize: 2
224:       test:
225:         requires: cuda
226:         args: -mat_type mpidensecuda
227:         suffix: mpiio_cuda_3
228:         nsize: 3
229:       test:
230:         requires: cuda
231:         args: -mat_type mpidensecuda
232:         suffix: mpiio_cuda_4
233:         nsize: 4
234:       test:
235:         requires: cuda
236:         args: -mat_type mpidensecuda
237:         suffix: mpiio_cuda_5
238:         nsize: 5

240: TEST*/