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: `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] = "";
35: PetscFunctionBegin;
39: PetscCall(PetscViewerHDF5GetGroup_Internal(viewer, &group));
40: PetscCall(PetscStrlen(path, &len));
41: relative = (PetscBool)(!len || path[0] != '/');
42: if (relative) {
43: PetscCall(PetscStrncpy(buf, group, sizeof(buf)));
44: if (!group || len) PetscCall(PetscStrlcat(buf, "/", sizeof(buf)));
45: PetscCall(PetscStrlcat(buf, path, sizeof(buf)));
46: PetscCall(PetscStrallocpy(buf, abspath));
47: } else {
48: PetscCall(PetscStrallocpy(path, abspath));
49: }
50: PetscFunctionReturn(PETSC_SUCCESS);
51: }
53: static PetscErrorCode PetscViewerHDF5CheckNamedObject_Internal(PetscViewer viewer, PetscObject obj)
54: {
55: PetscBool has;
57: PetscFunctionBegin;
58: PetscCall(PetscViewerHDF5HasObject(viewer, obj, &has));
59: if (!has) {
60: char *group;
61: PetscCall(PetscViewerHDF5GetGroup(viewer, NULL, &group));
62: SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Object (dataset) \"%s\" not stored in group %s", obj->name, group);
63: }
64: PetscFunctionReturn(PETSC_SUCCESS);
65: }
67: static PetscErrorCode PetscViewerSetFromOptions_HDF5(PetscViewer v, PetscOptionItems *PetscOptionsObject)
68: {
69: PetscBool flg = PETSC_FALSE, set;
70: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)v->data;
72: PetscFunctionBegin;
73: PetscOptionsHeadBegin(PetscOptionsObject, "HDF5 PetscViewer Options");
74: PetscCall(PetscOptionsBool("-viewer_hdf5_base_dimension2", "1d Vectors get 2 dimensions in HDF5", "PetscViewerHDF5SetBaseDimension2", hdf5->basedimension2, &hdf5->basedimension2, NULL));
75: PetscCall(PetscOptionsBool("-viewer_hdf5_sp_output", "Force data to be written in single precision", "PetscViewerHDF5SetSPOutput", hdf5->spoutput, &hdf5->spoutput, NULL));
76: PetscCall(PetscOptionsBool("-viewer_hdf5_collective", "Enable collective transfer mode", "PetscViewerHDF5SetCollective", flg, &flg, &set));
77: if (set) PetscCall(PetscViewerHDF5SetCollective(v, flg));
78: flg = PETSC_FALSE;
79: PetscCall(PetscOptionsBool("-viewer_hdf5_default_timestepping", "Set default timestepping state", "PetscViewerHDF5SetDefaultTimestepping", flg, &flg, &set));
80: if (set) PetscCall(PetscViewerHDF5SetDefaultTimestepping(v, flg));
81: PetscOptionsHeadEnd();
82: PetscFunctionReturn(PETSC_SUCCESS);
83: }
85: static PetscErrorCode PetscViewerView_HDF5(PetscViewer v, PetscViewer viewer)
86: {
87: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)v->data;
88: PetscBool flg;
90: PetscFunctionBegin;
91: if (hdf5->filename) PetscCall(PetscViewerASCIIPrintf(viewer, "Filename: %s\n", hdf5->filename));
92: PetscCall(PetscViewerASCIIPrintf(viewer, "Vectors with blocksize 1 saved as 2D datasets: %s\n", PetscBools[hdf5->basedimension2]));
93: PetscCall(PetscViewerASCIIPrintf(viewer, "Enforce single precision storage: %s\n", PetscBools[hdf5->spoutput]));
94: PetscCall(PetscViewerHDF5GetCollective(v, &flg));
95: PetscCall(PetscViewerASCIIPrintf(viewer, "MPI-IO transfer mode: %s\n", flg ? "collective" : "independent"));
96: PetscCall(PetscViewerASCIIPrintf(viewer, "Default timestepping: %s\n", PetscBools[hdf5->defTimestepping]));
97: PetscFunctionReturn(PETSC_SUCCESS);
98: }
100: static PetscErrorCode PetscViewerFileClose_HDF5(PetscViewer viewer)
101: {
102: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
104: PetscFunctionBegin;
105: PetscCall(PetscFree(hdf5->filename));
106: if (hdf5->file_id) PetscCallHDF5(H5Fclose, (hdf5->file_id));
107: PetscFunctionReturn(PETSC_SUCCESS);
108: }
110: static PetscErrorCode PetscViewerFlush_HDF5(PetscViewer viewer)
111: {
112: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
114: PetscFunctionBegin;
115: if (hdf5->file_id) PetscCallHDF5(H5Fflush, (hdf5->file_id, H5F_SCOPE_LOCAL));
116: PetscFunctionReturn(PETSC_SUCCESS);
117: }
119: static PetscErrorCode PetscViewerDestroy_HDF5(PetscViewer viewer)
120: {
121: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
123: PetscFunctionBegin;
124: PetscCallHDF5(H5Pclose, (hdf5->dxpl_id));
125: PetscCall(PetscViewerFileClose_HDF5(viewer));
126: while (hdf5->groups) {
127: PetscViewerHDF5GroupList *tmp = hdf5->groups->next;
129: PetscCall(PetscFree(hdf5->groups->name));
130: PetscCall(PetscFree(hdf5->groups));
131: hdf5->groups = tmp;
132: }
133: PetscCall(PetscFree(hdf5));
134: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileSetName_C", NULL));
135: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileGetName_C", NULL));
136: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileSetMode_C", NULL));
137: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileGetMode_C", NULL));
138: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5SetBaseDimension2_C", NULL));
139: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5SetSPOutput_C", NULL));
140: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5SetCollective_C", NULL));
141: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5GetCollective_C", NULL));
142: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5GetDefaultTimestepping_C", NULL));
143: PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerHDF5SetDefaultTimestepping_C", NULL));
144: PetscFunctionReturn(PETSC_SUCCESS);
145: }
147: static PetscErrorCode PetscViewerFileSetMode_HDF5(PetscViewer viewer, PetscFileMode type)
148: {
149: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
151: PetscFunctionBegin;
152: hdf5->btype = type;
153: PetscFunctionReturn(PETSC_SUCCESS);
154: }
156: static PetscErrorCode PetscViewerFileGetMode_HDF5(PetscViewer viewer, PetscFileMode *type)
157: {
158: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
160: PetscFunctionBegin;
161: *type = hdf5->btype;
162: PetscFunctionReturn(PETSC_SUCCESS);
163: }
165: static PetscErrorCode PetscViewerHDF5SetBaseDimension2_HDF5(PetscViewer viewer, PetscBool flg)
166: {
167: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
169: PetscFunctionBegin;
170: hdf5->basedimension2 = flg;
171: PetscFunctionReturn(PETSC_SUCCESS);
172: }
174: /*@
175: PetscViewerHDF5SetBaseDimension2 - Vectors of 1 dimension (i.e. bs/dof is 1) will be saved in the HDF5 file with a
176: dimension of 2.
178: Logically Collective
180: Input Parameters:
181: + viewer - the `PetscViewer`; if it is a `PETSCVIEWERHDF5` then this command is ignored
182: - flg - if `PETSC_TRUE` the vector will always have at least a dimension of 2 even if that first dimension is of size 1
184: Options Database Key:
185: . -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
187: Level: intermediate
189: Note:
190: Setting this option allegedly makes code that reads the HDF5 in easier since they do not have a "special case" of a bs/dof
191: of one when the dimension is lower. Others think the option is crazy.
193: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, PetscViewerFileSetMode()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerBinaryOpen()`
194: @*/
195: PetscErrorCode PetscViewerHDF5SetBaseDimension2(PetscViewer viewer, PetscBool flg)
196: {
197: PetscFunctionBegin;
199: PetscTryMethod(viewer, "PetscViewerHDF5SetBaseDimension2_C", (PetscViewer, PetscBool), (viewer, flg));
200: PetscFunctionReturn(PETSC_SUCCESS);
201: }
203: /*@
204: PetscViewerHDF5GetBaseDimension2 - Vectors of 1 dimension (i.e. bs/dof is 1) will be saved in the HDF5 file with a
205: dimension of 2.
207: Logically Collective
209: Input Parameter:
210: . viewer - the `PetscViewer`, must be `PETSCVIEWERHDF5`
212: Output Parameter:
213: . flg - if `PETSC_TRUE` the vector will always have at least a dimension of 2 even if that first dimension is of size 1
215: Level: intermediate
217: Note:
218: Setting this option allegedly makes code that reads the HDF5 in easier since they do not have a "special case" of a bs/dof
219: of one when the dimension is lower. Others think the option is crazy.
221: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerFileSetMode()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerBinaryOpen()`
222: @*/
223: PetscErrorCode PetscViewerHDF5GetBaseDimension2(PetscViewer viewer, PetscBool *flg)
224: {
225: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
227: PetscFunctionBegin;
229: *flg = hdf5->basedimension2;
230: PetscFunctionReturn(PETSC_SUCCESS);
231: }
233: static PetscErrorCode PetscViewerHDF5SetSPOutput_HDF5(PetscViewer viewer, PetscBool flg)
234: {
235: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
237: PetscFunctionBegin;
238: hdf5->spoutput = flg;
239: PetscFunctionReturn(PETSC_SUCCESS);
240: }
242: /*@
243: PetscViewerHDF5SetSPOutput - Data is written to disk in single precision even if PETSc is
244: compiled with double precision `PetscReal`.
246: Logically Collective
248: Input Parameters:
249: + viewer - the PetscViewer; if it is a `PETSCVIEWERHDF5` then this command is ignored
250: - flg - if `PETSC_TRUE` the data will be written to disk with single precision
252: Options Database Key:
253: . -viewer_hdf5_sp_output - turns on (true) or off (false) output in single precision
255: Level: intermediate
257: Note:
258: Setting this option does not make any difference if PETSc is compiled with single precision
259: in the first place. It does not affect reading datasets (HDF5 handle this internally).
261: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerFileSetMode()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerBinaryOpen()`,
262: `PetscReal`, `PetscViewerHDF5GetSPOutput()`
263: @*/
264: PetscErrorCode PetscViewerHDF5SetSPOutput(PetscViewer viewer, PetscBool flg)
265: {
266: PetscFunctionBegin;
268: PetscTryMethod(viewer, "PetscViewerHDF5SetSPOutput_C", (PetscViewer, PetscBool), (viewer, flg));
269: PetscFunctionReturn(PETSC_SUCCESS);
270: }
272: /*@
273: PetscViewerHDF5GetSPOutput - Data is written to disk in single precision even if PETSc is
274: compiled with double precision `PetscReal`.
276: Logically Collective
278: Input Parameter:
279: . viewer - the PetscViewer, must be of type `PETSCVIEWERHDF5`
281: Output Parameter:
282: . flg - if `PETSC_TRUE` the data will be written to disk with single precision
284: Level: intermediate
286: .seealso: [](sec_viewers), `PetscViewerFileSetMode()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerBinaryOpen()`,
287: `PetscReal`, `PetscViewerHDF5SetSPOutput()`
288: @*/
289: PetscErrorCode PetscViewerHDF5GetSPOutput(PetscViewer viewer, PetscBool *flg)
290: {
291: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
293: PetscFunctionBegin;
295: *flg = hdf5->spoutput;
296: PetscFunctionReturn(PETSC_SUCCESS);
297: }
299: static PetscErrorCode PetscViewerHDF5SetCollective_HDF5(PetscViewer viewer, PetscBool flg)
300: {
301: PetscFunctionBegin;
302: /* H5FD_MPIO_COLLECTIVE is wrong in hdf5 1.10.2, and is the same as H5FD_MPIO_INDEPENDENT in earlier versions
303: - see e.g. https://gitlab.cosma.dur.ac.uk/swift/swiftsim/issues/431 */
304: #if H5_VERSION_GE(1, 10, 3) && defined(H5_HAVE_PARALLEL)
305: {
306: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
307: PetscCallHDF5(H5Pset_dxpl_mpio, (hdf5->dxpl_id, flg ? H5FD_MPIO_COLLECTIVE : H5FD_MPIO_INDEPENDENT));
308: }
309: #else
310: if (flg) PetscCall(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"));
311: #endif
312: PetscFunctionReturn(PETSC_SUCCESS);
313: }
315: /*@
316: PetscViewerHDF5SetCollective - Use collective MPI-IO transfer mode for HDF5 reads and writes.
318: Logically Collective; flg must contain common value
320: Input Parameters:
321: + viewer - the `PetscViewer`; if it is not `PETSCVIEWERHDF5` then this command is ignored
322: - flg - `PETSC_TRUE` for collective mode; `PETSC_FALSE` for independent mode (default)
324: Options Database Key:
325: . -viewer_hdf5_collective - turns on (true) or off (false) collective transfers
327: Level: intermediate
329: Note:
330: Collective mode gives the MPI-IO layer underneath HDF5 a chance to do some additional collective optimizations and hence can perform better.
331: 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.
333: Developer Note:
334: In the HDF5 layer, `PETSC_TRUE` / `PETSC_FALSE` means `H5Pset_dxpl_mpio()` is called with `H5FD_MPIO_COLLECTIVE` / `H5FD_MPIO_INDEPENDENT`, respectively.
335: This in turn means use of MPI_File_{read,write}_all / MPI_File_{read,write} in the MPI-IO layer, respectively.
336: See HDF5 documentation and MPI-IO documentation for details.
338: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5GetCollective()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerHDF5Open()`
339: @*/
340: PetscErrorCode PetscViewerHDF5SetCollective(PetscViewer viewer, PetscBool flg)
341: {
342: PetscFunctionBegin;
345: PetscTryMethod(viewer, "PetscViewerHDF5SetCollective_C", (PetscViewer, PetscBool), (viewer, flg));
346: PetscFunctionReturn(PETSC_SUCCESS);
347: }
349: static PetscErrorCode PetscViewerHDF5GetCollective_HDF5(PetscViewer viewer, PetscBool *flg)
350: {
351: #if defined(H5_HAVE_PARALLEL)
352: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
353: H5FD_mpio_xfer_t mode;
354: #endif
356: PetscFunctionBegin;
357: #if !defined(H5_HAVE_PARALLEL)
358: *flg = PETSC_FALSE;
359: #else
360: PetscCallHDF5(H5Pget_dxpl_mpio, (hdf5->dxpl_id, &mode));
361: *flg = (mode == H5FD_MPIO_COLLECTIVE) ? PETSC_TRUE : PETSC_FALSE;
362: #endif
363: PetscFunctionReturn(PETSC_SUCCESS);
364: }
366: /*@
367: PetscViewerHDF5GetCollective - Return flag whether collective MPI-IO transfer mode is used for HDF5 reads and writes.
369: Not Collective
371: Input Parameter:
372: . viewer - the `PETSCVIEWERHDF5` `PetscViewer`
374: Output Parameter:
375: . flg - the flag
377: Level: intermediate
379: Note:
380: 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.
381: For more details, see `PetscViewerHDF5SetCollective()`.
383: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5SetCollective()`, `PetscViewerCreate()`, `PetscViewerSetType()`, `PetscViewerHDF5Open()`
384: @*/
385: PetscErrorCode PetscViewerHDF5GetCollective(PetscViewer viewer, PetscBool *flg)
386: {
387: PetscFunctionBegin;
391: PetscUseMethod(viewer, "PetscViewerHDF5GetCollective_C", (PetscViewer, PetscBool *), (viewer, flg));
392: PetscFunctionReturn(PETSC_SUCCESS);
393: }
395: static PetscErrorCode PetscViewerFileSetName_HDF5(PetscViewer viewer, const char name[])
396: {
397: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
398: hid_t plist_id;
400: PetscFunctionBegin;
401: if (hdf5->file_id) PetscCallHDF5(H5Fclose, (hdf5->file_id));
402: if (hdf5->filename) PetscCall(PetscFree(hdf5->filename));
403: PetscCall(PetscStrallocpy(name, &hdf5->filename));
404: /* Set up file access property list with parallel I/O access */
405: PetscCallHDF5Return(plist_id, H5Pcreate, (H5P_FILE_ACCESS));
406: #if defined(H5_HAVE_PARALLEL)
407: PetscCallHDF5(H5Pset_fapl_mpio, (plist_id, PetscObjectComm((PetscObject)viewer), MPI_INFO_NULL));
408: #endif
409: /* Create or open the file collectively */
410: switch (hdf5->btype) {
411: case FILE_MODE_READ:
412: if (PetscDefined(USE_DEBUG)) {
413: PetscMPIInt rank;
414: PetscBool flg;
416: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank));
417: if (rank == 0) {
418: PetscCall(PetscTestFile(hdf5->filename, 'r', &flg));
419: PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "File %s requested for reading does not exist", hdf5->filename);
420: }
421: PetscCallMPI(MPI_Barrier(PetscObjectComm((PetscObject)viewer)));
422: }
423: PetscCallHDF5Return(hdf5->file_id, H5Fopen, (name, H5F_ACC_RDONLY, plist_id));
424: break;
425: case FILE_MODE_APPEND:
426: case FILE_MODE_UPDATE: {
427: PetscBool flg;
428: PetscCall(PetscTestFile(hdf5->filename, 'r', &flg));
429: if (flg) PetscCallHDF5Return(hdf5->file_id, H5Fopen, (name, H5F_ACC_RDWR, plist_id));
430: else PetscCallHDF5Return(hdf5->file_id, H5Fcreate, (name, H5F_ACC_EXCL, H5P_DEFAULT, plist_id));
431: break;
432: }
433: case FILE_MODE_WRITE:
434: PetscCallHDF5Return(hdf5->file_id, H5Fcreate, (name, H5F_ACC_TRUNC, H5P_DEFAULT, plist_id));
435: break;
436: case FILE_MODE_UNDEFINED:
437: SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
438: default:
439: SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Unsupported file mode %s", PetscFileModes[hdf5->btype]);
440: }
441: PetscCheck(hdf5->file_id >= 0, PETSC_COMM_SELF, PETSC_ERR_LIB, "H5Fcreate failed for %s", name);
442: PetscCallHDF5(H5Pclose, (plist_id));
443: PetscCall(PetscViewerHDF5ResetAttachedDMPlexStorageVersion(viewer));
444: PetscFunctionReturn(PETSC_SUCCESS);
445: }
447: static PetscErrorCode PetscViewerFileGetName_HDF5(PetscViewer viewer, const char **name)
448: {
449: PetscViewer_HDF5 *vhdf5 = (PetscViewer_HDF5 *)viewer->data;
451: PetscFunctionBegin;
452: *name = vhdf5->filename;
453: PetscFunctionReturn(PETSC_SUCCESS);
454: }
456: static PetscErrorCode PetscViewerSetUp_HDF5(PetscViewer viewer)
457: {
458: /*
459: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
460: */
462: PetscFunctionBegin;
463: PetscFunctionReturn(PETSC_SUCCESS);
464: }
466: static PetscErrorCode PetscViewerHDF5SetDefaultTimestepping_HDF5(PetscViewer viewer, PetscBool flg)
467: {
468: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
470: PetscFunctionBegin;
471: hdf5->defTimestepping = flg;
472: PetscFunctionReturn(PETSC_SUCCESS);
473: }
475: static PetscErrorCode PetscViewerHDF5GetDefaultTimestepping_HDF5(PetscViewer viewer, PetscBool *flg)
476: {
477: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
479: PetscFunctionBegin;
480: *flg = hdf5->defTimestepping;
481: PetscFunctionReturn(PETSC_SUCCESS);
482: }
484: /*@
485: PetscViewerHDF5SetDefaultTimestepping - Set the flag for default timestepping
487: Logically Collective
489: Input Parameters:
490: + viewer - the `PetscViewer`; if it is not `PETSCVIEWERHDF5` then this command is ignored
491: - flg - if `PETSC_TRUE` we will assume that timestepping is on
493: Options Database Key:
494: . -viewer_hdf5_default_timestepping - turns on (true) or off (false) default timestepping
496: Level: intermediate
498: Note:
499: If the timestepping attribute is not found for an object, then the default timestepping is used
501: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5GetDefaultTimestepping()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5GetTimestep()`
502: @*/
503: PetscErrorCode PetscViewerHDF5SetDefaultTimestepping(PetscViewer viewer, PetscBool flg)
504: {
505: PetscFunctionBegin;
507: PetscTryMethod(viewer, "PetscViewerHDF5SetDefaultTimestepping_C", (PetscViewer, PetscBool), (viewer, flg));
508: PetscFunctionReturn(PETSC_SUCCESS);
509: }
511: /*@
512: PetscViewerHDF5GetDefaultTimestepping - Get the flag for default timestepping
514: Not Collective
516: Input Parameter:
517: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
519: Output Parameter:
520: . flg - if `PETSC_TRUE` we will assume that timestepping is on
522: Level: intermediate
524: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5SetDefaultTimestepping()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5GetTimestep()`
525: @*/
526: PetscErrorCode PetscViewerHDF5GetDefaultTimestepping(PetscViewer viewer, PetscBool *flg)
527: {
528: PetscFunctionBegin;
530: PetscUseMethod(viewer, "PetscViewerHDF5GetDefaultTimestepping_C", (PetscViewer, PetscBool *), (viewer, flg));
531: PetscFunctionReturn(PETSC_SUCCESS);
532: }
534: /*MC
535: PETSCVIEWERHDF5 - A viewer that writes to an HDF5 file
537: Level: beginner
539: .seealso: [](sec_viewers), `PetscViewerHDF5Open()`, `PetscViewerStringSPrintf()`, `PetscViewerSocketOpen()`, `PetscViewerDrawOpen()`, `PETSCVIEWERSOCKET`,
540: `PetscViewerCreate()`, `PetscViewerASCIIOpen()`, `PetscViewerBinaryOpen()`, `PETSCVIEWERBINARY`, `PETSCVIEWERDRAW`, `PETSCVIEWERSTRING`,
541: `PetscViewerMatlabOpen()`, `VecView()`, `DMView()`, `PetscViewerMatlabPutArray()`, `PETSCVIEWERASCII`, `PETSCVIEWERMATLAB`,
542: `PetscViewerFileSetName()`, `PetscViewerFileSetMode()`, `PetscViewerFormat`, `PetscViewerType`, `PetscViewerSetType()`
543: M*/
545: PETSC_EXTERN PetscErrorCode PetscViewerCreate_HDF5(PetscViewer v)
546: {
547: PetscViewer_HDF5 *hdf5;
549: PetscFunctionBegin;
550: #if !defined(H5_HAVE_PARALLEL)
551: {
552: PetscMPIInt size;
553: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)v), &size));
554: PetscCheck(size <= 1, PetscObjectComm((PetscObject)v), PETSC_ERR_SUP, "Cannot use parallel HDF5 viewer since the given HDF5 does not support parallel I/O (H5_HAVE_PARALLEL is unset)");
555: }
556: #endif
558: PetscCall(PetscNew(&hdf5));
560: v->data = (void *)hdf5;
561: v->ops->destroy = PetscViewerDestroy_HDF5;
562: v->ops->setfromoptions = PetscViewerSetFromOptions_HDF5;
563: v->ops->setup = PetscViewerSetUp_HDF5;
564: v->ops->view = PetscViewerView_HDF5;
565: v->ops->flush = PetscViewerFlush_HDF5;
566: hdf5->btype = FILE_MODE_UNDEFINED;
567: hdf5->filename = NULL;
568: hdf5->timestep = -1;
569: hdf5->groups = NULL;
571: PetscCallHDF5Return(hdf5->dxpl_id, H5Pcreate, (H5P_DATASET_XFER));
573: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileSetName_C", PetscViewerFileSetName_HDF5));
574: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileGetName_C", PetscViewerFileGetName_HDF5));
575: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileSetMode_C", PetscViewerFileSetMode_HDF5));
576: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileGetMode_C", PetscViewerFileGetMode_HDF5));
577: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5SetBaseDimension2_C", PetscViewerHDF5SetBaseDimension2_HDF5));
578: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5SetSPOutput_C", PetscViewerHDF5SetSPOutput_HDF5));
579: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5SetCollective_C", PetscViewerHDF5SetCollective_HDF5));
580: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5GetCollective_C", PetscViewerHDF5GetCollective_HDF5));
581: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5GetDefaultTimestepping_C", PetscViewerHDF5GetDefaultTimestepping_HDF5));
582: PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerHDF5SetDefaultTimestepping_C", PetscViewerHDF5SetDefaultTimestepping_HDF5));
583: PetscFunctionReturn(PETSC_SUCCESS);
584: }
586: /*@C
587: PetscViewerHDF5Open - Opens a file for HDF5 input/output as a `PETSCVIEWERHDF5` `PetscViewer`
589: Collective
591: Input Parameters:
592: + comm - MPI communicator
593: . name - name of file
594: - type - type of file
596: Output Parameter:
597: . hdf5v - `PetscViewer` for HDF5 input/output to use with the specified file
599: Options Database Keys:
600: + -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
601: - -viewer_hdf5_sp_output - forces (if true) the viewer to write data in single precision independent on the precision of PetscReal
603: Level: beginner
605: Notes:
606: Reading is always available, regardless of the mode. Available modes are
607: .vb
608: FILE_MODE_READ - open existing HDF5 file for read only access, fail if file does not exist [H5Fopen() with H5F_ACC_RDONLY]
609: FILE_MODE_WRITE - if file exists, fully overwrite it, else create new HDF5 file [H5FcreateH5Fcreate() with H5F_ACC_TRUNC]
610: FILE_MODE_APPEND - if file exists, keep existing contents [H5Fopen() with H5F_ACC_RDWR], else create new HDF5 file [H5FcreateH5Fcreate() with H5F_ACC_EXCL]
611: FILE_MODE_UPDATE - same as FILE_MODE_APPEND
612: .ve
614: 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.
616: This PetscViewer should be destroyed with PetscViewerDestroy().
618: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerASCIIOpen()`, `PetscViewerPushFormat()`, `PetscViewerDestroy()`, `PetscViewerHDF5SetBaseDimension2()`,
619: `PetscViewerHDF5SetSPOutput()`, `PetscViewerHDF5GetBaseDimension2()`, `VecView()`, `MatView()`, `VecLoad()`,
620: `MatLoad()`, `PetscFileMode`, `PetscViewer`, `PetscViewerSetType()`, `PetscViewerFileSetMode()`, `PetscViewerFileSetName()`
621: @*/
622: PetscErrorCode PetscViewerHDF5Open(MPI_Comm comm, const char name[], PetscFileMode type, PetscViewer *hdf5v)
623: {
624: PetscFunctionBegin;
625: PetscCall(PetscViewerCreate(comm, hdf5v));
626: PetscCall(PetscViewerSetType(*hdf5v, PETSCVIEWERHDF5));
627: PetscCall(PetscViewerFileSetMode(*hdf5v, type));
628: PetscCall(PetscViewerFileSetName(*hdf5v, name));
629: PetscCall(PetscViewerSetFromOptions(*hdf5v));
630: PetscFunctionReturn(PETSC_SUCCESS);
631: }
633: /*@C
634: PetscViewerHDF5GetFileId - Retrieve the file id, this file ID then can be used in direct HDF5 calls
636: Not Collective
638: Input Parameter:
639: . viewer - the `PetscViewer`
641: Output Parameter:
642: . file_id - The file id
644: Level: intermediate
646: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`
647: @*/
648: PetscErrorCode PetscViewerHDF5GetFileId(PetscViewer viewer, hid_t *file_id)
649: {
650: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
652: PetscFunctionBegin;
654: if (file_id) *file_id = hdf5->file_id;
655: PetscFunctionReturn(PETSC_SUCCESS);
656: }
658: /*@C
659: PetscViewerHDF5PushGroup - Set the current HDF5 group for output
661: Not Collective
663: Input Parameters:
664: + viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
665: - name - The group name
667: Level: intermediate
669: Notes:
670: This is designed to mnemonically resemble the Unix cd command.
671: .vb
672: 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.
673: `NULL`, empty string, or any sequence of all slashes (e.g. "///") is interpreted as the root group "/".
674: "." means the current group is pushed again.
675: .ve
677: Example:
678: Suppose the current group is "/a".
679: .vb
680: If name is `NULL`, empty string, or a sequence of all slashes (e.g. "///"), then the new group will be "/".
681: If name is ".", then the new group will be "/a".
682: If name is "b", then the new group will be "/a/b".
683: If name is "/b", then the new group will be "/b".
684: .ve
686: Developer Note:
687: The root group "/" is internally stored as `NULL`.
689: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5OpenGroup()`, `PetscViewerHDF5WriteGroup()`
690: @*/
691: PetscErrorCode PetscViewerHDF5PushGroup(PetscViewer viewer, const char name[])
692: {
693: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
694: PetscViewerHDF5GroupList *groupNode;
695: size_t i, len;
696: char buf[PETSC_MAX_PATH_LEN];
697: const char *gname;
699: PetscFunctionBegin;
702: PetscCall(PetscStrlen(name, &len));
703: gname = NULL;
704: if (len) {
705: if (len == 1 && name[0] == '.') {
706: /* use current name */
707: gname = (hdf5->groups && hdf5->groups->name) ? hdf5->groups->name : NULL;
708: } else if (name[0] == '/') {
709: /* absolute */
710: for (i = 1; i < len; i++) {
711: if (name[i] != '/') {
712: gname = name;
713: break;
714: }
715: }
716: } else {
717: /* relative */
718: const char *parent = (hdf5->groups && hdf5->groups->name) ? hdf5->groups->name : "";
719: PetscCall(PetscSNPrintf(buf, sizeof(buf), "%s/%s", parent, name));
720: gname = buf;
721: }
722: }
723: PetscCall(PetscNew(&groupNode));
724: PetscCall(PetscStrallocpy(gname, (char **)&groupNode->name));
725: groupNode->next = hdf5->groups;
726: hdf5->groups = groupNode;
727: PetscFunctionReturn(PETSC_SUCCESS);
728: }
730: /*@
731: PetscViewerHDF5PopGroup - Return the current HDF5 group for output to the previous value
733: Not Collective
735: Input Parameter:
736: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
738: Level: intermediate
740: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5OpenGroup()`, `PetscViewerHDF5WriteGroup()`
741: @*/
742: PetscErrorCode PetscViewerHDF5PopGroup(PetscViewer viewer)
743: {
744: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
745: PetscViewerHDF5GroupList *groupNode;
747: PetscFunctionBegin;
749: PetscCheck(hdf5->groups, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "HDF5 group stack is empty, cannot pop");
750: groupNode = hdf5->groups;
751: hdf5->groups = hdf5->groups->next;
752: PetscCall(PetscFree(groupNode->name));
753: PetscCall(PetscFree(groupNode));
754: PetscFunctionReturn(PETSC_SUCCESS);
755: }
757: PetscErrorCode PetscViewerHDF5GetGroup_Internal(PetscViewer viewer, const char *name[])
758: {
759: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
761: PetscFunctionBegin;
764: if (hdf5->groups) *name = hdf5->groups->name;
765: else *name = NULL;
766: PetscFunctionReturn(PETSC_SUCCESS);
767: }
769: /*@C
770: PetscViewerHDF5OpenGroup - Open the HDF5 group with the name (full path) returned by `PetscViewerHDF5GetGroup()`,
771: and return this group's ID and file ID.
772: If `PetscViewerHDF5GetGroup()` yields NULL, then group ID is file ID.
774: Not Collective
776: Input Parameters:
777: + viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
778: - path - (Optional) The path relative to the pushed group
780: Output Parameters:
781: + fileId - The HDF5 file ID
782: - groupId - The HDF5 group ID
784: Level: intermediate
786: Note:
787: If path starts with '/', it is taken as an absolute path overriding currently pushed group, else path is relative to the current pushed group.
788: `NULL` or empty path means the current pushed group.
790: If the viewer is writable, the group is created if it doesn't exist yet.
792: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5WriteGroup()`
793: @*/
794: PetscErrorCode PetscViewerHDF5OpenGroup(PetscViewer viewer, const char path[], hid_t *fileId, hid_t *groupId)
795: {
796: hid_t file_id;
797: H5O_type_t type;
798: const char *fileName = NULL;
799: char *groupName = NULL;
800: PetscBool writable, has;
802: PetscFunctionBegin;
803: PetscCall(PetscViewerWritable(viewer, &writable));
804: PetscCall(PetscViewerHDF5GetFileId(viewer, &file_id));
805: PetscCall(PetscViewerFileGetName(viewer, &fileName));
806: PetscCall(PetscViewerHDF5GetGroup(viewer, path, &groupName));
807: PetscCall(PetscViewerHDF5Traverse_Internal(viewer, groupName, writable, &has, &type));
808: if (!has) {
809: PetscCheck(writable, PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Group %s does not exist and file %s is not open for writing", groupName, fileName);
810: SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_LIB, "HDF5 failed to create group %s although file %s is open for writing", groupName, fileName);
811: }
812: PetscCheck(type == H5O_TYPE_GROUP, PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Path %s in file %s resolves to something which is not a group", groupName, fileName);
813: PetscCallHDF5Return(*groupId, H5Gopen2, (file_id, groupName, H5P_DEFAULT));
814: PetscCall(PetscFree(groupName));
815: *fileId = file_id;
816: PetscFunctionReturn(PETSC_SUCCESS);
817: }
819: /*@C
820: PetscViewerHDF5WriteGroup - Ensure the HDF5 group exists in the HDF5 file
822: Not Collective
824: Input Parameters:
825: + viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
826: - path - (Optional) The path relative to the pushed group
828: Level: intermediate
830: Note:
831: If path starts with '/', it is taken as an absolute path overriding currently pushed group, else path is relative to the current pushed group.
832: `NULL` or empty path means the current pushed group.
834: This will fail if the viewer is not writable.
836: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5OpenGroup()`
837: @*/
838: PetscErrorCode PetscViewerHDF5WriteGroup(PetscViewer viewer, const char path[])
839: {
840: hid_t fileId, groupId;
842: PetscFunctionBegin;
843: PetscCall(PetscViewerHDF5OpenGroup(viewer, path, &fileId, &groupId)); // make sure group is actually created
844: PetscCallHDF5(H5Gclose, (groupId));
845: PetscFunctionReturn(PETSC_SUCCESS);
846: }
848: /*@
849: PetscViewerHDF5PushTimestepping - Activate 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: Notes:
859: On first `PetscViewerHDF5PushTimestepping()`, the initial time step is set to 0.
860: Next timesteps can then be set using `PetscViewerHDF5IncrementTimestep()` or `PetscViewerHDF5SetTimestep()`.
861: Current timestep value determines which timestep is read from or written to any dataset on the next HDF5 I/O operation [e.g. `VecView()`].
862: Use `PetscViewerHDF5PopTimestepping()` to deactivate timestepping mode; calling it by the end of the program is NOT mandatory.
863: Current timestep is remembered between `PetscViewerHDF5PopTimestepping()` and the next `PetscViewerHDF5PushTimestepping()`.
865: If a dataset was stored with timestepping, it can be loaded only in the timestepping mode again.
866: Loading a timestepped dataset with timestepping disabled, or vice-versa results in an error.
868: Developer note:
869: Timestepped HDF5 dataset has an extra dimension and attribute "timestepping" set to true.
871: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PopTimestepping()`, `PetscViewerHDF5IsTimestepping()`, `PetscViewerHDF5SetTimestep()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5GetTimestep()`
872: @*/
873: PetscErrorCode PetscViewerHDF5PushTimestepping(PetscViewer viewer)
874: {
875: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
877: PetscFunctionBegin;
879: PetscCheck(!hdf5->timestepping, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "Timestepping is already pushed");
880: hdf5->timestepping = PETSC_TRUE;
881: if (hdf5->timestep < 0) hdf5->timestep = 0;
882: PetscFunctionReturn(PETSC_SUCCESS);
883: }
885: /*@
886: PetscViewerHDF5PopTimestepping - Deactivate timestepping mode for subsequent HDF5 reading and writing.
888: Not Collective
890: Input Parameter:
891: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
893: Level: intermediate
895: Note:
896: See `PetscViewerHDF5PushTimestepping()` for details.
898: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5IsTimestepping()`, `PetscViewerHDF5SetTimestep()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5GetTimestep()`
899: @*/
900: PetscErrorCode PetscViewerHDF5PopTimestepping(PetscViewer viewer)
901: {
902: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
904: PetscFunctionBegin;
906: PetscCheck(hdf5->timestepping, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "Timestepping has not been pushed yet. Call PetscViewerHDF5PushTimestepping() first");
907: hdf5->timestepping = PETSC_FALSE;
908: PetscFunctionReturn(PETSC_SUCCESS);
909: }
911: /*@
912: PetscViewerHDF5IsTimestepping - Ask the viewer whether it is in timestepping mode currently.
914: Not Collective
916: Input Parameter:
917: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
919: Output Parameter:
920: . flg - is timestepping active?
922: Level: intermediate
924: Note:
925: See `PetscViewerHDF5PushTimestepping()` for details.
927: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5PopTimestepping()`, `PetscViewerHDF5SetTimestep()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5GetTimestep()`
928: @*/
929: PetscErrorCode PetscViewerHDF5IsTimestepping(PetscViewer viewer, PetscBool *flg)
930: {
931: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
933: PetscFunctionBegin;
935: *flg = hdf5->timestepping;
936: PetscFunctionReturn(PETSC_SUCCESS);
937: }
939: /*@
940: PetscViewerHDF5IncrementTimestep - Increments current timestep for the HDF5 output. Fields are stacked in time.
942: Not Collective
944: Input Parameter:
945: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
947: Level: intermediate
949: Note:
950: This can be called only if the viewer is in timestepping mode. See `PetscViewerHDF5PushTimestepping()` for details.
952: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5SetTimestep()`, `PetscViewerHDF5GetTimestep()`
953: @*/
954: PetscErrorCode PetscViewerHDF5IncrementTimestep(PetscViewer viewer)
955: {
956: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
958: PetscFunctionBegin;
960: PetscCheck(hdf5->timestepping, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "Timestepping has not been pushed yet. Call PetscViewerHDF5PushTimestepping() first");
961: ++hdf5->timestep;
962: PetscFunctionReturn(PETSC_SUCCESS);
963: }
965: /*@
966: PetscViewerHDF5SetTimestep - Set the current timestep for the HDF5 output. Fields are stacked in time.
968: Logically Collective
970: Input Parameters:
971: + viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
972: - timestep - The timestep
974: Level: intermediate
976: Note:
977: This can be called only if the viewer is in timestepping mode. See `PetscViewerHDF5PushTimestepping()` for details.
979: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5GetTimestep()`
980: @*/
981: PetscErrorCode PetscViewerHDF5SetTimestep(PetscViewer viewer, PetscInt timestep)
982: {
983: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
985: PetscFunctionBegin;
988: PetscCheck(timestep >= 0, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "Timestep %" PetscInt_FMT " is negative", timestep);
989: PetscCheck(hdf5->timestepping, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "Timestepping has not been pushed yet. Call PetscViewerHDF5PushTimestepping() first");
990: hdf5->timestep = timestep;
991: PetscFunctionReturn(PETSC_SUCCESS);
992: }
994: /*@
995: PetscViewerHDF5GetTimestep - Get the current timestep for the HDF5 output. Fields are stacked in time.
997: Not Collective
999: Input Parameter:
1000: . viewer - the `PetscViewer` of type `PETSCVIEWERHDF5`
1002: Output Parameter:
1003: . timestep - The timestep
1005: Level: intermediate
1007: Note:
1008: This can be called only if the viewer is in the timestepping mode. See `PetscViewerHDF5PushTimestepping()` for details.
1010: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5PushTimestepping()`, `PetscViewerHDF5IncrementTimestep()`, `PetscViewerHDF5SetTimestep()`
1011: @*/
1012: PetscErrorCode PetscViewerHDF5GetTimestep(PetscViewer viewer, PetscInt *timestep)
1013: {
1014: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
1016: PetscFunctionBegin;
1019: PetscCheck(hdf5->timestepping, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "Timestepping has not been pushed yet. Call PetscViewerHDF5PushTimestepping() first");
1020: *timestep = hdf5->timestep;
1021: PetscFunctionReturn(PETSC_SUCCESS);
1022: }
1024: /*@C
1025: PetscDataTypeToHDF5DataType - Converts the PETSc name of a datatype to its HDF5 name.
1027: Not Collective
1029: Input Parameter:
1030: . ptype - the PETSc datatype name (for example `PETSC_DOUBLE`)
1032: Output Parameter:
1033: . mtype - the HDF5 datatype
1035: Level: advanced
1037: .seealso: [](sec_viewers), `PetscDataType`, `PetscHDF5DataTypeToPetscDataType()`
1038: @*/
1039: PetscErrorCode PetscDataTypeToHDF5DataType(PetscDataType ptype, hid_t *htype)
1040: {
1041: PetscFunctionBegin;
1042: if (ptype == PETSC_INT)
1043: #if defined(PETSC_USE_64BIT_INDICES)
1044: *htype = H5T_NATIVE_LLONG;
1045: #else
1046: *htype = H5T_NATIVE_INT;
1047: #endif
1048: else if (ptype == PETSC_DOUBLE) *htype = H5T_NATIVE_DOUBLE;
1049: else if (ptype == PETSC_LONG) *htype = H5T_NATIVE_LONG;
1050: else if (ptype == PETSC_SHORT) *htype = H5T_NATIVE_SHORT;
1051: else if (ptype == PETSC_ENUM) *htype = H5T_NATIVE_INT;
1052: else if (ptype == PETSC_BOOL) *htype = H5T_NATIVE_INT;
1053: else if (ptype == PETSC_FLOAT) *htype = H5T_NATIVE_FLOAT;
1054: else if (ptype == PETSC_CHAR) *htype = H5T_NATIVE_CHAR;
1055: else if (ptype == PETSC_BIT_LOGICAL) *htype = H5T_NATIVE_UCHAR;
1056: else if (ptype == PETSC_STRING) *htype = H5Tcopy(H5T_C_S1);
1057: else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported PETSc datatype");
1058: PetscFunctionReturn(PETSC_SUCCESS);
1059: }
1061: /*@C
1062: PetscHDF5DataTypeToPetscDataType - Finds the PETSc name of a datatype from its HDF5 name
1064: Not Collective
1066: Input Parameter:
1067: . htype - the HDF5 datatype (for example `H5T_NATIVE_DOUBLE`, ...)
1069: Output Parameter:
1070: . ptype - the PETSc datatype name (for example `PETSC_DOUBLE`)
1072: Level: advanced
1074: .seealso: [](sec_viewers), `PetscDataType`, `PetscHDF5DataTypeToPetscDataType()`
1075: @*/
1076: PetscErrorCode PetscHDF5DataTypeToPetscDataType(hid_t htype, PetscDataType *ptype)
1077: {
1078: PetscFunctionBegin;
1079: #if defined(PETSC_USE_64BIT_INDICES)
1080: if (htype == H5T_NATIVE_INT) *ptype = PETSC_LONG;
1081: else if (htype == H5T_NATIVE_LLONG) *ptype = PETSC_INT;
1082: #else
1083: if (htype == H5T_NATIVE_INT) *ptype = PETSC_INT;
1084: #endif
1085: else if (htype == H5T_NATIVE_DOUBLE) *ptype = PETSC_DOUBLE;
1086: else if (htype == H5T_NATIVE_LONG) *ptype = PETSC_LONG;
1087: else if (htype == H5T_NATIVE_SHORT) *ptype = PETSC_SHORT;
1088: else if (htype == H5T_NATIVE_FLOAT) *ptype = PETSC_FLOAT;
1089: else if (htype == H5T_NATIVE_CHAR) *ptype = PETSC_CHAR;
1090: else if (htype == H5T_NATIVE_UCHAR) *ptype = PETSC_CHAR;
1091: else if (htype == H5T_C_S1) *ptype = PETSC_STRING;
1092: else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported HDF5 datatype");
1093: PetscFunctionReturn(PETSC_SUCCESS);
1094: }
1096: /*@C
1097: PetscViewerHDF5WriteAttribute - Write an attribute
1099: Collective
1101: Input Parameters:
1102: + viewer - The `PETSCVIEWERHDF5` viewer
1103: . parent - The parent dataset/group name
1104: . name - The attribute name
1105: . datatype - The attribute type
1106: - value - The attribute value
1108: Level: advanced
1110: Note:
1111: 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.
1113: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5WriteObjectAttribute()`, `PetscViewerHDF5ReadAttribute()`, `PetscViewerHDF5HasAttribute()`,
1114: `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1115: @*/
1116: PetscErrorCode PetscViewerHDF5WriteAttribute(PetscViewer viewer, const char parent[], const char name[], PetscDataType datatype, const void *value)
1117: {
1118: char *parentAbsPath;
1119: hid_t h5, dataspace, obj, attribute, dtype;
1120: PetscBool has;
1122: PetscFunctionBegin;
1128: PetscCall(PetscViewerHDF5GetGroup(viewer, parent, &parentAbsPath));
1129: PetscCall(PetscViewerHDF5Traverse_Internal(viewer, parentAbsPath, PETSC_TRUE, NULL, NULL));
1130: PetscCall(PetscViewerHDF5HasAttribute_Internal(viewer, parentAbsPath, name, &has));
1131: PetscCall(PetscDataTypeToHDF5DataType(datatype, &dtype));
1132: if (datatype == PETSC_STRING) {
1133: size_t len;
1134: PetscCall(PetscStrlen((const char *)value, &len));
1135: PetscCallHDF5(H5Tset_size, (dtype, len + 1));
1136: }
1137: PetscCall(PetscViewerHDF5GetFileId(viewer, &h5));
1138: PetscCallHDF5Return(dataspace, H5Screate, (H5S_SCALAR));
1139: PetscCallHDF5Return(obj, H5Oopen, (h5, parentAbsPath, H5P_DEFAULT));
1140: if (has) {
1141: PetscCallHDF5Return(attribute, H5Aopen_name, (obj, name));
1142: } else {
1143: PetscCallHDF5Return(attribute, H5Acreate2, (obj, name, dtype, dataspace, H5P_DEFAULT, H5P_DEFAULT));
1144: }
1145: PetscCallHDF5(H5Awrite, (attribute, dtype, value));
1146: if (datatype == PETSC_STRING) PetscCallHDF5(H5Tclose, (dtype));
1147: PetscCallHDF5(H5Aclose, (attribute));
1148: PetscCallHDF5(H5Oclose, (obj));
1149: PetscCallHDF5(H5Sclose, (dataspace));
1150: PetscCall(PetscFree(parentAbsPath));
1151: PetscFunctionReturn(PETSC_SUCCESS);
1152: }
1154: /*@C
1155: PetscViewerHDF5WriteObjectAttribute - Write an attribute to the dataset matching the given `PetscObject` by name
1157: Collective
1159: Input Parameters:
1160: + viewer - The `PETSCVIEWERHDF5` viewer
1161: . obj - The object whose name is used to lookup the parent dataset, relative to the current group.
1162: . name - The attribute name
1163: . datatype - The attribute type
1164: - value - The attribute value
1166: Level: advanced
1168: Note:
1169: 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).
1170: You might want to check first if it does using `PetscViewerHDF5HasObject()`.
1172: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5WriteAttribute()`, `PetscViewerHDF5ReadObjectAttribute()`, `PetscViewerHDF5HasObjectAttribute()`,
1173: `PetscViewerHDF5HasObject()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1174: @*/
1175: PetscErrorCode PetscViewerHDF5WriteObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscDataType datatype, const void *value)
1176: {
1177: PetscFunctionBegin;
1182: PetscCall(PetscViewerHDF5CheckNamedObject_Internal(viewer, obj));
1183: PetscCall(PetscViewerHDF5WriteAttribute(viewer, obj->name, name, datatype, value));
1184: PetscFunctionReturn(PETSC_SUCCESS);
1185: }
1187: /*@C
1188: PetscViewerHDF5ReadAttribute - Read an attribute
1190: Collective
1192: Input Parameters:
1193: + viewer - The `PETSCVIEWERHDF5` viewer
1194: . parent - The parent dataset/group name
1195: . name - The attribute name
1196: . datatype - The attribute type
1197: - defaultValue - The pointer to the default value
1199: Output Parameter:
1200: . value - The pointer to the read HDF5 attribute value
1202: Level: advanced
1204: Notes:
1205: If defaultValue is `NULL` and the attribute is not found, an error occurs.
1207: If defaultValue is not `NULL` and the attribute is not found, `defaultValue` is copied to value.
1209: The pointers `defaultValue` and value can be the same; for instance
1210: .vb
1211: flg = PETSC_FALSE;
1212: PetscCall(`PetscViewerHDF5ReadAttribute`(viewer,name,"attr",PETSC_BOOL,&flg,&flg));
1213: .ve
1214: is valid, but make sure the default value is initialized.
1216: If the datatype is `PETSC_STRING`, the output string is newly allocated so one must `PetscFree()` it when no longer needed.
1218: 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.
1220: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5ReadObjectAttribute()`, `PetscViewerHDF5WriteAttribute()`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5HasObject()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1221: @*/
1222: PetscErrorCode PetscViewerHDF5ReadAttribute(PetscViewer viewer, const char parent[], const char name[], PetscDataType datatype, const void *defaultValue, void *value)
1223: {
1224: char *parentAbsPath;
1225: hid_t h5, obj, attribute, dtype;
1226: PetscBool has;
1228: PetscFunctionBegin;
1234: PetscCall(PetscDataTypeToHDF5DataType(datatype, &dtype));
1235: PetscCall(PetscViewerHDF5GetGroup(viewer, parent, &parentAbsPath));
1236: PetscCall(PetscViewerHDF5Traverse_Internal(viewer, parentAbsPath, PETSC_FALSE, &has, NULL));
1237: if (has) PetscCall(PetscViewerHDF5HasAttribute_Internal(viewer, parentAbsPath, name, &has));
1238: if (!has) {
1239: if (defaultValue) {
1240: if (defaultValue != value) {
1241: if (datatype == PETSC_STRING) {
1242: PetscCall(PetscStrallocpy(*(char **)defaultValue, (char **)value));
1243: } else {
1244: size_t len;
1245: PetscCallHDF5ReturnNoCheck(len, H5Tget_size, (dtype));
1246: PetscCall(PetscMemcpy(value, defaultValue, len));
1247: }
1248: }
1249: PetscCall(PetscFree(parentAbsPath));
1250: PetscFunctionReturn(PETSC_SUCCESS);
1251: } else SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Attribute %s/%s does not exist and default value not provided", parentAbsPath, name);
1252: }
1253: PetscCall(PetscViewerHDF5GetFileId(viewer, &h5));
1254: PetscCallHDF5Return(obj, H5Oopen, (h5, parentAbsPath, H5P_DEFAULT));
1255: PetscCallHDF5Return(attribute, H5Aopen_name, (obj, name));
1256: if (datatype == PETSC_STRING) {
1257: size_t len;
1258: hid_t atype;
1259: PetscCallHDF5Return(atype, H5Aget_type, (attribute));
1260: PetscCallHDF5ReturnNoCheck(len, H5Tget_size, (atype));
1261: PetscCall(PetscMalloc((len + 1) * sizeof(char), value));
1262: PetscCallHDF5(H5Tset_size, (dtype, len + 1));
1263: PetscCallHDF5(H5Aread, (attribute, dtype, *(char **)value));
1264: } else {
1265: PetscCallHDF5(H5Aread, (attribute, dtype, value));
1266: }
1267: PetscCallHDF5(H5Aclose, (attribute));
1268: /* H5Oclose can be used to close groups, datasets, or committed datatypes */
1269: PetscCallHDF5(H5Oclose, (obj));
1270: PetscCall(PetscFree(parentAbsPath));
1271: PetscFunctionReturn(PETSC_SUCCESS);
1272: }
1274: /*@C
1275: PetscViewerHDF5ReadObjectAttribute - Read an attribute from the dataset matching the given `PetscObject` by name
1277: Collective
1279: Input Parameters:
1280: + viewer - The `PETSCVIEWERHDF5` viewer
1281: . obj - The object whose name is used to lookup the parent dataset, relative to the current group.
1282: . name - The attribute name
1283: - datatype - The attribute type
1285: Output Parameter:
1286: . value - The attribute value
1288: Level: advanced
1290: Note:
1291: This fails if current_group/object_name doesn't resolve to a dataset (the path doesn't exist or is not a dataset).
1292: You might want to check first if it does using `PetscViewerHDF5HasObject()`.
1294: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5ReadAttribute()` `PetscViewerHDF5WriteObjectAttribute()`, `PetscViewerHDF5HasObjectAttribute()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1295: @*/
1296: PetscErrorCode PetscViewerHDF5ReadObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscDataType datatype, void *defaultValue, void *value)
1297: {
1298: PetscFunctionBegin;
1303: PetscCall(PetscViewerHDF5CheckNamedObject_Internal(viewer, obj));
1304: PetscCall(PetscViewerHDF5ReadAttribute(viewer, obj->name, name, datatype, defaultValue, value));
1305: PetscFunctionReturn(PETSC_SUCCESS);
1306: }
1308: static inline PetscErrorCode PetscViewerHDF5Traverse_Inner_Internal(hid_t h5, const char name[], PetscBool createGroup, PetscBool *exists_)
1309: {
1310: htri_t exists;
1311: hid_t group;
1313: PetscFunctionBegin;
1314: PetscCallHDF5Return(exists, H5Lexists, (h5, name, H5P_DEFAULT));
1315: if (exists) PetscCallHDF5Return(exists, H5Oexists_by_name, (h5, name, H5P_DEFAULT));
1316: if (!exists && createGroup) {
1317: PetscCallHDF5Return(group, H5Gcreate2, (h5, name, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT));
1318: PetscCallHDF5(H5Gclose, (group));
1319: exists = PETSC_TRUE;
1320: }
1321: *exists_ = (PetscBool)exists;
1322: PetscFunctionReturn(PETSC_SUCCESS);
1323: }
1325: static PetscErrorCode PetscViewerHDF5Traverse_Internal(PetscViewer viewer, const char name[], PetscBool createGroup, PetscBool *has, H5O_type_t *otype)
1326: {
1327: const char rootGroupName[] = "/";
1328: hid_t h5;
1329: PetscBool exists = PETSC_FALSE;
1330: PetscInt i;
1331: int n;
1332: char **hierarchy;
1333: char buf[PETSC_MAX_PATH_LEN] = "";
1335: PetscFunctionBegin;
1338: else name = rootGroupName;
1339: if (has) {
1341: *has = PETSC_FALSE;
1342: }
1343: if (otype) {
1345: *otype = H5O_TYPE_UNKNOWN;
1346: }
1347: PetscCall(PetscViewerHDF5GetFileId(viewer, &h5));
1349: /*
1350: Unfortunately, H5Oexists_by_name() fails if any object in hierarchy is missing.
1351: Hence, each of them needs to be tested separately:
1352: 1) whether it's a valid link
1353: 2) whether this link resolves to an object
1354: See H5Oexists_by_name() documentation.
1355: */
1356: PetscCall(PetscStrToArray(name, '/', &n, &hierarchy));
1357: if (!n) {
1358: /* Assume group "/" always exists in accordance with HDF5 >= 1.10.0. See H5Lexists() documentation. */
1359: if (has) *has = PETSC_TRUE;
1360: if (otype) *otype = H5O_TYPE_GROUP;
1361: PetscCall(PetscStrToArrayDestroy(n, hierarchy));
1362: PetscFunctionReturn(PETSC_SUCCESS);
1363: }
1364: for (i = 0; i < n; i++) {
1365: PetscCall(PetscStrlcat(buf, "/", sizeof(buf)));
1366: PetscCall(PetscStrlcat(buf, hierarchy[i], sizeof(buf)));
1367: PetscCall(PetscViewerHDF5Traverse_Inner_Internal(h5, buf, createGroup, &exists));
1368: if (!exists) break;
1369: }
1370: PetscCall(PetscStrToArrayDestroy(n, hierarchy));
1372: /* If the object exists, get its type */
1373: if (exists && otype) {
1374: H5O_info_t info;
1376: /* We could use H5Iget_type() here but that would require opening the object. This way we only need its name. */
1377: PetscCallHDF5(H5Oget_info_by_name, (h5, name, &info, H5P_DEFAULT));
1378: *otype = info.type;
1379: }
1380: if (has) *has = exists;
1381: PetscFunctionReturn(PETSC_SUCCESS);
1382: }
1384: /*@C
1385: PetscViewerHDF5HasGroup - Check whether the current (pushed) group exists in the HDF5 file
1387: Collective
1389: Input Parameters:
1390: + viewer - The `PETSCVIEWERHDF5` viewer
1391: - path - (Optional) The path relative to the pushed group
1393: Output Parameter:
1394: . has - Flag for group existence
1396: Level: advanced
1398: Notes:
1399: If path starts with '/', it is taken as an absolute path overriding currently pushed group, else path is relative to the current pushed group.
1400: `NULL` or empty path means the current pushed group.
1402: If path exists but is not a group, `PETSC_FALSE` is returned.
1404: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5HasDataset()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`, `PetscViewerHDF5OpenGroup()`
1405: @*/
1406: PetscErrorCode PetscViewerHDF5HasGroup(PetscViewer viewer, const char path[], PetscBool *has)
1407: {
1408: H5O_type_t type;
1409: char *abspath;
1411: PetscFunctionBegin;
1415: PetscCall(PetscViewerHDF5GetGroup(viewer, path, &abspath));
1416: PetscCall(PetscViewerHDF5Traverse_Internal(viewer, abspath, PETSC_FALSE, NULL, &type));
1417: *has = (PetscBool)(type == H5O_TYPE_GROUP);
1418: PetscCall(PetscFree(abspath));
1419: PetscFunctionReturn(PETSC_SUCCESS);
1420: }
1422: /*@C
1423: PetscViewerHDF5HasDataset - Check whether a given dataset exists in the HDF5 file
1425: Collective
1427: Input Parameters:
1428: + viewer - The `PETSCVIEWERHDF5` viewer
1429: - path - The dataset path
1431: Output Parameter:
1432: . has - Flag whether dataset exists
1434: Level: advanced
1436: Notes:
1437: If path starts with '/', it is taken as an absolute path overriding currently pushed group, else path is relative to the current pushed group.
1439: If `path` is `NULL` or empty, has is set to `PETSC_FALSE`.
1441: If `path` exists but is not a dataset, has is set to `PETSC_FALSE` as well.
1443: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5HasObject()`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5HasGroup()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1444: @*/
1445: PetscErrorCode PetscViewerHDF5HasDataset(PetscViewer viewer, const char path[], PetscBool *has)
1446: {
1447: H5O_type_t type;
1448: char *abspath;
1450: PetscFunctionBegin;
1454: PetscCall(PetscViewerHDF5GetGroup(viewer, path, &abspath));
1455: PetscCall(PetscViewerHDF5Traverse_Internal(viewer, abspath, PETSC_FALSE, NULL, &type));
1456: *has = (PetscBool)(type == H5O_TYPE_DATASET);
1457: PetscCall(PetscFree(abspath));
1458: PetscFunctionReturn(PETSC_SUCCESS);
1459: }
1461: /*@
1462: PetscViewerHDF5HasObject - Check whether a dataset with the same name as given object exists in the HDF5 file under current group
1464: Collective
1466: Input Parameters:
1467: + viewer - The `PETSCVIEWERHDF5` viewer
1468: - obj - The named object
1470: Output Parameter:
1471: . has - Flag for dataset existence
1473: Level: advanced
1475: Notes:
1476: If the object is unnamed, an error occurs.
1478: If the path current_group/object_name exists but is not a dataset, has is set to `PETSC_FALSE` as well.
1480: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5HasDataset()`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1481: @*/
1482: PetscErrorCode PetscViewerHDF5HasObject(PetscViewer viewer, PetscObject obj, PetscBool *has)
1483: {
1484: size_t len;
1486: PetscFunctionBegin;
1490: PetscCall(PetscStrlen(obj->name, &len));
1491: PetscCheck(len, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONG, "Object must be named");
1492: PetscCall(PetscViewerHDF5HasDataset(viewer, obj->name, has));
1493: PetscFunctionReturn(PETSC_SUCCESS);
1494: }
1496: /*@C
1497: PetscViewerHDF5HasAttribute - Check whether an attribute exists
1499: Collective
1501: Input Parameters:
1502: + viewer - The `PETSCVIEWERHDF5` viewer
1503: . parent - The parent dataset/group name
1504: - name - The attribute name
1506: Output Parameter:
1507: . has - Flag for attribute existence
1509: Level: advanced
1511: Note:
1512: 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.
1514: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5HasObjectAttribute()`, `PetscViewerHDF5WriteAttribute()`, `PetscViewerHDF5ReadAttribute()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1515: @*/
1516: PetscErrorCode PetscViewerHDF5HasAttribute(PetscViewer viewer, const char parent[], const char name[], PetscBool *has)
1517: {
1518: char *parentAbsPath;
1520: PetscFunctionBegin;
1525: PetscCall(PetscViewerHDF5GetGroup(viewer, parent, &parentAbsPath));
1526: PetscCall(PetscViewerHDF5Traverse_Internal(viewer, parentAbsPath, PETSC_FALSE, has, NULL));
1527: if (*has) PetscCall(PetscViewerHDF5HasAttribute_Internal(viewer, parentAbsPath, name, has));
1528: PetscCall(PetscFree(parentAbsPath));
1529: PetscFunctionReturn(PETSC_SUCCESS);
1530: }
1532: /*@C
1533: PetscViewerHDF5HasObjectAttribute - Check whether an attribute is attached to the dataset matching the given `PetscObject` by name
1535: Collective
1537: Input Parameters:
1538: + viewer - The `PETSCVIEWERHDF5` viewer
1539: . obj - The object whose name is used to lookup the parent dataset, relative to the current group.
1540: - name - The attribute name
1542: Output Parameter:
1543: . has - Flag for attribute existence
1545: Level: advanced
1547: Note:
1548: This fails if current_group/object_name doesn't resolve to a dataset (the path doesn't exist or is not a dataset).
1549: You might want to check first if it does using `PetscViewerHDF5HasObject()`.
1551: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerHDF5HasAttribute()`, `PetscViewerHDF5WriteObjectAttribute()`, `PetscViewerHDF5ReadObjectAttribute()`, `PetscViewerHDF5HasObject()`, `PetscViewerHDF5PushGroup()`, `PetscViewerHDF5PopGroup()`, `PetscViewerHDF5GetGroup()`
1552: @*/
1553: PetscErrorCode PetscViewerHDF5HasObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscBool *has)
1554: {
1555: PetscFunctionBegin;
1560: PetscCall(PetscViewerHDF5CheckNamedObject_Internal(viewer, obj));
1561: PetscCall(PetscViewerHDF5HasAttribute(viewer, obj->name, name, has));
1562: PetscFunctionReturn(PETSC_SUCCESS);
1563: }
1565: static PetscErrorCode PetscViewerHDF5HasAttribute_Internal(PetscViewer viewer, const char parent[], const char name[], PetscBool *has)
1566: {
1567: hid_t h5;
1568: htri_t hhas;
1570: PetscFunctionBegin;
1571: PetscCall(PetscViewerHDF5GetFileId(viewer, &h5));
1572: PetscCallHDF5Return(hhas, H5Aexists_by_name, (h5, parent, name, H5P_DEFAULT));
1573: *has = hhas ? PETSC_TRUE : PETSC_FALSE;
1574: PetscFunctionReturn(PETSC_SUCCESS);
1575: }
1577: /*
1578: The variable Petsc_Viewer_HDF5_keyval is used to indicate an MPI attribute that
1579: is attached to a communicator, in this case the attribute is a PetscViewer.
1580: */
1581: PetscMPIInt Petsc_Viewer_HDF5_keyval = MPI_KEYVAL_INVALID;
1583: /*@C
1584: PETSC_VIEWER_HDF5_ - Creates an `PETSCVIEWERHDF5` `PetscViewer` shared by all processors in a communicator.
1586: Collective
1588: Input Parameter:
1589: . comm - the MPI communicator to share the `PETSCVIEWERHDF5` `PetscViewer`
1591: Options Database Key:
1592: . -viewer_hdf5_filename <name> - name of the HDF5 file
1594: Environmental variable:
1595: . `PETSC_VIEWER_HDF5_FILENAME` - name of the HDF5 file
1597: Level: intermediate
1599: Note:
1600: Unlike almost all other PETSc routines, `PETSC_VIEWER_HDF5_()` does not return
1601: an error code. The HDF5 `PetscViewer` is usually used in the form
1602: $ XXXView(XXX object, PETSC_VIEWER_HDF5_(comm));
1604: .seealso: [](sec_viewers), `PETSCVIEWERHDF5`, `PetscViewerHDF5Open()`, `PetscViewerCreate()`, `PetscViewerDestroy()`
1605: @*/
1606: PetscViewer PETSC_VIEWER_HDF5_(MPI_Comm comm)
1607: {
1608: PetscErrorCode ierr;
1609: PetscMPIInt mpi_ierr;
1610: PetscBool flg;
1611: PetscViewer viewer;
1612: char fname[PETSC_MAX_PATH_LEN];
1613: MPI_Comm ncomm;
1615: PetscFunctionBegin;
1616: ierr = PetscCommDuplicate(comm, &ncomm, NULL);
1617: if (ierr) {
1618: ierr = PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_INITIAL, " ");
1619: PetscFunctionReturn(NULL);
1620: }
1621: if (Petsc_Viewer_HDF5_keyval == MPI_KEYVAL_INVALID) {
1622: mpi_ierr = MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN, MPI_COMM_NULL_DELETE_FN, &Petsc_Viewer_HDF5_keyval, NULL);
1623: if (mpi_ierr) {
1624: ierr = PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_INITIAL, " ");
1625: PetscFunctionReturn(NULL);
1626: }
1627: }
1628: mpi_ierr = MPI_Comm_get_attr(ncomm, Petsc_Viewer_HDF5_keyval, (void **)&viewer, (int *)&flg);
1629: if (mpi_ierr) {
1630: ierr = PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_INITIAL, " ");
1631: PetscFunctionReturn(NULL);
1632: }
1633: if (!flg) { /* PetscViewer not yet created */
1634: ierr = PetscOptionsGetenv(ncomm, "PETSC_VIEWER_HDF5_FILENAME", fname, PETSC_MAX_PATH_LEN, &flg);
1635: if (ierr) {
1636: ierr = PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1637: PetscFunctionReturn(NULL);
1638: }
1639: if (!flg) {
1640: ierr = PetscStrncpy(fname, "output.h5", sizeof(fname));
1641: if (ierr) {
1642: ierr = PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1643: PetscFunctionReturn(NULL);
1644: }
1645: }
1646: ierr = PetscViewerHDF5Open(ncomm, fname, FILE_MODE_WRITE, &viewer);
1647: if (ierr) {
1648: ierr = PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1649: PetscFunctionReturn(NULL);
1650: }
1651: ierr = PetscObjectRegisterDestroy((PetscObject)viewer);
1652: if (ierr) {
1653: ierr = PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1654: PetscFunctionReturn(NULL);
1655: }
1656: mpi_ierr = MPI_Comm_set_attr(ncomm, Petsc_Viewer_HDF5_keyval, (void *)viewer);
1657: if (mpi_ierr) {
1658: ierr = PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_INITIAL, " ");
1659: PetscFunctionReturn(NULL);
1660: }
1661: }
1662: ierr = PetscCommDestroy(&ncomm);
1663: if (ierr) {
1664: ierr = PetscError(PETSC_COMM_SELF, __LINE__, "PETSC_VIEWER_HDF5_", __FILE__, PETSC_ERR_PLIB, PETSC_ERROR_REPEAT, " ");
1665: PetscFunctionReturn(NULL);
1666: }
1667: PetscFunctionReturn(viewer);
1668: }