Actual source code: ex46.c
petsc-3.12.4 2020-02-04
Report Typos and Errors

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

      4: /*T
      5:  Concepts: viewers^skipheader
      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;

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

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

     31:   VecView(x,viewer);

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

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

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

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

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

     59:   VecLoad(x,viewer);

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

     66:   PetscViewerDestroy(&viewer);
     67:   return(0);
     68: }

     70: PetscErrorCode VecFill(Vec x)
     71: {
     73:   PetscInt       i,s,e;

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

     85: PetscErrorCode VecCompare(Vec a,Vec b)
     86: {
     87:   PetscInt       locmin[2],locmax[2];
     88:   PetscReal      min[2],max[2];
     89:   Vec            ref;

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

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

     99:   PetscPrintf(PETSC_COMM_WORLD,"VecCompare\n");
    100:   PetscPrintf(PETSC_COMM_WORLD,"  min(a)   = %+1.2e [loc %D]\n",(double)min[0],locmin[0]);
    101:   PetscPrintf(PETSC_COMM_WORLD,"  max(a)   = %+1.2e [loc %D]\n",(double)max[0],locmax[0]);

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

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

    120: PetscErrorCode HeaderlessBinaryRead(const char name[])
    121: {
    123:   int            fdes;
    124:   PetscScalar    buffer[VEC_LEN];
    125:   PetscInt       i;
    126:   PetscMPIInt    rank;
    127:   PetscBool      dataverified = PETSC_TRUE;

    130:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
    131:   if (!rank) {
    132:     PetscBinaryOpen(name,FILE_MODE_READ,&fdes);
    133:     PetscBinaryRead(fdes,buffer,VEC_LEN,NULL,PETSC_SCALAR);
    134:     PetscBinaryClose(fdes);

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

    158: PetscErrorCode TestBinary(void)
    159: {
    161:   Vec            x,y;
    162:   PetscBool      skipheader = PETSC_TRUE;
    163:   PetscBool      usempiio = PETSC_FALSE;
    164:
    166:   VecCreate(PETSC_COMM_WORLD,&x);
    167:   VecSetSizes(x,PETSC_DECIDE,VEC_LEN);
    168:   VecSetFromOptions(x);
    169:   VecFill(x);
    170:   MyVecDump("xH.pbvec",skipheader,usempiio,x);

    172:   VecCreate(PETSC_COMM_WORLD,&y);
    173:   VecSetSizes(y,PETSC_DECIDE,VEC_LEN);
    174:   VecSetFromOptions(y);

    176:   MyVecLoad("xH.pbvec",skipheader,usempiio,y);
    177:   VecCompare(x,y);

    179:   VecDestroy(&y);
    180:   VecDestroy(&x);

    182:   HeaderlessBinaryRead("xH.pbvec");
    183:   return(0);
    184: }

    186: #if defined(PETSC_HAVE_MPIIO)
    187: PetscErrorCode TestBinaryMPIIO(void)
    188: {
    190:   Vec            x,y;
    191:   PetscBool      skipheader = PETSC_TRUE;
    192:   PetscBool      usempiio = PETSC_TRUE;

    195:   VecCreate(PETSC_COMM_WORLD,&x);
    196:   VecSetSizes(x,PETSC_DECIDE,VEC_LEN);
    197:   VecSetFromOptions(x);
    198:   VecFill(x);
    199:   MyVecDump("xHmpi.pbvec",skipheader,usempiio,x);

    201:   VecCreate(PETSC_COMM_WORLD,&y);
    202:   VecSetSizes(y,PETSC_DECIDE,VEC_LEN);
    203:   VecSetFromOptions(y);

    205:   MyVecLoad("xHmpi.pbvec",skipheader,usempiio,y);
    206:   VecCompare(x,y);

    208:   VecDestroy(&y);
    209:   VecDestroy(&x);

    211:   HeaderlessBinaryRead("xHmpi.pbvec");
    212:   return(0);
    213: }
    214: #endif

    216: int main(int argc,char **args)
    217: {
    219:   PetscBool      usempiio = PETSC_FALSE;

    221:   PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
    222:   PetscOptionsGetBool(NULL,NULL,"-usempiio",&usempiio,NULL);
    223:   if (!usempiio) {
    224:     TestBinary();
    225:   } else {
    226: #if defined(PETSC_HAVE_MPIIO)
    227:     TestBinaryMPIIO();
    228: #else
    229:     PetscPrintf(PETSC_COMM_WORLD,"Warning: Executing TestBinaryMPIIO() requires a working MPI-2 implementation\n");
    230: #endif
    231:   }
    232:   PetscFinalize();
    233:   return ierr;
    234: }

    237: /*TEST

    239:    test:
    240:       output_file: output/ex46_1_p1.out

    242:    test:
    243:       suffix: 2
    244:       nsize: 6
    245:       output_file: output/ex46_1_p6.out

    247:    test:
    248:       suffix: 3
    249:       nsize: 12
    250:       output_file: output/ex46_1_p12.out

    252:    test:
    253:       suffix: mpiio
    254:       nsize: 6
    255:       args: -usempiio
    256:       output_file: output/ex46_2_p6.out
    257:       requires: mpiio

    259: TEST*/