1: static char help[] = "Parallel HDF5 Vec Viewing.\n\n";

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

  6: int main(int argc, char **argv)
  7: {
  8:   Vec         x1, x2, *x3ts, *x4ts;
  9:   Vec         x1r, x2r, x3r, x4r;
 10:   PetscViewer viewer;
 11:   PetscRandom rand;
 12:   PetscMPIInt rank;
 13:   PetscInt    i, n = 6, n_timesteps = 5;
 14:   PetscBool   equal;
 15:   MPI_Comm    comm;

 17:   PetscFunctionBegin;
 18:   PetscFunctionBeginUser;
 19:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 20:   comm = PETSC_COMM_WORLD;
 21:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
 22:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
 23:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-n_timesteps", &n_timesteps, NULL));
 24:   PetscCheck(n_timesteps >= 0, comm, PETSC_ERR_USER_INPUT, "-n_timesteps must be nonnegative");

 26:   /* create, initialize and write vectors */
 27:   PetscCall(PetscRandomCreate(comm, &rand));
 28:   PetscCall(PetscRandomSetFromOptions(rand));
 29:   PetscCall(PetscViewerHDF5Open(comm, "ex19.h5", FILE_MODE_WRITE, &viewer));

 31:   PetscCall(VecCreate(comm, &x1));
 32:   PetscCall(PetscObjectSetName((PetscObject)x1, "x1"));
 33:   PetscCall(VecSetSizes(x1, PETSC_DECIDE, n));
 34:   PetscCall(VecSetFromOptions(x1));
 35:   PetscCall(VecSetRandom(x1, rand));
 36:   PetscCall(VecView(x1, viewer));

 38:   PetscCall(PetscViewerHDF5PushGroup(viewer, "/testBlockSize"));
 39:   PetscCall(VecCreate(comm, &x2));
 40:   PetscCall(PetscObjectSetName((PetscObject)x2, "x2"));
 41:   PetscCall(VecSetSizes(x2, PETSC_DECIDE, n));
 42:   PetscCall(VecSetBlockSize(x2, 2));
 43:   PetscCall(VecSetFromOptions(x2));
 44:   PetscCall(VecSetRandom(x2, rand));
 45:   PetscCall(VecView(x2, viewer));
 46:   PetscCall(PetscViewerHDF5PopGroup(viewer));

 48:   PetscCall(PetscViewerHDF5PushGroup(viewer, "/testTimestep"));
 49:   PetscCall(PetscViewerHDF5PushTimestepping(viewer));

 51:   PetscCall(VecDuplicateVecs(x1, n_timesteps, &x3ts));
 52:   for (i = 0; i < n_timesteps; i++) {
 53:     PetscCall(PetscObjectSetName((PetscObject)x3ts[i], "x3"));
 54:     PetscCall(VecSetRandom(x3ts[i], rand));
 55:     PetscCall(VecView(x3ts[i], viewer));
 56:     PetscCall(PetscViewerHDF5IncrementTimestep(viewer));
 57:   }

 59:   PetscCall(PetscViewerHDF5PushGroup(viewer, "testBlockSize"));
 60:   PetscCall(VecDuplicateVecs(x2, n_timesteps, &x4ts));
 61:   for (i = 0; i < n_timesteps; i++) {
 62:     PetscCall(PetscObjectSetName((PetscObject)x4ts[i], "x4"));
 63:     PetscCall(VecSetRandom(x4ts[i], rand));
 64:     PetscCall(PetscViewerHDF5SetTimestep(viewer, i));
 65:     PetscCall(VecView(x4ts[i], viewer));
 66:   }
 67:   PetscCall(PetscViewerHDF5PopGroup(viewer));

 69:   PetscCall(PetscViewerHDF5PopTimestepping(viewer));
 70:   PetscCall(PetscViewerHDF5PopGroup(viewer));

 72:   PetscCall(PetscViewerDestroy(&viewer));
 73:   PetscCall(PetscRandomDestroy(&rand));

 75:   /* read and compare */
 76:   PetscCall(PetscViewerHDF5Open(comm, "ex19.h5", FILE_MODE_READ, &viewer));

 78:   PetscCall(VecDuplicate(x1, &x1r));
 79:   PetscCall(PetscObjectSetName((PetscObject)x1r, "x1"));
 80:   PetscCall(VecLoad(x1r, viewer));
 81:   PetscCall(VecEqual(x1, x1r, &equal));
 82:   if (!equal) {
 83:     PetscCall(VecView(x1, PETSC_VIEWER_STDOUT_WORLD));
 84:     PetscCall(VecView(x1r, PETSC_VIEWER_STDOUT_WORLD));
 85:     SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x1 != x1r");
 86:   }

 88:   PetscCall(PetscViewerHDF5PushGroup(viewer, "/testBlockSize"));
 89:   PetscCall(VecDuplicate(x2, &x2r));
 90:   PetscCall(PetscObjectSetName((PetscObject)x2r, "x2"));
 91:   PetscCall(VecLoad(x2r, viewer));
 92:   PetscCall(VecEqual(x2, x2r, &equal));
 93:   if (!equal) {
 94:     PetscCall(VecView(x2, PETSC_VIEWER_STDOUT_WORLD));
 95:     PetscCall(VecView(x2r, PETSC_VIEWER_STDOUT_WORLD));
 96:     SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x2 != x2r");
 97:   }
 98:   PetscCall(PetscViewerHDF5PopGroup(viewer));

100:   PetscCall(PetscViewerHDF5PushGroup(viewer, "/testTimestep"));
101:   PetscCall(PetscViewerHDF5PushTimestepping(viewer));

103:   PetscCall(VecDuplicate(x1, &x3r));
104:   PetscCall(PetscObjectSetName((PetscObject)x3r, "x3"));
105:   for (i = 0; i < n_timesteps; i++) {
106:     PetscCall(PetscViewerHDF5SetTimestep(viewer, i));
107:     PetscCall(VecLoad(x3r, viewer));
108:     PetscCall(VecEqual(x3r, x3ts[i], &equal));
109:     if (!equal) {
110:       PetscCall(VecView(x3r, PETSC_VIEWER_STDOUT_WORLD));
111:       PetscCall(VecView(x3ts[i], PETSC_VIEWER_STDOUT_WORLD));
112:       SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x3ts[%" PetscInt_FMT "] != x3r", i);
113:     }
114:   }

116:   PetscCall(PetscViewerHDF5PushGroup(viewer, "testBlockSize"));
117:   PetscCall(VecDuplicate(x2, &x4r));
118:   PetscCall(PetscObjectSetName((PetscObject)x4r, "x4"));
119:   PetscCall(PetscViewerHDF5SetTimestep(viewer, 0));
120:   for (i = 0; i < n_timesteps; i++) {
121:     PetscCall(VecLoad(x4r, viewer));
122:     PetscCall(VecEqual(x4r, x4ts[i], &equal));
123:     if (!equal) {
124:       PetscCall(VecView(x4r, PETSC_VIEWER_STDOUT_WORLD));
125:       PetscCall(VecView(x4ts[i], PETSC_VIEWER_STDOUT_WORLD));
126:       SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x4ts[%" PetscInt_FMT "] != x4r", i);
127:     }
128:     PetscCall(PetscViewerHDF5IncrementTimestep(viewer));
129:   }
130:   PetscCall(PetscViewerHDF5PopGroup(viewer));

132:   PetscCall(PetscViewerHDF5PopTimestepping(viewer));
133:   PetscCall(PetscViewerHDF5PopGroup(viewer));

135:   /* cleanup */
136:   PetscCall(PetscViewerDestroy(&viewer));
137:   PetscCall(VecDestroy(&x1));
138:   PetscCall(VecDestroy(&x2));
139:   PetscCall(VecDestroyVecs(n_timesteps, &x3ts));
140:   PetscCall(VecDestroyVecs(n_timesteps, &x4ts));
141:   PetscCall(VecDestroy(&x1r));
142:   PetscCall(VecDestroy(&x2r));
143:   PetscCall(VecDestroy(&x3r));
144:   PetscCall(VecDestroy(&x4r));
145:   PetscCall(PetscFinalize());
146:   return 0;
147: }

149: /*TEST

151:      build:
152:        requires: hdf5

154:      test:
155:        nsize: {{1 2 3 4}}

157: TEST*/