Actual source code: ex46.c


  2: static char help[] = "Tests PetscViewerBinary VecView()/VecLoad() function correctly when binary header is skipped.\n\n";

  4: /*T
  5:  Concepts: viewers^skipheader^mpiio
  6: T*/

  8: #include <petscviewer.h>
  9: #include <petscvec.h>

 11: #define VEC_LEN 10
 12: const PetscReal test_values[] = { 0.311256, 88.068, 11.077444, 9953.62, 7.345, 64.8943, 3.1458, 6699.95, 0.00084, 0.0647 };

 14: PetscErrorCode MyVecDump(const char fname[],PetscBool skippheader,PetscBool usempiio,Vec x)
 15: {
 16:   MPI_Comm       comm;
 17:   PetscViewer    viewer;
 18:   PetscBool      ismpiio,isskip;

 21:   PetscObjectGetComm((PetscObject)x,&comm);

 23:   PetscViewerCreate(comm,&viewer);
 24:   PetscViewerSetType(viewer,PETSCVIEWERBINARY);
 25:   if (skippheader) PetscViewerBinarySetSkipHeader(viewer,PETSC_TRUE);
 26:   PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
 27:   if (usempiio) PetscViewerBinarySetUseMPIIO(viewer,PETSC_TRUE);
 28:   PetscViewerFileSetName(viewer,fname);

 30:   VecView(x,viewer);

 32:   PetscViewerBinaryGetUseMPIIO(viewer,&ismpiio);
 33:   if (ismpiio) PetscPrintf(comm,"*** PetscViewer[write] using MPI-IO ***\n");
 34:   PetscViewerBinaryGetSkipHeader(viewer,&isskip);
 35:   if (isskip) PetscPrintf(comm,"*** PetscViewer[write] skipping header ***\n");

 37:   PetscViewerDestroy(&viewer);
 38:   return 0;
 39: }

 41: PetscErrorCode MyVecLoad(const char fname[],PetscBool skippheader,PetscBool usempiio,Vec x)
 42: {
 43:   MPI_Comm       comm;
 44:   PetscViewer    viewer;
 45:   PetscBool      ismpiio,isskip;

 48:   PetscObjectGetComm((PetscObject)x,&comm);

 50:   PetscViewerCreate(comm,&viewer);
 51:   PetscViewerSetType(viewer,PETSCVIEWERBINARY);
 52:   if (skippheader) PetscViewerBinarySetSkipHeader(viewer,PETSC_TRUE);
 53:   PetscViewerFileSetMode(viewer,FILE_MODE_READ);
 54:   if (usempiio) PetscViewerBinarySetUseMPIIO(viewer,PETSC_TRUE);
 55:   PetscViewerFileSetName(viewer,fname);

 57:   VecLoad(x,viewer);

 59:   PetscViewerBinaryGetSkipHeader(viewer,&isskip);
 60:   if (isskip) PetscPrintf(comm,"*** PetscViewer[load] skipping header ***\n");
 61:   PetscViewerBinaryGetUseMPIIO(viewer,&ismpiio);
 62:   if (ismpiio) PetscPrintf(comm,"*** PetscViewer[load] using MPI-IO ***\n");

 64:   PetscViewerDestroy(&viewer);
 65:   return 0;
 66: }

 68: PetscErrorCode VecFill(Vec x)
 69: {
 70:   PetscInt       i,s,e;

 73:   VecGetOwnershipRange(x,&s,&e);
 74:   for (i=s; i<e; i++) {
 75:     VecSetValue(x,i,(PetscScalar)test_values[i],INSERT_VALUES);
 76:   }
 77:   VecAssemblyBegin(x);
 78:   VecAssemblyEnd(x);
 79:   return 0;
 80: }

 82: PetscErrorCode VecCompare(Vec a,Vec b)
 83: {
 84:   PetscInt       locmin[2],locmax[2];
 85:   PetscReal      min[2],max[2];
 86:   Vec            ref;

 89:   VecMin(a,&locmin[0],&min[0]);
 90:   VecMax(a,&locmax[0],&max[0]);

 92:   VecMin(b,&locmin[1],&min[1]);
 93:   VecMax(b,&locmax[1],&max[1]);

 95:   PetscPrintf(PETSC_COMM_WORLD,"VecCompare\n");
 96:   PetscPrintf(PETSC_COMM_WORLD,"  min(a)   = %+1.2e [loc %" PetscInt_FMT "]\n",(double)min[0],locmin[0]);
 97:   PetscPrintf(PETSC_COMM_WORLD,"  max(a)   = %+1.2e [loc %" PetscInt_FMT "]\n",(double)max[0],locmax[0]);

 99:   PetscPrintf(PETSC_COMM_WORLD,"  min(b)   = %+1.2e [loc %" PetscInt_FMT "]\n",(double)min[1],locmin[1]);
100:   PetscPrintf(PETSC_COMM_WORLD,"  max(b)   = %+1.2e [loc %" PetscInt_FMT "]\n",(double)max[1],locmax[1]);

102:   VecDuplicate(a,&ref);
103:   VecCopy(a,ref);
104:   VecAXPY(ref,-1.0,b);
105:   VecMin(ref,&locmin[0],&min[0]);
106:   if (PetscAbsReal(min[0]) > 1.0e-10) {
107:     PetscPrintf(PETSC_COMM_WORLD,"  ERROR: min(a-b) > 1.0e-10\n");
108:     PetscPrintf(PETSC_COMM_WORLD,"  min(a-b) = %+1.10e\n",(double)PetscAbsReal(min[0]));
109:   } else {
110:     PetscPrintf(PETSC_COMM_WORLD,"  min(a-b) < 1.0e-10\n");
111:   }
112:   VecDestroy(&ref);
113:   return 0;
114: }

116: PetscErrorCode HeaderlessBinaryRead(const char name[])
117: {
118:   int            fdes;
119:   PetscScalar    buffer[VEC_LEN];
120:   PetscInt       i;
121:   PetscMPIInt    rank;
122:   PetscBool      dataverified = PETSC_TRUE;

125:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
126:   if (rank == 0) {
127:     PetscBinaryOpen(name,FILE_MODE_READ,&fdes);
128:     PetscBinaryRead(fdes,buffer,VEC_LEN,NULL,PETSC_SCALAR);
129:     PetscBinaryClose(fdes);

131:     for (i=0; i<VEC_LEN; i++) {
132:       PetscScalar v;
133:       v = PetscAbsScalar(test_values[i]-buffer[i]);
134: #if defined(PETSC_USE_COMPLEX)
135:       if ((PetscRealPart(v) > 1.0e-10) || (PetscImaginaryPart(v) > 1.0e-10)) {
136:         PetscPrintf(PETSC_COMM_SELF,"ERROR: Difference > 1.0e-10 occurred (delta = (%+1.12e,%+1.12e) [loc %" PetscInt_FMT "])\n",(double)PetscRealPart(buffer[i]),(double)PetscImaginaryPart(buffer[i]),i);
137:         dataverified = PETSC_FALSE;
138:       }
139: #else
140:       if (PetscRealPart(v) > 1.0e-10) {
141:         PetscPrintf(PETSC_COMM_SELF,"ERROR: Difference > 1.0e-10 occurred (delta = %+1.12e [loc %" PetscInt_FMT "])\n",(double)PetscRealPart(buffer[i]),i);
142:         dataverified = PETSC_FALSE;
143:       }
144: #endif
145:     }
146:     if (dataverified) {
147:       PetscPrintf(PETSC_COMM_SELF,"Headerless read of data verified\n");
148:     }
149:   }
150:   return 0;
151: }

153: PetscErrorCode TestBinary(void)
154: {
155:   Vec            x,y;
156:   PetscBool      skipheader = PETSC_TRUE;
157:   PetscBool      usempiio = PETSC_FALSE;

160:   VecCreate(PETSC_COMM_WORLD,&x);
161:   VecSetSizes(x,PETSC_DECIDE,VEC_LEN);
162:   VecSetFromOptions(x);
163:   VecFill(x);
164:   MyVecDump("xH.pbvec",skipheader,usempiio,x);

166:   VecCreate(PETSC_COMM_WORLD,&y);
167:   VecSetSizes(y,PETSC_DECIDE,VEC_LEN);
168:   VecSetFromOptions(y);

170:   MyVecLoad("xH.pbvec",skipheader,usempiio,y);
171:   VecCompare(x,y);

173:   VecDestroy(&y);
174:   VecDestroy(&x);

176:   HeaderlessBinaryRead("xH.pbvec");
177:   return 0;
178: }

180: #if defined(PETSC_HAVE_MPIIO)
181: PetscErrorCode TestBinaryMPIIO(void)
182: {
183:   Vec            x,y;
184:   PetscBool      skipheader = PETSC_TRUE;
185:   PetscBool      usempiio = PETSC_TRUE;

188:   VecCreate(PETSC_COMM_WORLD,&x);
189:   VecSetSizes(x,PETSC_DECIDE,VEC_LEN);
190:   VecSetFromOptions(x);
191:   VecFill(x);
192:   MyVecDump("xHmpi.pbvec",skipheader,usempiio,x);

194:   VecCreate(PETSC_COMM_WORLD,&y);
195:   VecSetSizes(y,PETSC_DECIDE,VEC_LEN);
196:   VecSetFromOptions(y);

198:   MyVecLoad("xHmpi.pbvec",skipheader,usempiio,y);
199:   VecCompare(x,y);

201:   VecDestroy(&y);
202:   VecDestroy(&x);

204:   HeaderlessBinaryRead("xHmpi.pbvec");
205:   return 0;
206: }
207: #endif

209: int main(int argc,char **args)
210: {
211:   PetscBool      usempiio = PETSC_FALSE;

213:   PetscInitialize(&argc,&args,(char*)0,help);
214:   PetscOptionsGetBool(NULL,NULL,"-usempiio",&usempiio,NULL);
215:   if (!usempiio) {
216:     TestBinary();
217:   } else {
218: #if defined(PETSC_HAVE_MPIIO)
219:     TestBinaryMPIIO();
220: #else
221:     PetscPrintf(PETSC_COMM_WORLD,"Warning: Executing TestBinaryMPIIO() requires a working MPI-2 implementation\n");
222: #endif
223:   }
224:   PetscFinalize();
225:   return 0;
226: }

228: /*TEST

230:    test:
231:       output_file: output/ex46_1_p1.out

233:    test:
234:       suffix: 2
235:       nsize: 6
236:       output_file: output/ex46_1_p6.out

238:    test:
239:       suffix: 3
240:       nsize: 12
241:       output_file: output/ex46_1_p12.out

243:    testset:
244:       requires: mpiio
245:       args: -usempiio
246:       test:
247:          suffix: mpiio_1
248:          output_file: output/ex46_2_p1.out
249:       test:
250:          suffix: mpiio_2
251:          nsize: 6
252:          output_file: output/ex46_2_p6.out
253:       test:
254:          suffix: mpiio_3
255:          nsize: 12
256:          output_file: output/ex46_2_p12.out

258: TEST*/