Actual source code: ex47.c

  1: static char help[] = "Tests PetscViewerHDF5 VecView()/VecLoad() function.\n\n";

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

  7: int main(int argc, char **args)
  8: {
  9:   Vec         x, y;
 10:   PetscReal   norm, dnorm;
 11:   PetscViewer H5viewer;
 12:   char        filename[PETSC_MAX_PATH_LEN];
 13:   PetscBool   flg;

 15:   PetscFunctionBeginUser;
 16:   PetscCall(PetscInitialize(&argc, &args, NULL, help));
 17:   PetscCall(PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), &flg));
 18:   if (!flg) PetscCall(PetscStrncpy(filename, "x.h5", sizeof(filename)));
 19:   PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
 20:   PetscCall(VecSetFromOptions(x));
 21:   PetscCall(VecSetSizes(x, 11, PETSC_DETERMINE));
 22:   PetscCall(VecSet(x, 22.3));

 24:   PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, filename, FILE_MODE_WRITE, &H5viewer));
 25:   PetscCall(PetscViewerSetFromOptions(H5viewer));

 27:   /* Write the Vec without one extra dimension for BS */
 28:   PetscCall(PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_FALSE));
 29:   PetscCall(PetscObjectSetName((PetscObject)x, "noBsDim"));
 30:   PetscCall(VecView(x, H5viewer));

 32:   /* Write the Vec with one extra, 1-sized, dimension for BS */
 33:   PetscCall(PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_TRUE));
 34:   PetscCall(PetscObjectSetName((PetscObject)x, "bsDim"));
 35:   PetscCall(VecView(x, H5viewer));

 37:   PetscCall(PetscViewerDestroy(&H5viewer));
 38:   PetscCall(VecDuplicate(x, &y));

 40:   /* Create the HDF5 viewer for reading */
 41:   PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, filename, FILE_MODE_READ, &H5viewer));
 42:   PetscCall(PetscViewerSetFromOptions(H5viewer));

 44:   /* Load the Vec without the BS dim and compare */
 45:   PetscCall(PetscObjectSetName((PetscObject)y, "noBsDim"));
 46:   PetscCall(VecLoad(y, H5viewer));
 47:   PetscCall(VecAXPY(y, -1.0, x));
 48:   PetscCall(VecNorm(y, NORM_2, &norm));
 49:   PetscCall(VecNorm(x, NORM_2, &dnorm));
 50:   PetscCheck(norm / dnorm <= 1.e-6, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Vec read in 'noBsDim' does not match vector written out %g", (double)(norm / dnorm));

 52:   /* Load the Vec with one extra, 1-sized, BS dim and compare */
 53:   PetscCall(PetscObjectSetName((PetscObject)y, "bsDim"));
 54:   PetscCall(VecLoad(y, H5viewer));
 55:   PetscCall(VecAXPY(y, -1.0, x));
 56:   PetscCall(VecNorm(y, NORM_2, &norm));
 57:   PetscCall(VecNorm(x, NORM_2, &dnorm));
 58:   PetscCheck(norm / dnorm <= 1.e-6, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Vec read in 'bsDim' does not match vector written out %g", (double)(norm / dnorm));

 60:   PetscCall(PetscViewerDestroy(&H5viewer));
 61:   PetscCall(VecDestroy(&y));
 62:   PetscCall(VecDestroy(&x));
 63:   PetscCall(PetscFinalize());
 64:   return 0;
 65: }

 67: /*TEST

 69:    build:
 70:      requires: hdf5

 72:    test:
 73:      requires: hdf5

 75:    test:
 76:      suffix: 2
 77:      nsize: 4

 79:    test:
 80:      suffix: 3
 81:      nsize: 4
 82:      args: -viewer_hdf5_base_dimension2

 84:    test:
 85:      suffix: 4
 86:      nsize: 4
 87:      args: -viewer_hdf5_sp_output

 89: TEST*/