Actual source code: ex46.c


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

  4: #include <petscviewer.h>
  5: #include <petscvec.h>

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

 10: PetscErrorCode MyVecDump(const char fname[],PetscBool skippheader,PetscBool usempiio,Vec x)
 11: {
 12:   MPI_Comm       comm;
 13:   PetscViewer    viewer;
 14:   PetscBool      ismpiio,isskip;

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

 19:   PetscViewerCreate(comm,&viewer);
 20:   PetscViewerSetType(viewer,PETSCVIEWERBINARY);
 21:   if (skippheader) PetscViewerBinarySetSkipHeader(viewer,PETSC_TRUE);
 22:   PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
 23:   if (usempiio) PetscViewerBinarySetUseMPIIO(viewer,PETSC_TRUE);
 24:   PetscViewerFileSetName(viewer,fname);

 26:   VecView(x,viewer);

 28:   PetscViewerBinaryGetUseMPIIO(viewer,&ismpiio);
 29:   if (ismpiio) PetscPrintf(comm,"*** PetscViewer[write] using MPI-IO ***\n");
 30:   PetscViewerBinaryGetSkipHeader(viewer,&isskip);
 31:   if (isskip) PetscPrintf(comm,"*** PetscViewer[write] skipping header ***\n");

 33:   PetscViewerDestroy(&viewer);
 34:   return 0;
 35: }

 37: PetscErrorCode MyVecLoad(const char fname[],PetscBool skippheader,PetscBool usempiio,Vec x)
 38: {
 39:   MPI_Comm       comm;
 40:   PetscViewer    viewer;
 41:   PetscBool      ismpiio,isskip;

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

 46:   PetscViewerCreate(comm,&viewer);
 47:   PetscViewerSetType(viewer,PETSCVIEWERBINARY);
 48:   if (skippheader) PetscViewerBinarySetSkipHeader(viewer,PETSC_TRUE);
 49:   PetscViewerFileSetMode(viewer,FILE_MODE_READ);
 50:   if (usempiio) PetscViewerBinarySetUseMPIIO(viewer,PETSC_TRUE);
 51:   PetscViewerFileSetName(viewer,fname);

 53:   VecLoad(x,viewer);

 55:   PetscViewerBinaryGetSkipHeader(viewer,&isskip);
 56:   if (isskip) PetscPrintf(comm,"*** PetscViewer[load] skipping header ***\n");
 57:   PetscViewerBinaryGetUseMPIIO(viewer,&ismpiio);
 58:   if (ismpiio) PetscPrintf(comm,"*** PetscViewer[load] using MPI-IO ***\n");

 60:   PetscViewerDestroy(&viewer);
 61:   return 0;
 62: }

 64: PetscErrorCode VecFill(Vec x)
 65: {
 66:   PetscInt       i,s,e;

 69:   VecGetOwnershipRange(x,&s,&e);
 70:   for (i=s; i<e; i++) {
 71:     VecSetValue(x,i,(PetscScalar)test_values[i],INSERT_VALUES);
 72:   }
 73:   VecAssemblyBegin(x);
 74:   VecAssemblyEnd(x);
 75:   return 0;
 76: }

 78: PetscErrorCode VecCompare(Vec a,Vec b)
 79: {
 80:   PetscInt       locmin[2],locmax[2];
 81:   PetscReal      min[2],max[2];
 82:   Vec            ref;

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

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

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

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

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

112: PetscErrorCode HeaderlessBinaryRead(const char name[])
113: {
114:   int            fdes;
115:   PetscScalar    buffer[VEC_LEN];
116:   PetscInt       i;
117:   PetscMPIInt    rank;
118:   PetscBool      dataverified = PETSC_TRUE;

121:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
122:   if (rank == 0) {
123:     PetscBinaryOpen(name,FILE_MODE_READ,&fdes);
124:     PetscBinaryRead(fdes,buffer,VEC_LEN,NULL,PETSC_SCALAR);
125:     PetscBinaryClose(fdes);

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

149: PetscErrorCode TestBinary(void)
150: {
151:   Vec            x,y;
152:   PetscBool      skipheader = PETSC_TRUE;
153:   PetscBool      usempiio = PETSC_FALSE;

156:   VecCreate(PETSC_COMM_WORLD,&x);
157:   VecSetSizes(x,PETSC_DECIDE,VEC_LEN);
158:   VecSetFromOptions(x);
159:   VecFill(x);
160:   MyVecDump("xH.pbvec",skipheader,usempiio,x);

162:   VecCreate(PETSC_COMM_WORLD,&y);
163:   VecSetSizes(y,PETSC_DECIDE,VEC_LEN);
164:   VecSetFromOptions(y);

166:   MyVecLoad("xH.pbvec",skipheader,usempiio,y);
167:   VecCompare(x,y);

169:   VecDestroy(&y);
170:   VecDestroy(&x);

172:   HeaderlessBinaryRead("xH.pbvec");
173:   return 0;
174: }

176: #if defined(PETSC_HAVE_MPIIO)
177: PetscErrorCode TestBinaryMPIIO(void)
178: {
179:   Vec            x,y;
180:   PetscBool      skipheader = PETSC_TRUE;
181:   PetscBool      usempiio = PETSC_TRUE;

184:   VecCreate(PETSC_COMM_WORLD,&x);
185:   VecSetSizes(x,PETSC_DECIDE,VEC_LEN);
186:   VecSetFromOptions(x);
187:   VecFill(x);
188:   MyVecDump("xHmpi.pbvec",skipheader,usempiio,x);

190:   VecCreate(PETSC_COMM_WORLD,&y);
191:   VecSetSizes(y,PETSC_DECIDE,VEC_LEN);
192:   VecSetFromOptions(y);

194:   MyVecLoad("xHmpi.pbvec",skipheader,usempiio,y);
195:   VecCompare(x,y);

197:   VecDestroy(&y);
198:   VecDestroy(&x);

200:   HeaderlessBinaryRead("xHmpi.pbvec");
201:   return 0;
202: }
203: #endif

205: int main(int argc,char **args)
206: {
207:   PetscBool      usempiio = PETSC_FALSE;

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

224: /*TEST

226:    test:
227:       output_file: output/ex46_1_p1.out

229:    test:
230:       suffix: 2
231:       nsize: 6
232:       output_file: output/ex46_1_p6.out

234:    test:
235:       suffix: 3
236:       nsize: 12
237:       output_file: output/ex46_1_p12.out

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

254: TEST*/