Actual source code: hdf5v.c

  1: #include <petsc/private/viewerhdf5impl.h>

  3: static PetscErrorCode PetscViewerHDF5Traverse_Internal(PetscViewer, const char[], PetscBool, PetscBool *, H5O_type_t *);
  4: static PetscErrorCode PetscViewerHDF5HasAttribute_Internal(PetscViewer, const char[], const char[], PetscBool *);

  6: /*@C
  7:   PetscViewerHDF5GetGroup - Get the current HDF5 group name (full path), set with `PetscViewerHDF5PushGroup()`/`PetscViewerHDF5PopGroup()`.

  9:   Not collective

 11:   Input Parameters:
 12: + viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
 13: - path - (Optional) The path relative to the pushed group

 15:   Output Parameter:
 16: . abspath - The absolute HDF5 path (group)

 18:   Level: intermediate

 20:   Notes:
 21:   If path starts with '/', it is taken as an absolute path overriding currently pushed group, else path is relative to the current pushed group.
 22:   So NULL or empty path means the current pushed group.

 24:   The output abspath is newly allocated so needs to be freed.

 26: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5OpenGroup()`, `PetscViewerHDF5WriteGroup()`
 27: @*/
 28: PetscErrorCode PetscViewerHDF5GetGroup(PetscViewer viewer, const char path[], char *abspath[])
 29: {
 30:   size_t      len;
 31:   PetscBool   relative = PETSC_FALSE;
 32:   const char *group;
 33:   char        buf[PETSC_MAX_PATH_LEN] = "";

 38:   PetscViewerHDF5GetGroup_Internal(viewer, &group);
 39:   PetscStrlen(path, &len);
 40:   relative = (PetscBool)(!len || path[0] != '/');
 41:   if (relative) {
 42:     PetscStrcpy(buf, group);
 43:     if (!group || len) PetscStrcat(buf, "/");
 44:     PetscStrcat(buf, path);
 45:     PetscStrallocpy(buf, abspath);
 46:   } else {
 47:     PetscStrallocpy(path, abspath);
 48:   }
 49:   return 0;
 50: }

 52: static PetscErrorCode PetscViewerHDF5CheckNamedObject_Internal(PetscViewer viewer, PetscObject obj)
 53: {
 54:   PetscBool has;

 56:   PetscViewerHDF5HasObject(viewer, obj, &has);
 57:   if (!has) {
 58:     char *group;
 59:     PetscViewerHDF5GetGroup(viewer, NULL, &group);
 60:     SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Object (dataset) \"%s\" not stored in group %s", obj->name, group);
 61:   }
 62:   return 0;
 63: }

 65: static PetscErrorCode PetscViewerSetFromOptions_HDF5(PetscViewer v, PetscOptionItems *PetscOptionsObject)
 66: {
 67:   PetscBool         flg  = PETSC_FALSE, set;
 68:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)v->data;

 70:   PetscOptionsHeadBegin(PetscOptionsObject, "HDF5 PetscViewer Options");
 71:   PetscOptionsBool("-viewer_hdf5_base_dimension2", "1d Vectors get 2 dimensions in HDF5", "PetscViewerHDF5SetBaseDimension2", hdf5->basedimension2, &hdf5->basedimension2, NULL);
 72:   PetscOptionsBool("-viewer_hdf5_sp_output", "Force data to be written in single precision", "PetscViewerHDF5SetSPOutput", hdf5->spoutput, &hdf5->spoutput, NULL);
 73:   PetscOptionsBool("-viewer_hdf5_collective", "Enable collective transfer mode", "PetscViewerHDF5SetCollective", flg, &flg, &set);
 74:   if (set) PetscViewerHDF5SetCollective(v, flg);
 75:   flg = PETSC_FALSE;
 76:   PetscOptionsBool("-viewer_hdf5_default_timestepping", "Set default timestepping state", "PetscViewerHDF5SetDefaultTimestepping", flg, &flg, &set);
 77:   if (set) PetscViewerHDF5SetDefaultTimestepping(v, flg);
 78:   PetscOptionsHeadEnd();
 79:   return 0;
 80: }

 82: static PetscErrorCode PetscViewerView_HDF5(PetscViewer v, PetscViewer viewer)
 83: {
 84:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)v->data;
 85:   PetscBool         flg;

 87:   if (hdf5->filename) PetscViewerASCIIPrintf(viewer, "Filename: %s\n", hdf5->filename);
 88:   PetscViewerASCIIPrintf(viewer, "Vectors with blocksize 1 saved as 2D datasets: %s\n", PetscBools[hdf5->basedimension2]);
 89:   PetscViewerASCIIPrintf(viewer, "Enforce single precision storage: %s\n", PetscBools[hdf5->spoutput]);
 90:   PetscViewerHDF5GetCollective(v, &flg);
 91:   PetscViewerASCIIPrintf(viewer, "MPI-IO transfer mode: %s\n", flg ? "collective" : "independent");
 92:   PetscViewerASCIIPrintf(viewer, "Default timestepping: %s\n", PetscBools[hdf5->defTimestepping]);
 93:   return 0;
 94: }

 96: static PetscErrorCode PetscViewerFileClose_HDF5(PetscViewer viewer)
 97: {
 98:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

100:   PetscFree(hdf5->filename);
101:   if (hdf5->file_id) PetscCallHDF5(H5Fclose, (hdf5->file_id));
102:   return 0;
103: }

105: static PetscErrorCode PetscViewerFlush_HDF5(PetscViewer viewer)
106: {
107:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

109:   if (hdf5->file_id) PetscCallHDF5(H5Fflush, (hdf5->file_id, H5F_SCOPE_LOCAL));
110:   return 0;
111: }

113: static PetscErrorCode PetscViewerDestroy_HDF5(PetscViewer viewer)
114: {
115:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

117:   PetscCallHDF5(H5Pclose, (hdf5->dxpl_id));
118:   PetscViewerFileClose_HDF5(viewer);
119:   while (hdf5->groups) {
120:     PetscViewerHDF5GroupList *tmp = hdf5->groups->next;

122:     PetscFree(hdf5->groups->name);
123:     PetscFree(hdf5->groups);
124:     hdf5->groups = tmp;
125:   }
126:   PetscFree(hdf5);
127:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileSetName_C", NULL);
128:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileGetName_C", NULL);
129:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileSetMode_C", NULL);
130:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileGetMode_C", NULL);
131:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5SetBaseDimension2_C", NULL);
132:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5SetSPOutput_C", NULL);
133:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5SetCollective_C", NULL);
134:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5GetCollective_C", NULL);
135:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5GetDefaultTimestepping_C", NULL);
136:   PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5SetDefaultTimestepping_C", NULL);
137:   return 0;
138: }

140: static PetscErrorCode PetscViewerFileSetMode_HDF5(PetscViewer viewer, PetscFileMode type)
141: {
142:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

144:   hdf5->btype = type;
145:   return 0;
146: }

148: static PetscErrorCode PetscViewerFileGetMode_HDF5(PetscViewer viewer, PetscFileMode *type)
149: {
150:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

152:   *type = hdf5->btype;
153:   return 0;
154: }

156: static PetscErrorCode PetscViewerHDF5SetBaseDimension2_HDF5(PetscViewer viewer, PetscBool flg)
157: {
158:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

160:   hdf5->basedimension2 = flg;
161:   return 0;
162: }

164: /*@
165:      PetscViewerHDF5SetBaseDimension2 - Vectors of 1 dimension (i.e. bs/dof is 1) will be saved in the HDF5 file with a
166:        dimension of 2.

168:     Logically Collective

170:   Input Parameters:
171: +  viewer - the `PetscViewer`; if it is a `PETSCVIEWERHDF5` then this command is ignored
172: -  flg - if `PETSC_TRUE` the vector will always have at least a dimension of 2 even if that first dimension is of size 1

174:   Options Database Key:
175: .  -viewer_hdf5_base_dimension2 - turns on (true) or off (false) using a dimension of 2 in the HDF5 file even if the bs/dof of the vector is 1

177:   Note:
178:   Setting this option allegedly makes code that reads the HDF5 in easier since they do not have a "special case" of a bs/dof
179:   of one when the dimension is lower. Others think the option is crazy.

181:   Level: intermediate

183: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, PetscViewerFileSetMode()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerBinaryOpen()`
184: @*/
185: PetscErrorCode PetscViewerHDF5SetBaseDimension2(PetscViewer viewer, PetscBool flg)
186: {
188:   PetscTryMethod(viewer, "PetscViewerHDF5SetBaseDimension2_C", (PetscViewer, PetscBool), (viewer, flg));
189:   return 0;
190: }

192: /*@
193:      PetscViewerHDF5GetBaseDimension2 - Vectors of 1 dimension (i.e. bs/dof is 1) will be saved in the HDF5 file with a
194:        dimension of 2.

196:     Logically Collective

198:   Input Parameter:
199: .  viewer - the `PetscViewer`, must be `PETSCVIEWERHDF5`

201:   Output Parameter:
202: .  flg - if `PETSC_TRUE` the vector will always have at least a dimension of 2 even if that first dimension is of size 1

204:   Note:
205:   Setting this option allegedly makes code that reads the HDF5 in easier since they do not have a "special case" of a bs/dof
206:   of one when the dimension is lower. Others think the option is crazy.

208:   Level: intermediate

210: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerFileSetMode()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerBinaryOpen()`
211: @*/
212: PetscErrorCode PetscViewerHDF5GetBaseDimension2(PetscViewer viewer, PetscBool *flg)
213: {
214:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

217:   *flg = hdf5->basedimension2;
218:   return 0;
219: }

221: static PetscErrorCode PetscViewerHDF5SetSPOutput_HDF5(PetscViewer viewer, PetscBool flg)
222: {
223:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

225:   hdf5->spoutput = flg;
226:   return 0;
227: }

229: /*@
230:      PetscViewerHDF5SetSPOutput - Data is written to disk in single precision even if PETSc is
231:        compiled with double precision `PetscReal`.

233:     Logically Collective

235:   Input Parameters:
236: +  viewer - the PetscViewer; if it is a `PETSCVIEWERHDF5` then this command is ignored
237: -  flg - if `PETSC_TRUE` the data will be written to disk with single precision

239:   Options Database Key:
240: .  -viewer_hdf5_sp_output - turns on (true) or off (false) output in single precision

242:   Note:
243:     Setting this option does not make any difference if PETSc is compiled with single precision
244:          in the first place. It does not affect reading datasets (HDF5 handle this internally).

246:   Level: intermediate

248: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerFileSetMode()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerBinaryOpen()`,
249:           `PetscReal`, `PetscViewerHDF5GetSPOutput()`
250: @*/
251: PetscErrorCode PetscViewerHDF5SetSPOutput(PetscViewer viewer, PetscBool flg)
252: {
254:   PetscTryMethod(viewer, "PetscViewerHDF5SetSPOutput_C", (PetscViewer, PetscBool), (viewer, flg));
255:   return 0;
256: }

258: /*@
259:      PetscViewerHDF5GetSPOutput - Data is written to disk in single precision even if PETSc is
260:        compiled with double precision `PetscReal`.

262:     Logically Collective

264:   Input Parameter:
265: .  viewer - the PetscViewer, must be of type `PETSCVIEWERHDF5`

267:   Output Parameter:
268: .  flg - if `PETSC_TRUE` the data will be written to disk with single precision

270:   Level: intermediate

272: .seealso: [](sec_viewers), `PetscViewerFileSetMode()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerBinaryOpen()`,
273:           `PetscReal`, `PetscViewerHDF5SetSPOutput()`
274: @*/
275: PetscErrorCode PetscViewerHDF5GetSPOutput(PetscViewer viewer, PetscBool *flg)
276: {
277:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

280:   *flg = hdf5->spoutput;
281:   return 0;
282: }

284: static PetscErrorCode PetscViewerHDF5SetCollective_HDF5(PetscViewer viewer, PetscBool flg)
285: {
286:   /* H5FD_MPIO_COLLECTIVE is wrong in hdf5 1.10.2, and is the same as H5FD_MPIO_INDEPENDENT in earlier versions
287:      - see e.g. https://gitlab.cosma.dur.ac.uk/swift/swiftsim/issues/431 */
288: #if H5_VERSION_GE(1, 10, 3) && defined(H5_HAVE_PARALLEL)
289:   {
290:     PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
291:     PetscCallHDF5(H5Pset_dxpl_mpio, (hdf5->dxpl_id, flg ? H5FD_MPIO_COLLECTIVE : H5FD_MPIO_INDEPENDENT));
292:   }
293: #else
294:   if (flg) PetscPrintf(PetscObjectComm((PetscObject)viewer), "Warning: PetscViewerHDF5SetCollective(viewer,PETSC_TRUE) is ignored for HDF5 versions prior to 1.10.3 or if built without MPI support\n");
295: #endif
296:   return 0;
297: }

299: /*@
300:   PetscViewerHDF5SetCollective - Use collective MPI-IO transfer mode for HDF5 reads and writes.

302:   Logically Collective; flg must contain common value

304:   Input Parameters:
305: + viewer - the `PetscViewer`; if it is not `PETSCVIEWERHDF5` then this command is ignored
306: - flg - `PETSC_TRUE` for collective mode; `PETSC_FALSE` for independent mode (default)

308:   Options Database Key:
309: . -viewer_hdf5_collective - turns on (true) or off (false) collective transfers

311:   Note:
312:   Collective mode gives the MPI-IO layer underneath HDF5 a chance to do some additional collective optimizations and hence can perform better.
313:   However, this works correctly only since HDF5 1.10.3 and if HDF5 is installed for MPI; hence, we ignore this setting for older versions.

315:   Developer Note:
316:   In the HDF5 layer, `PETSC_TRUE` / `PETSC_FALSE` means `H5Pset_dxpl_mpio()` is called with `H5FD_MPIO_COLLECTIVE` / `H5FD_MPIO_INDEPENDENT`, respectively.
317:   This in turn means use of MPI_File_{read,write}_all /  MPI_File_{read,write} in the MPI-IO layer, respectively.
318:   See HDF5 documentation and MPI-IO documentation for details.

320:   Level: intermediate

322: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5GetCollective()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerHDF5Open()`
323: @*/
324: PetscErrorCode PetscViewerHDF5SetCollective(PetscViewer viewer, PetscBool flg)
325: {
328:   PetscTryMethod(viewer, "PetscViewerHDF5SetCollective_C", (PetscViewer, PetscBool), (viewer, flg));
329:   return 0;
330: }

332: static PetscErrorCode PetscViewerHDF5GetCollective_HDF5(PetscViewer viewer, PetscBool *flg)
333: {
334: #if defined(H5_HAVE_PARALLEL)
335:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
336:   H5FD_mpio_xfer_t  mode;
337: #endif

339: #if !defined(H5_HAVE_PARALLEL)
340:   *flg = PETSC_FALSE;
341: #else
342:   PetscCallHDF5(H5Pget_dxpl_mpio, (hdf5->dxpl_id, &mode));
343:   *flg = (mode == H5FD_MPIO_COLLECTIVE) ? PETSC_TRUE : PETSC_FALSE;
344: #endif
345:   return 0;
346: }

348: /*@
349:   PetscViewerHDF5GetCollective - Return flag whether collective MPI-IO transfer mode is used for HDF5 reads and writes.

351:   Not Collective

353:   Input Parameters:
354: . viewer - the `PETSCVIEWERHDF5` `PetscViewer`

356:   Output Parameters:
357: . flg - the flag

359:   Level: intermediate

361:   Note:
362:   This setting works correctly only since HDF5 1.10.3 and if HDF5 was installed for MPI. For older versions, `PETSC_FALSE` will be always returned.
363:   For more details, see `PetscViewerHDF5SetCollective()`.

365: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5SetCollective()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerHDF5Open()`
366: @*/
367: PetscErrorCode PetscViewerHDF5GetCollective(PetscViewer viewer, PetscBool *flg)
368: {

372:   PetscUseMethod(viewer, "PetscViewerHDF5GetCollective_C", (PetscViewer, PetscBool *), (viewer, flg));
373:   return 0;
374: }

376: static PetscErrorCode PetscViewerFileSetName_HDF5(PetscViewer viewer, const char name[])
377: {
378:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
379:   hid_t             plist_id;

381:   if (hdf5->file_id) PetscCallHDF5(H5Fclose, (hdf5->file_id));
382:   if (hdf5->filename) PetscFree(hdf5->filename);
383:   PetscStrallocpy(name, &hdf5->filename);
384:   /* Set up file access property list with parallel I/O access */
385:   PetscCallHDF5Return(plist_id, H5Pcreate, (H5P_FILE_ACCESS));
386: #if defined(H5_HAVE_PARALLEL)
387:   PetscCallHDF5(H5Pset_fapl_mpio, (plist_id, PetscObjectComm((PetscObject)viewer), MPI_INFO_NULL));
388: #endif
389:   /* Create or open the file collectively */
390:   switch (hdf5->btype) {
391:   case FILE_MODE_READ:
392:     if (PetscDefined(USE_DEBUG)) {
393:       PetscMPIInt rank;
394:       PetscBool   flg;

396:       MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank);
397:       if (rank == 0) {
398:         PetscTestFile(hdf5->filename, 'r', &flg);
400:       }
401:       MPI_Barrier(PetscObjectComm((PetscObject)viewer));
402:     }
403:     PetscCallHDF5Return(hdf5->file_id, H5Fopen, (name, H5F_ACC_RDONLY, plist_id));
404:     break;
405:   case FILE_MODE_APPEND:
406:   case FILE_MODE_UPDATE: {
407:     PetscBool flg;
408:     PetscTestFile(hdf5->filename, 'r', &flg);
409:     if (flg) PetscCallHDF5Return(hdf5->file_id, H5Fopen, (name, H5F_ACC_RDWR, plist_id));
410:     else PetscCallHDF5Return(hdf5->file_id, H5Fcreate, (name, H5F_ACC_EXCL, H5P_DEFAULT, plist_id));
411:     break;
412:   }
413:   case FILE_MODE_WRITE:
414:     PetscCallHDF5Return(hdf5->file_id, H5Fcreate, (name, H5F_ACC_TRUNC, H5P_DEFAULT, plist_id));
415:     break;
416:   case FILE_MODE_UNDEFINED:
417:     SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
418:   default:
419:     SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Unsupported file mode %s", PetscFileModes[hdf5->btype]);
420:   }
422:   PetscCallHDF5(H5Pclose, (plist_id));
423:   PetscViewerHDF5ResetAttachedDMPlexStorageVersion(viewer);
424:   return 0;
425: }

427: static PetscErrorCode PetscViewerFileGetName_HDF5(PetscViewer viewer, const char **name)
428: {
429:   PetscViewer_HDF5 *vhdf5 = (PetscViewer_HDF5 *)viewer->data;

431:   *name = vhdf5->filename;
432:   return 0;
433: }

435: static PetscErrorCode PetscViewerSetUp_HDF5(PetscViewer viewer)
436: {
437:   /*
438:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
439:   */

441:   return 0;
442: }

444: static PetscErrorCode PetscViewerHDF5SetDefaultTimestepping_HDF5(PetscViewer viewer, PetscBool flg)
445: {
446:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

448:   hdf5->defTimestepping = flg;
449:   return 0;
450: }

452: static PetscErrorCode PetscViewerHDF5GetDefaultTimestepping_HDF5(PetscViewer viewer, PetscBool *flg)
453: {
454:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

456:   *flg = hdf5->defTimestepping;
457:   return 0;
458: }

460: /*@
461:   PetscViewerHDF5SetDefaultTimestepping - Set the flag for default timestepping

463:   Logically Collective

465:   Input Parameters:
466: + viewer - the `PetscViewer`; if it is not `PETSCVIEWERHDF5` then this command is ignored
467: - flg    - if `PETSC_TRUE` we will assume that timestepping is on

469:   Options Database Key:
470: . -viewer_hdf5_default_timestepping - turns on (true) or off (false) default timestepping

472:   Note:
473:   If the timestepping attribute is not found for an object, then the default timestepping is used

475:   Level: intermediate

477: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5GetDefaultTimestepping()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5GetTimestep()`
478: @*/
479: PetscErrorCode PetscViewerHDF5SetDefaultTimestepping(PetscViewer viewer, PetscBool flg)
480: {
482:   PetscTryMethod(viewer, "PetscViewerHDF5SetDefaultTimestepping_C", (PetscViewer, PetscBool), (viewer, flg));
483:   return 0;
484: }

486: /*@
487:   PetscViewerHDF5GetDefaultTimestepping - Get the flag for default timestepping

489:   Not collective

491:   Input Parameter:
492: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`

494:   Output Parameter:
495: . flg    - if `PETSC_TRUE` we will assume that timestepping is on

497:   Level: intermediate

499: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5SetDefaultTimestepping()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5GetTimestep()`
500: @*/
501: PetscErrorCode PetscViewerHDF5GetDefaultTimestepping(PetscViewer viewer, PetscBool *flg)
502: {
504:   PetscUseMethod(viewer, "PetscViewerHDF5GetDefaultTimestepping_C", (PetscViewer, PetscBool *), (viewer, flg));
505:   return 0;
506: }

508: /*MC
509:    PETSCVIEWERHDF5 - A viewer that writes to an HDF5 file

511:   Level: beginner

513: .seealso: [](sec_viewers), `PetscViewerHDF5Open()`, `PetscViewerStringSPrintf()`, `PetscViewerSocketOpen()`, `PetscViewerDrawOpen()`, `PETSCVIEWERSOCKET`,
514:           `PetscViewerCreate()`, `PetscViewerASCIIOpen()`, `PetscViewerBinaryOpen()`, `PETSCVIEWERBINARY`, `PETSCVIEWERDRAW`, `PETSCVIEWERSTRING`,
515:           `PetscViewerMatlabOpen()`, `VecView()`, `DMView()`, `PetscViewerMatlabPutArray()`, `PETSCVIEWERASCII`, `PETSCVIEWERMATLAB`,
516:           `PetscViewerFileSetName()`, `PetscViewerFileSetMode()`, `PetscViewerFormat`, `PetscViewerType`, `PetscViewerSetType()`
517: M*/

519: PETSC_EXTERN PetscErrorCode PetscViewerCreate_HDF5(PetscViewer v)
520: {
521:   PetscViewer_HDF5 *hdf5;

523: #if !defined(H5_HAVE_PARALLEL)
524:   {
525:     PetscMPIInt size;
526:     MPI_Comm_size(PetscObjectComm((PetscObject)v), &size);
528:   }
529: #endif

531:   PetscNew(&hdf5);

533:   v->data                = (void *)hdf5;
534:   v->ops->destroy        = PetscViewerDestroy_HDF5;
535:   v->ops->setfromoptions = PetscViewerSetFromOptions_HDF5;
536:   v->ops->setup          = PetscViewerSetUp_HDF5;
537:   v->ops->view           = PetscViewerView_HDF5;
538:   v->ops->flush          = PetscViewerFlush_HDF5;
539:   hdf5->btype            = FILE_MODE_UNDEFINED;
540:   hdf5->filename         = NULL;
541:   hdf5->timestep         = -1;
542:   hdf5->groups           = NULL;

544:   PetscCallHDF5Return(hdf5->dxpl_id, H5Pcreate, (H5P_DATASET_XFER));

546:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileSetName_C", PetscViewerFileSetName_HDF5);
547:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileGetName_C", PetscViewerFileGetName_HDF5);
548:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileSetMode_C", PetscViewerFileSetMode_HDF5);
549:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileGetMode_C", PetscViewerFileGetMode_HDF5);
550:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5SetBaseDimension2_C", PetscViewerHDF5SetBaseDimension2_HDF5);
551:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5SetSPOutput_C", PetscViewerHDF5SetSPOutput_HDF5);
552:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5SetCollective_C", PetscViewerHDF5SetCollective_HDF5);
553:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5GetCollective_C", PetscViewerHDF5GetCollective_HDF5);
554:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5GetDefaultTimestepping_C", PetscViewerHDF5GetDefaultTimestepping_HDF5);
555:   PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5SetDefaultTimestepping_C", PetscViewerHDF5SetDefaultTimestepping_HDF5);
556:   return 0;
557: }

559: /*@C
560:    PetscViewerHDF5Open - Opens a file for HDF5 input/output as a `PETSCVIEWERHDF5` `PetscViewer`

562:    Collective

564:    Input Parameters:
565: +  comm - MPI communicator
566: .  name - name of file
567: -  type - type of file

569:    Output Parameter:
570: .  hdf5v - `PetscViewer` for HDF5 input/output to use with the specified file

572:   Options Database Keys:
573: +  -viewer_hdf5_base_dimension2 - turns on (true) or off (false) using a dimension of 2 in the HDF5 file even if the bs/dof of the vector is 1
574: -  -viewer_hdf5_sp_output - forces (if true) the viewer to write data in single precision independent on the precision of PetscReal

576:    Level: beginner

578:    Notes:
579:    Reading is always available, regardless of the mode. Available modes are
580: .vb
581:   FILE_MODE_READ - open existing HDF5 file for read only access, fail if file does not exist [H5Fopen() with H5F_ACC_RDONLY]
582:   FILE_MODE_WRITE - if file exists, fully overwrite it, else create new HDF5 file [H5FcreateH5Fcreate() with H5F_ACC_TRUNC]
583:   FILE_MODE_APPEND - if file exists, keep existing contents [H5Fopen() with H5F_ACC_RDWR], else create new HDF5 file [H5FcreateH5Fcreate() with H5F_ACC_EXCL]
584:   FILE_MODE_UPDATE - same as FILE_MODE_APPEND
585: .ve

587:    In case of `FILE_MODE_APPEND` / `FILE_MODE_UPDATE`, any stored object (dataset, attribute) can be selectively overwritten if the same fully qualified name (/group/path/to/object) is specified.

589:    This PetscViewer should be destroyed with PetscViewerDestroy().

591: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerASCIIOpen()`, `PetscViewerPushFormat()`, `PetscViewerDestroy()`, `PetscViewerHDF5SetBaseDimension2()`,
592:           `PetscViewerHDF5SetSPOutput()`, `PetscViewerHDF5GetBaseDimension2()`, `VecView()`, `MatView()`, `VecLoad()`,
593:           `MatLoad()`, `PetscFileMode`, `PetscViewer`, `PetscViewerSetType()`, `PetscViewerFileSetMode()`, `PetscViewerFileSetName()`
594: @*/
595: PetscErrorCode PetscViewerHDF5Open(MPI_Comm comm, const char name[], PetscFileMode type, PetscViewer *hdf5v)
596: {
597:   PetscViewerCreate(comm, hdf5v);
598:   PetscViewerSetType(*hdf5v, PETSCVIEWERHDF5);
599:   PetscViewerFileSetMode(*hdf5v, type);
600:   PetscViewerFileSetName(*hdf5v, name);
601:   PetscViewerSetFromOptions(*hdf5v);
602:   return 0;
603: }

605: /*@C
606:   PetscViewerHDF5GetFileId - Retrieve the file id, this file ID then can be used in direct HDF5 calls

608:   Not collective

610:   Input Parameter:
611: . viewer - the `PetscViewer`

613:   Output Parameter:
614: . file_id - The file id

616:   Level: intermediate

618: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`
619: @*/
620: PetscErrorCode PetscViewerHDF5GetFileId(PetscViewer viewer, hid_t *file_id)
621: {
622:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

625:   if (file_id) *file_id = hdf5->file_id;
626:   return 0;
627: }

629: /*@C
630:   PetscViewerHDF5PushGroup - Set the current HDF5 group for output

632:   Not collective

634:   Input Parameters:
635: + viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
636: - name - The group name

638:   Level: intermediate

640:   Notes:
641:   This is designed to mnemonically resemble the Unix cd command.
642: .vb
643:   If name begins with '/', it is interpreted as an absolute path fully replacing current group, otherwise it is taken as relative to the current group.
644:   NULL, empty string, or any sequence of all slashes (e.g. "///") is interpreted as the root group "/".
645:   "." means the current group is pushed again.
646: .ve

648:   Example:
649:   Suppose the current group is "/a".
650:   + If name is NULL, empty string, or a sequence of all slashes (e.g. "///"), then the new group will be "/".
651:   . If name is ".", then the new group will be "/a".
652:   . If name is "b", then the new group will be "/a/b".
653:   - If name is "/b", then the new group will be "/b".

655:   Developer Note:
656:   The root group "/" is internally stored as NULL.

658: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5OpenGroup()`, `PetscViewerHDF5WriteGroup()`
659: @*/
660: PetscErrorCode PetscViewerHDF5PushGroup(PetscViewer viewer, const char name[])
661: {
662:   PetscViewer_HDF5         *hdf5 = (PetscViewer_HDF5 *)viewer->data;
663:   PetscViewerHDF5GroupList *groupNode;
664:   size_t                    i, len;
665:   char                      buf[PETSC_MAX_PATH_LEN];
666:   const char               *gname;

670:   PetscStrlen(name, &len);
671:   gname = NULL;
672:   if (len) {
673:     if (len == 1 && name[0] == '.') {
674:       /* use current name */
675:       gname = (hdf5->groups && hdf5->groups->name) ? hdf5->groups->name : NULL;
676:     } else if (name[0] == '/') {
677:       /* absolute */
678:       for (i = 1; i < len; i++) {
679:         if (name[i] != '/') {
680:           gname = name;
681:           break;
682:         }
683:       }
684:     } else {
685:       /* relative */
686:       const char *parent = (hdf5->groups && hdf5->groups->name) ? hdf5->groups->name : "";
687:       PetscSNPrintf(buf, sizeof(buf), "%s/%s", parent, name);
688:       gname = buf;
689:     }
690:   }
691:   PetscNew(&groupNode);
692:   PetscStrallocpy(gname, (char **)&groupNode->name);
693:   groupNode->next = hdf5->groups;
694:   hdf5->groups    = groupNode;
695:   return 0;
696: }

698: /*@
699:   PetscViewerHDF5PopGroup - Return the current HDF5 group for output to the previous value

701:   Not collective

703:   Input Parameter:
704: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`

706:   Level: intermediate

708: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5OpenGroup()`, `PetscViewerHDF5WriteGroup()`
709: @*/
710: PetscErrorCode PetscViewerHDF5PopGroup(PetscViewer viewer)
711: {
712:   PetscViewer_HDF5         *hdf5 = (PetscViewer_HDF5 *)viewer->data;
713:   PetscViewerHDF5GroupList *groupNode;

717:   groupNode    = hdf5->groups;
718:   hdf5->groups = hdf5->groups->next;
719:   PetscFree(groupNode->name);
720:   PetscFree(groupNode);
721:   return 0;
722: }

724: PetscErrorCode PetscViewerHDF5GetGroup_Internal(PetscViewer viewer, const char *name[])
725: {
726:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

730:   if (hdf5->groups) *name = hdf5->groups->name;
731:   else *name = NULL;
732:   return 0;
733: }

735: /*@C
736:   PetscViewerHDF5OpenGroup - Open the HDF5 group with the name (full path) returned by `PetscViewerHDF5GetGroup()`,
737:   and return this group's ID and file ID.
738:   If `PetscViewerHDF5GetGroup()` yields NULL, then group ID is file ID.

740:   Not collective

742:   Input Parameters:
743: + viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
744: - path - (Optional) The path relative to the pushed group

746:   Output Parameters:
747: + fileId - The HDF5 file ID
748: - groupId - The HDF5 group ID

750:   Note:
751:   If path starts with '/', it is taken as an absolute path overriding currently pushed group, else path is relative to the current pushed group.
752:   So NULL or empty path means the current pushed group.

754:   If the viewer is writable, the group is created if it doesn't exist yet.

756:   Level: intermediate

758: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5WriteGroup()`
759: @*/
760: PetscErrorCode PetscViewerHDF5OpenGroup(PetscViewer viewer, const char path[], hid_t *fileId, hid_t *groupId)
761: {
762:   hid_t       file_id;
763:   H5O_type_t  type;
764:   const char *fileName  = NULL;
765:   char       *groupName = NULL;
766:   PetscBool   writable, has;

768:   PetscViewerWritable(viewer, &writable);
769:   PetscViewerHDF5GetFileId(viewer, &file_id);
770:   PetscViewerFileGetName(viewer, &fileName);
771:   PetscViewerHDF5GetGroup(viewer, path, &groupName);
772:   PetscViewerHDF5Traverse_Internal(viewer, groupName, writable, &has, &type);
773:   if (!has) {
775:     SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_LIB, "HDF5 failed to create group %s although file %s is open for writing", groupName, fileName);
776:   }
778:   PetscCallHDF5Return(*groupId, H5Gopen2, (file_id, groupName, H5P_DEFAULT));
779:   PetscFree(groupName);
780:   *fileId = file_id;
781:   return 0;
782: }

784: /*@C
785:   PetscViewerHDF5WriteGroup - Ensure the HDF5 group exists in the HDF5 file

787:   Not collective

789:   Input Parameters:
790: + viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
791: - path - (Optional) The path relative to the pushed group

793:   Note:
794:   If path starts with '/', it is taken as an absolute path overriding currently pushed group, else path is relative to the current pushed group.
795:   So NULL or empty path means the current pushed group.

797:   This will fail if the viewer is not writable.

799:   Level: intermediate

801: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5OpenGroup()`
802: @*/
803: PetscErrorCode PetscViewerHDF5WriteGroup(PetscViewer viewer, const char path[])
804: {
805:   hid_t fileId, groupId;

807:   PetscViewerHDF5OpenGroup(viewer, path, &fileId, &groupId); // make sure group is actually created
808:   PetscCallHDF5(H5Gclose, (groupId));
809:   return 0;
810: }

812: /*@
813:   PetscViewerHDF5PushTimestepping - Activate timestepping mode for subsequent HDF5 reading and writing.

815:   Not collective

817:   Input Parameter:
818: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`

820:   Level: intermediate

822:   Notes:
823:   On first `PetscViewerHDF5PushTimestepping()`, the initial time step is set to 0.
824:   Next timesteps can then be set using `PetscViewerHDF5IncrementTimestep()` or `PetscViewerHDF5SetTimestep()`.
825:   Current timestep value determines which timestep is read from or written to any dataset on the next HDF5 I/O operation [e.g. `VecView()`].
826:   Use `PetscViewerHDF5PopTimestepping()` to deactivate timestepping mode; calling it by the end of the program is NOT mandatory.
827:   Current timestep is remembered between `PetscViewerHDF5PopTimestepping()` and the next `PetscViewerHDF5PushTimestepping()`.

829:   If a dataset was stored with timestepping, it can be loaded only in the timestepping mode again.
830:   Loading a timestepped dataset with timestepping disabled, or vice-versa results in an error.

832:   Developer note:
833:   Timestepped HDF5 dataset has an extra dimension and attribute "timestepping" set to true.

835: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PopTimestepping()`, `PetscViewerHDF5IsTimestepping()`, `PetscViewerHDF5SetTimestep()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5GetTimestep()`
836: @*/
837: PetscErrorCode PetscViewerHDF5PushTimestepping(PetscViewer viewer)
838: {
839:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

843:   hdf5->timestepping = PETSC_TRUE;
844:   if (hdf5->timestep < 0) hdf5->timestep = 0;
845:   return 0;
846: }

848: /*@
849:   PetscViewerHDF5PopTimestepping - Deactivate timestepping mode for subsequent HDF5 reading and writing.

851:   Not collective

853:   Input Parameter:
854: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`

856:   Level: intermediate

858:   Note:
859:   See `PetscViewerHDF5PushTimestepping()` for details.

861: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5IsTimestepping()`, `PetscViewerHDF5SetTimestep()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5GetTimestep()`
862: @*/
863: PetscErrorCode PetscViewerHDF5PopTimestepping(PetscViewer viewer)
864: {
865:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

869:   hdf5->timestepping = PETSC_FALSE;
870:   return 0;
871: }

873: /*@
874:   PetscViewerHDF5IsTimestepping - Ask the viewer whether it is in timestepping mode currently.

876:   Not collective

878:   Input Parameter:
879: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`

881:   Output Parameter:
882: . flg - is timestepping active?

884:   Level: intermediate

886:   Note:
887:   See `PetscViewerHDF5PushTimestepping()` for details.

889: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5PopTimestepping()`, `PetscViewerHDF5SetTimestep()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5GetTimestep()`
890: @*/
891: PetscErrorCode PetscViewerHDF5IsTimestepping(PetscViewer viewer, PetscBool *flg)
892: {
893:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

896:   *flg = hdf5->timestepping;
897:   return 0;
898: }

900: /*@
901:   PetscViewerHDF5IncrementTimestep - Increments current timestep for the HDF5 output. Fields are stacked in time.

903:   Not collective

905:   Input Parameter:
906: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`

908:   Level: intermediate

910:   Note:
911:   This can be called only if the viewer is in timestepping mode. See `PetscViewerHDF5PushTimestepping()` for details.

913: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5SetTimestep()`, `PetscViewerHDF5GetTimestep()`
914: @*/
915: PetscErrorCode PetscViewerHDF5IncrementTimestep(PetscViewer viewer)
916: {
917:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

921:   ++hdf5->timestep;
922:   return 0;
923: }

925: /*@
926:   PetscViewerHDF5SetTimestep - Set the current timestep for the HDF5 output. Fields are stacked in time.

928:   Logically collective

930:   Input Parameters:
931: + viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
932: - timestep - The timestep

934:   Level: intermediate

936:   Note:
937:   This can be called only if the viewer is in timestepping mode. See `PetscViewerHDF5PushTimestepping()` for details.

939: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5GetTimestep()`
940: @*/
941: PetscErrorCode PetscViewerHDF5SetTimestep(PetscViewer viewer, PetscInt timestep)
942: {
943:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

949:   hdf5->timestep = timestep;
950:   return 0;
951: }

953: /*@
954:   PetscViewerHDF5GetTimestep - Get the current timestep for the HDF5 output. Fields are stacked in time.

956:   Not collective

958:   Input Parameter:
959: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`

961:   Output Parameter:
962: . timestep - The timestep

964:   Level: intermediate

966:   Note:
967:   This can be called only if the viewer is in the timestepping mode. See `PetscViewerHDF5PushTimestepping()` for details.

969: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5SetTimestep()`
970: @*/
971: PetscErrorCode PetscViewerHDF5GetTimestep(PetscViewer viewer, PetscInt *timestep)
972: {
973:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;

978:   *timestep = hdf5->timestep;
979:   return 0;
980: }

982: /*@C
983:   PetscDataTypeToHDF5DataType - Converts the PETSc name of a datatype to its HDF5 name.

985:   Not collective

987:   Input Parameter:
988: . ptype - the PETSc datatype name (for example `PETSC_DOUBLE`)

990:   Output Parameter:
991: . mtype - the HDF5  datatype

993:   Level: advanced

995: .seealso: [](sec_viewers), `PetscDataType`, `PetscHDF5DataTypeToPetscDataType()`
996: @*/
997: PetscErrorCode PetscDataTypeToHDF5DataType(PetscDataType ptype, hid_t *htype)
998: {
999:   if (ptype == PETSC_INT)
1000: #if defined(PETSC_USE_64BIT_INDICES)
1001:     *htype = H5T_NATIVE_LLONG;
1002: #else
1003:     *htype = H5T_NATIVE_INT;
1004: #endif
1005:   else if (ptype == PETSC_DOUBLE) *htype = H5T_NATIVE_DOUBLE;
1006:   else if (ptype == PETSC_LONG) *htype = H5T_NATIVE_LONG;
1007:   else if (ptype == PETSC_SHORT) *htype = H5T_NATIVE_SHORT;
1008:   else if (ptype == PETSC_ENUM) *htype = H5T_NATIVE_INT;
1009:   else if (ptype == PETSC_BOOL) *htype = H5T_NATIVE_INT;
1010:   else if (ptype == PETSC_FLOAT) *htype = H5T_NATIVE_FLOAT;
1011:   else if (ptype == PETSC_CHAR) *htype = H5T_NATIVE_CHAR;
1012:   else if (ptype == PETSC_BIT_LOGICAL) *htype = H5T_NATIVE_UCHAR;
1013:   else if (ptype == PETSC_STRING) *htype = H5Tcopy(H5T_C_S1);
1014:   else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported PETSc datatype");
1015:   return 0;
1016: }

1018: /*@C
1019:   PetscHDF5DataTypeToPetscDataType - Finds the PETSc name of a datatype from its HDF5 name

1021:   Not collective

1023:   Input Parameter:
1024: . htype - the HDF5 datatype (for example `H5T_NATIVE_DOUBLE`, ...)

1026:   Output Parameter:
1027: . ptype - the PETSc datatype name (for example `PETSC_DOUBLE`)

1029:   Level: advanced

1031: .seealso: [](sec_viewers), `PetscDataType`, `PetscHDF5DataTypeToPetscDataType()`
1032: @*/
1033: PetscErrorCode PetscHDF5DataTypeToPetscDataType(hid_t htype, PetscDataType *ptype)
1034: {
1035: #if defined(PETSC_USE_64BIT_INDICES)
1036:   if (htype == H5T_NATIVE_INT) *ptype = PETSC_LONG;
1037:   else if (htype == H5T_NATIVE_LLONG) *ptype = PETSC_INT;
1038: #else
1039:   if (htype == H5T_NATIVE_INT) *ptype = PETSC_INT;
1040: #endif
1041:   else if (htype == H5T_NATIVE_DOUBLE) *ptype = PETSC_DOUBLE;
1042:   else if (htype == H5T_NATIVE_LONG) *ptype = PETSC_LONG;
1043:   else if (htype == H5T_NATIVE_SHORT) *ptype = PETSC_SHORT;
1044:   else if (htype == H5T_NATIVE_FLOAT) *ptype = PETSC_FLOAT;
1045:   else if (htype == H5T_NATIVE_CHAR) *ptype = PETSC_CHAR;
1046:   else if (htype == H5T_NATIVE_UCHAR) *ptype = PETSC_CHAR;
1047:   else if (htype == H5T_C_S1) *ptype = PETSC_STRING;
1048:   else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported HDF5 datatype");
1049:   return 0;
1050: }

1052: /*@C
1053:  PetscViewerHDF5WriteAttribute - Write an attribute

1055:   Collective

1057:   Input Parameters:
1058: + viewer - The `PETSCVIEWERHDF5` viewer
1059: . parent - The parent dataset/group name
1060: . name   - The attribute name
1061: . datatype - The attribute type
1062: - value    - The attribute value

1064:   Level: advanced

1066:   Note:
1067:   If parent starts with '/', it is taken as an absolute path overriding currently pushed group, else parent is relative to the current pushed group. NULL means the current pushed group.

1069: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5WriteObjectAttribute()`, `PetscViewerHDF5ReadAttribute()`, `PetscViewerHDF5HasAttribute()`,
1070:           `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1071: @*/
1072: PetscErrorCode PetscViewerHDF5WriteAttribute(PetscViewer viewer, const char parent[], const char name[], PetscDataType datatype, const void *value)
1073: {
1074:   char     *parentAbsPath;
1075:   hid_t     h5, dataspace, obj, attribute, dtype;
1076:   PetscBool has;

1083:   PetscViewerHDF5GetGroup(viewer, parent, &parentAbsPath);
1084:   PetscViewerHDF5Traverse_Internal(viewer, parentAbsPath, PETSC_TRUE, NULL, NULL);
1085:   PetscViewerHDF5HasAttribute_Internal(viewer, parentAbsPath, name, &has);
1086:   PetscDataTypeToHDF5DataType(datatype, &dtype);
1087:   if (datatype == PETSC_STRING) {
1088:     size_t len;
1089:     PetscStrlen((const char *)value, &len);
1090:     PetscCallHDF5(H5Tset_size, (dtype, len + 1));
1091:   }
1092:   PetscViewerHDF5GetFileId(viewer, &h5);
1093:   PetscCallHDF5Return(dataspace, H5Screate, (H5S_SCALAR));
1094:   PetscCallHDF5Return(obj, H5Oopen, (h5, parentAbsPath, H5P_DEFAULT));
1095:   if (has) {
1096:     PetscCallHDF5Return(attribute, H5Aopen_name, (obj, name));
1097:   } else {
1098:     PetscCallHDF5Return(attribute, H5Acreate2, (obj, name, dtype, dataspace, H5P_DEFAULT, H5P_DEFAULT));
1099:   }
1100:   PetscCallHDF5(H5Awrite, (attribute, dtype, value));
1101:   if (datatype == PETSC_STRING) PetscCallHDF5(H5Tclose, (dtype));
1102:   PetscCallHDF5(H5Aclose, (attribute));
1103:   PetscCallHDF5(H5Oclose, (obj));
1104:   PetscCallHDF5(H5Sclose, (dataspace));
1105:   PetscFree(parentAbsPath);
1106:   return 0;
1107: }

1109: /*@C
1110:  PetscViewerHDF5WriteObjectAttribute - Write an attribute to the dataset matching the given `PetscObject` by name

1112:   Collective

1114:   Input Parameters:
1115: + viewer   - The `PETSCVIEWERHDF5` viewer
1116: . obj      - The object whose name is used to lookup the parent dataset, relative to the current group.
1117: . name     - The attribute name
1118: . datatype - The attribute type
1119: - value    - The attribute value

1121:   Note:
1122:   This fails if the path current_group/object_name doesn't resolve to a dataset (the path doesn't exist or is not a dataset).
1123:   You might want to check first if it does using `PetscViewerHDF5HasObject()`.

1125:   Level: advanced

1127: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5WriteAttribute()`, `PetscViewerHDF5ReadObjectAttribute()`, `PetscViewerHDF5HasObjectAttribute()`,
1128:           `PetscViewerHDF5HasObject()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1129: @*/
1130: PetscErrorCode PetscViewerHDF5WriteObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscDataType datatype, const void *value)
1131: {
1136:   PetscViewerHDF5CheckNamedObject_Internal(viewer, obj);
1137:   PetscViewerHDF5WriteAttribute(viewer, obj->name, name, datatype, value);
1138:   return 0;
1139: }

1141: /*@C
1142:  PetscViewerHDF5ReadAttribute - Read an attribute

1144:   Collective

1146:   Input Parameters:
1147: + viewer - The `PETSCVIEWERHDF5` viewer
1148: . parent - The parent dataset/group name
1149: . name   - The attribute name
1150: . datatype - The attribute type
1151: - defaultValue - The pointer to the default value

1153:   Output Parameter:
1154: . value    - The pointer to the read HDF5 attribute value

1156:   Notes:
1157:   If defaultValue is NULL and the attribute is not found, an error occurs.
1158:   If defaultValue is not NULL and the attribute is not found, defaultValue is copied to value.
1159:   The pointers defaultValue and value can be the same; for instance
1160: $  flg = `PETSC_FALSE`;
1161: $  `PetscViewerHDF5ReadAttribute`(viewer,name,"attr",PETSC_BOOL,&flg,&flg);
1162:   is valid, but make sure the default value is initialized.

1164:   If the datatype is `PETSC_STRING`, the output string is newly allocated so one must `PetscFree()` it when no longer needed.

1166:   If parent starts with '/', it is taken as an absolute path overriding currently pushed group, else parent is relative to the current pushed group. NULL means the current pushed group.

1168:   Level: advanced

1170: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5ReadObjectAttribute()`, `PetscViewerHDF5WriteAttribute()`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5HasObject()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1171: @*/
1172: PetscErrorCode PetscViewerHDF5ReadAttribute(PetscViewer viewer, const char parent[], const char name[], PetscDataType datatype, const void *defaultValue, void *value)
1173: {
1174:   char     *parentAbsPath;
1175:   hid_t     h5, obj, attribute, dtype;
1176:   PetscBool has;

1183:   PetscDataTypeToHDF5DataType(datatype, &dtype);
1184:   PetscViewerHDF5GetGroup(viewer, parent, &parentAbsPath);
1185:   PetscViewerHDF5Traverse_Internal(viewer, parentAbsPath, PETSC_FALSE, &has, NULL);
1186:   if (has) PetscViewerHDF5HasAttribute_Internal(viewer, parentAbsPath, name, &has);
1187:   if (!has) {
1188:     if (defaultValue) {
1189:       if (defaultValue != value) {
1190:         if (datatype == PETSC_STRING) {
1191:           PetscStrallocpy(*(char **)defaultValue, (char **)value);
1192:         } else {
1193:           size_t len;
1194:           PetscCallHDF5ReturnNoCheck(len, H5Tget_size, (dtype));
1195:           PetscMemcpy(value, defaultValue, len);
1196:         }
1197:       }
1198:       PetscFree(parentAbsPath);
1199:       return 0;
1200:     } else SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Attribute %s/%s does not exist and default value not provided", parentAbsPath, name);
1201:   }
1202:   PetscViewerHDF5GetFileId(viewer, &h5);
1203:   PetscCallHDF5Return(obj, H5Oopen, (h5, parentAbsPath, H5P_DEFAULT));
1204:   PetscCallHDF5Return(attribute, H5Aopen_name, (obj, name));
1205:   if (datatype == PETSC_STRING) {
1206:     size_t len;
1207:     hid_t  atype;
1208:     PetscCallHDF5Return(atype, H5Aget_type, (attribute));
1209:     PetscCallHDF5ReturnNoCheck(len, H5Tget_size, (atype));
1210:     PetscMalloc((len + 1) * sizeof(char), value);
1211:     PetscCallHDF5(H5Tset_size, (dtype, len + 1));
1212:     PetscCallHDF5(H5Aread, (attribute, dtype, *(char **)value));
1213:   } else {
1214:     PetscCallHDF5(H5Aread, (attribute, dtype, value));
1215:   }
1216:   PetscCallHDF5(H5Aclose, (attribute));
1217:   /* H5Oclose can be used to close groups, datasets, or committed datatypes */
1218:   PetscCallHDF5(H5Oclose, (obj));
1219:   PetscFree(parentAbsPath);
1220:   return 0;
1221: }

1223: /*@C
1224:  PetscViewerHDF5ReadObjectAttribute - Read an attribute from the dataset matching the given `PetscObject` by name

1226:   Collective

1228:   Input Parameters:
1229: + viewer   - The `PETSCVIEWERHDF5` viewer
1230: . obj      - The object whose name is used to lookup the parent dataset, relative to the current group.
1231: . name     - The attribute name
1232: - datatype - The attribute type

1234:   Output Parameter:
1235: . value    - The attribute value

1237:   Note:
1238:   This fails if current_group/object_name doesn't resolve to a dataset (the path doesn't exist or is not a dataset).
1239:   You might want to check first if it does using `PetscViewerHDF5HasObject()`.

1241:   Level: advanced

1243: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5ReadAttribute()` `PetscViewerHDF5WriteObjectAttribute()`, `PetscViewerHDF5HasObjectAttribute()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1244: @*/
1245: PetscErrorCode PetscViewerHDF5ReadObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscDataType datatype, void *defaultValue, void *value)
1246: {
1251:   PetscViewerHDF5CheckNamedObject_Internal(viewer, obj);
1252:   PetscViewerHDF5ReadAttribute(viewer, obj->name, name, datatype, defaultValue, value);
1253:   return 0;
1254: }

1256: static inline PetscErrorCode PetscViewerHDF5Traverse_Inner_Internal(hid_t h5, const char name[], PetscBool createGroup, PetscBool *exists_)
1257: {
1258:   htri_t exists;
1259:   hid_t  group;

1261:   PetscCallHDF5Return(exists, H5Lexists, (h5, name, H5P_DEFAULT));
1262:   if (exists) PetscCallHDF5Return(exists, H5Oexists_by_name, (h5, name, H5P_DEFAULT));
1263:   if (!exists && createGroup) {
1264:     PetscCallHDF5Return(group, H5Gcreate2, (h5, name, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT));
1265:     PetscCallHDF5(H5Gclose, (group));
1266:     exists = PETSC_TRUE;
1267:   }
1268:   *exists_ = (PetscBool)exists;
1269:   return 0;
1270: }

1272: static PetscErrorCode PetscViewerHDF5Traverse_Internal(PetscViewer viewer, const char name[], PetscBool createGroup, PetscBool *has, H5O_type_t *otype)
1273: {
1274:   const char rootGroupName[] = "/";
1275:   hid_t      h5;
1276:   PetscBool  exists = PETSC_FALSE;
1277:   PetscInt   i;
1278:   int        n;
1279:   char     **hierarchy;
1280:   char       buf[PETSC_MAX_PATH_LEN] = "";

1284:   else name = rootGroupName;
1285:   if (has) {
1287:     *has = PETSC_FALSE;
1288:   }
1289:   if (otype) {
1291:     *otype = H5O_TYPE_UNKNOWN;
1292:   }
1293:   PetscViewerHDF5GetFileId(viewer, &h5);

1295:   /*
1296:      Unfortunately, H5Oexists_by_name() fails if any object in hierarchy is missing.
1297:      Hence, each of them needs to be tested separately:
1298:      1) whether it's a valid link
1299:      2) whether this link resolves to an object
1300:      See H5Oexists_by_name() documentation.
1301:   */
1302:   PetscStrToArray(name, '/', &n, &hierarchy);
1303:   if (!n) {
1304:     /*  Assume group "/" always exists in accordance with HDF5 >= 1.10.0. See H5Lexists() documentation. */
1305:     if (has) *has = PETSC_TRUE;
1306:     if (otype) *otype = H5O_TYPE_GROUP;
1307:     PetscStrToArrayDestroy(n, hierarchy);
1308:     return 0;
1309:   }
1310:   for (i = 0; i < n; i++) {
1311:     PetscStrcat(buf, "/");
1312:     PetscStrcat(buf, hierarchy[i]);
1313:     PetscViewerHDF5Traverse_Inner_Internal(h5, buf, createGroup, &exists);
1314:     if (!exists) break;
1315:   }
1316:   PetscStrToArrayDestroy(n, hierarchy);

1318:   /* If the object exists, get its type */
1319:   if (exists && otype) {
1320:     H5O_info_t info;

1322:     /* We could use H5Iget_type() here but that would require opening the object. This way we only need its name. */
1323:     PetscCallHDF5(H5Oget_info_by_name, (h5, name, &info, H5P_DEFAULT));
1324:     *otype = info.type;
1325:   }
1326:   if (has) *has = exists;
1327:   return 0;
1328: }

1330: /*@C
1331:  PetscViewerHDF5HasGroup - Check whether the current (pushed) group exists in the HDF5 file

1333:   Collective

1335:   Input Parameters:
1336: + viewer - The `PETSCVIEWERHDF5` viewer
1337: - path - (Optional) The path relative to the pushed group

1339:   Output Parameter:
1340: . has - Flag for group existence

1342:   Level: advanced

1344:   Notes:
1345:   If path starts with '/', it is taken as an absolute path overriding currently pushed group, else path is relative to the current pushed group.
1346:   So NULL or empty path means the current pushed group.

1348:   If path exists but is not a group, `PETSC_FALSE` is returned.

1350: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5HasDataset()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5OpenGroup()`
1351: @*/
1352: PetscErrorCode PetscViewerHDF5HasGroup(PetscViewer viewer, const char path[], PetscBool *has)
1353: {
1354:   H5O_type_t type;
1355:   char      *abspath;

1360:   PetscViewerHDF5GetGroup(viewer, path, &abspath);
1361:   PetscViewerHDF5Traverse_Internal(viewer, abspath, PETSC_FALSE, NULL, &type);
1362:   *has = (PetscBool)(type == H5O_TYPE_GROUP);
1363:   PetscFree(abspath);
1364:   return 0;
1365: }

1367: /*@C
1368:  PetscViewerHDF5HasDataset - Check whether a given dataset exists in the HDF5 file

1370:   Collective

1372:   Input Parameters:
1373: + viewer - The `PETSCVIEWERHDF5` viewer
1374: - path - The dataset path

1376:   Output Parameter:
1377: . has - Flag whether dataset exists

1379:   Level: advanced

1381:   Notes:
1382:   If path starts with '/', it is taken as an absolute path overriding currently pushed group, else path is relative to the current pushed group.

1384:   If path is NULL or empty, has is set to `PETSC_FALSE`.

1386:   If path exists but is not a dataset, has is set to `PETSC_FALSE` as well.

1388: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5HasObject()`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5HasGroup()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1389: @*/
1390: PetscErrorCode PetscViewerHDF5HasDataset(PetscViewer viewer, const char path[], PetscBool *has)
1391: {
1392:   H5O_type_t type;
1393:   char      *abspath;

1398:   PetscViewerHDF5GetGroup(viewer, path, &abspath);
1399:   PetscViewerHDF5Traverse_Internal(viewer, abspath, PETSC_FALSE, NULL, &type);
1400:   *has = (PetscBool)(type == H5O_TYPE_DATASET);
1401:   PetscFree(abspath);
1402:   return 0;
1403: }

1405: /*@
1406:  PetscViewerHDF5HasObject - Check whether a dataset with the same name as given object exists in the HDF5 file under current group

1408:   Collective

1410:   Input Parameters:
1411: + viewer - The `PETSCVIEWERHDF5` viewer
1412: - obj    - The named object

1414:   Output Parameter:
1415: . has    - Flag for dataset existence

1417:   Notes:
1418:   If the object is unnamed, an error occurs.

1420:   If the path current_group/object_name exists but is not a dataset, has is set to `PETSC_FALSE` as well.

1422:   Level: advanced

1424: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5HasDataset()`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1425: @*/
1426: PetscErrorCode PetscViewerHDF5HasObject(PetscViewer viewer, PetscObject obj, PetscBool *has)
1427: {
1428:   size_t len;

1433:   PetscStrlen(obj->name, &len);
1435:   PetscViewerHDF5HasDataset(viewer, obj->name, has);
1436:   return 0;
1437: }

1439: /*@C
1440:  PetscViewerHDF5HasAttribute - Check whether an attribute exists

1442:   Collective

1444:   Input Parameters:
1445: + viewer - The `PETSCVIEWERHDF5` viewer
1446: . parent - The parent dataset/group name
1447: - name   - The attribute name

1449:   Output Parameter:
1450: . has    - Flag for attribute existence

1452:   Level: advanced

1454:   Note:
1455:   If parent starts with '/', it is taken as an absolute path overriding currently pushed group, else parent is relative to the current pushed group. NULL means the current pushed group.

1457: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5HasObjectAttribute()`, `PetscViewerHDF5WriteAttribute()`, `PetscViewerHDF5ReadAttribute()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1458: @*/
1459: PetscErrorCode PetscViewerHDF5HasAttribute(PetscViewer viewer, const char parent[], const char name[], PetscBool *has)
1460: {
1461:   char *parentAbsPath;

1467:   PetscViewerHDF5GetGroup(viewer, parent, &parentAbsPath);
1468:   PetscViewerHDF5Traverse_Internal(viewer, parentAbsPath, PETSC_FALSE, has, NULL);
1469:   if (*has) PetscViewerHDF5HasAttribute_Internal(viewer, parentAbsPath, name, has);
1470:   PetscFree(parentAbsPath);
1471:   return 0;
1472: }

1474: /*@C
1475:  PetscViewerHDF5HasObjectAttribute - Check whether an attribute is attached to the dataset matching the given `PetscObject` by name

1477:   Collective

1479:   Input Parameters:
1480: + viewer - The `PETSCVIEWERHDF5` viewer
1481: . obj    - The object whose name is used to lookup the parent dataset, relative to the current group.
1482: - name   - The attribute name

1484:   Output Parameter:
1485: . has    - Flag for attribute existence

1487:   Note:
1488:   This fails if current_group/object_name doesn't resolve to a dataset (the path doesn't exist or is not a dataset).
1489:   You might want to check first if it does using `PetscViewerHDF5HasObject()`.

1491:   Level: advanced

1493: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5WriteObjectAttribute()`, `PetscViewerHDF5ReadObjectAttribute()`, `PetscViewerHDF5HasObject()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1494: @*/
1495: PetscErrorCode PetscViewerHDF5HasObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscBool *has)
1496: {
1501:   PetscViewerHDF5CheckNamedObject_Internal(viewer, obj);
1502:   PetscViewerHDF5HasAttribute(viewer, obj->name, name, has);
1503:   return 0;
1504: }

1506: static PetscErrorCode PetscViewerHDF5HasAttribute_Internal(PetscViewer viewer, const char parent[], const char name[], PetscBool *has)
1507: {
1508:   hid_t  h5;
1509:   htri_t hhas;

1511:   PetscViewerHDF5GetFileId(viewer, &h5);
1512:   PetscCallHDF5Return(hhas, H5Aexists_by_name, (h5, parent, name, H5P_DEFAULT));
1513:   *has = hhas ? PETSC_TRUE : PETSC_FALSE;
1514:   return 0;
1515: }

1517: /*
1518:   The variable Petsc_Viewer_HDF5_keyval is used to indicate an MPI attribute that
1519:   is attached to a communicator, in this case the attribute is a PetscViewer.
1520: */
1521: PetscMPIInt Petsc_Viewer_HDF5_keyval = MPI_KEYVAL_INVALID;

1523: /*@C
1524:   PETSC_VIEWER_HDF5_ - Creates an `PETSCVIEWERHDF5` `PetscViewer` shared by all processors in a communicator.

1526:   Collective

1528:   Input Parameter:
1529: . comm - the MPI communicator to share the HDF5 `PetscViewer`

1531:   Level: intermediate

1533:   Options Database Key:
1534: . -viewer_hdf5_filename <name> - name of the HDF5 file

1536:   Environmental variable:
1537: . `PETSC_VIEWER_HDF5_FILENAME` - name of the HDF5 file

1539:   Note:
1540:   Unlike almost all other PETSc routines, `PETSC_VIEWER_HDF5_()` does not return
1541:   an error code.  The HDF5 `PetscViewer` is usually used in the form
1542: $       XXXView(XXX object, PETSC_VIEWER_HDF5_(comm));

1544: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerCreate()`, `PetscViewerDestroy()`
1545: @*/
1546: PetscViewer PETSC_VIEWER_HDF5_(MPI_Comm comm)
1547: {
1549:   PetscBool      flg;
1550:   PetscViewer    viewer;
1551:   char           fname[PETSC_MAX_PATH_LEN];
1552:   MPI_Comm       ncomm;

1554:   PetscCommDuplicate(comm, &ncomm, NULL);
1555:   if (ierr) {
1556:     PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_INITIAL, " ");
1557:     return NULL;
1558:   }
1559:   if (Petsc_Viewer_HDF5_keyval == MPI_KEYVAL_INVALID) {
1560:     MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN, MPI_COMM_NULL_DELETE_FN, &Petsc_Viewer_HDF5_keyval, NULL);
1561:     if (ierr) {
1562:       PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_INITIAL, " ");
1563:       return NULL;
1564:     }
1565:   }
1566:   MPI_Comm_get_attr(ncomm, Petsc_Viewer_HDF5_keyval, (void **)&viewer, (int *)&flg);
1567:   if (ierr) {
1568:     PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_INITIAL, " ");
1569:     return NULL;
1570:   }
1571:   if (!flg) { /* PetscViewer not yet created */
1572:     PetscOptionsGetenv(ncomm, "PETSC_VIEWER_HDF5_FILENAME", fname, PETSC_MAX_PATH_LEN, &flg);
1573:     if (ierr) {
1574:       PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1575:       return NULL;
1576:     }
1577:     if (!flg) {
1578:       PetscStrcpy(fname, "output.h5");
1579:       if (ierr) {
1580:         PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1581:         return NULL;
1582:       }
1583:     }
1584:     PetscViewerHDF5Open(ncomm, fname, FILE_MODE_WRITE, &viewer);
1585:     if (ierr) {
1586:       PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1587:       return NULL;
1588:     }
1589:     PetscObjectRegisterDestroy((PetscObject)viewer);
1590:     if (ierr) {
1591:       PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1592:       return NULL;
1593:     }
1594:     MPI_Comm_set_attr(ncomm, Petsc_Viewer_HDF5_keyval, (void *)viewer);
1595:     if (ierr) {
1596:       PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_INITIAL, " ");
1597:       return NULL;
1598:     }
1599:   }
1600:   PetscCommDestroy(&ncomm);
1601:   if (ierr) {
1602:     PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1603:     return NULL;
1604:   }
1605:   return viewer;
1606: }