Actual source code: hdf5v.c

petsc-main 2021-04-20
Report Typos and Errors
  1: #include <petsc/private/viewerimpl.h>
  2: #include <petsc/private/viewerhdf5impl.h>
  3: #include <petscviewerhdf5.h>

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

  8: static PetscErrorCode PetscViewerHDF5GetAbsolutePath_Internal(PetscViewer viewer, const char path[], char *abspath[])
  9: {
 10:   PetscBool      relative = PETSC_FALSE;
 11:   const char     *group;
 12:   char           buf[PETSC_MAX_PATH_LEN] = "";

 16:   PetscViewerHDF5GetGroup(viewer, &group);
 17:   PetscViewerHDF5PathIsRelative(path, PETSC_TRUE, &relative);
 18:   if (relative) {
 19:     PetscStrcpy(buf, group);
 20:     PetscStrcat(buf, "/");
 21:     PetscStrcat(buf, path);
 22:     PetscStrallocpy(buf, abspath);
 23:   } else {
 24:     PetscStrallocpy(path, abspath);
 25:   }
 26:   return(0);
 27: }

 29: static PetscErrorCode PetscViewerHDF5CheckNamedObject_Internal(PetscViewer viewer, PetscObject obj)
 30: {
 31:   PetscBool      has;

 35:   PetscViewerHDF5HasObject(viewer, obj, &has);
 36:   if (!has) {
 37:     const char *group;
 38:     PetscViewerHDF5GetGroup(viewer, &group);
 39:     SETERRQ2(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Object (dataset) \"%s\" not stored in group %s", obj->name, group ? group : "/");
 40:   }
 41:   return(0);
 42: }

 44: static PetscErrorCode PetscViewerSetFromOptions_HDF5(PetscOptionItems *PetscOptionsObject,PetscViewer v)
 45: {
 46:   PetscErrorCode   ierr;
 47:   PetscBool        flg = PETSC_FALSE, set;
 48:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*)v->data;

 51:   PetscOptionsHead(PetscOptionsObject,"HDF5 PetscViewer Options");
 52:   PetscOptionsBool("-viewer_hdf5_base_dimension2","1d Vectors get 2 dimensions in HDF5","PetscViewerHDF5SetBaseDimension2",hdf5->basedimension2,&hdf5->basedimension2,NULL);
 53:   PetscOptionsBool("-viewer_hdf5_sp_output","Force data to be written in single precision","PetscViewerHDF5SetSPOutput",hdf5->spoutput,&hdf5->spoutput,NULL);
 54:   PetscOptionsBool("-viewer_hdf5_collective","Enable collective transfer mode","PetscViewerHDF5SetCollective",flg,&flg,&set);
 55:   if (set) {PetscViewerHDF5SetCollective(v,flg);}
 56:   PetscOptionsTail();
 57:   return(0);
 58: }

 60: static PetscErrorCode PetscViewerView_HDF5(PetscViewer v,PetscViewer viewer)
 61: {
 62:   PetscViewer_HDF5  *hdf5 = (PetscViewer_HDF5*)v->data;
 63:   PetscBool         flg;
 64:   PetscErrorCode    ierr;

 67:   if (hdf5->filename) {
 68:     PetscViewerASCIIPrintf(viewer,"Filename: %s\n",hdf5->filename);
 69:   }
 70:   PetscViewerASCIIPrintf(viewer,"Vectors with blocksize 1 saved as 2D datasets: %s\n",PetscBools[hdf5->basedimension2]);
 71:   PetscViewerASCIIPrintf(viewer,"Enforce single precision storage: %s\n",PetscBools[hdf5->spoutput]);
 72:   PetscViewerHDF5GetCollective(v,&flg);
 73:   PetscViewerASCIIPrintf(viewer,"MPI-IO transfer mode: %s\n",flg ? "collective" : "independent");
 74:   return(0);
 75: }

 77: static PetscErrorCode PetscViewerFileClose_HDF5(PetscViewer viewer)
 78: {
 79:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*)viewer->data;
 80:   PetscErrorCode   ierr;

 83:   PetscFree(hdf5->filename);
 84:   if (hdf5->file_id) PetscStackCallHDF5(H5Fclose,(hdf5->file_id));
 85:   return(0);
 86: }

 88: static PetscErrorCode PetscViewerDestroy_HDF5(PetscViewer viewer)
 89: {
 90:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
 91:   PetscErrorCode   ierr;

 94:   PetscStackCallHDF5(H5Pclose,(hdf5->dxpl_id));
 95:   PetscViewerFileClose_HDF5(viewer);
 96:   while (hdf5->groups) {
 97:     PetscViewerHDF5GroupList *tmp = hdf5->groups->next;

 99:     PetscFree(hdf5->groups->name);
100:     PetscFree(hdf5->groups);
101:     hdf5->groups = tmp;
102:   }
103:   PetscFree(hdf5);
104:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetName_C",NULL);
105:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileGetName_C",NULL);
106:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetMode_C",NULL);
107:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerHDF5SetBaseDimension2_C",NULL);
108:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerHDF5SetSPOutput_C",NULL);
109:   return(0);
110: }

112: static PetscErrorCode  PetscViewerFileSetMode_HDF5(PetscViewer viewer, PetscFileMode type)
113: {
114:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

117:   hdf5->btype = type;
118:   return(0);
119: }

121: static PetscErrorCode  PetscViewerFileGetMode_HDF5(PetscViewer viewer, PetscFileMode *type)
122: {
123:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

126:   *type = hdf5->btype;
127:   return(0);
128: }

130: static PetscErrorCode  PetscViewerHDF5SetBaseDimension2_HDF5(PetscViewer viewer, PetscBool flg)
131: {
132:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

135:   hdf5->basedimension2 = flg;
136:   return(0);
137: }

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

143:     Logically Collective on PetscViewer

145:   Input Parameters:
146: +  viewer - the PetscViewer; if it is not hdf5 then this command is ignored
147: -  flg - if PETSC_TRUE the vector will always have at least a dimension of 2 even if that first dimension is of size 1

149:   Options Database:
150: .  -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


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

157:   Level: intermediate

159: .seealso: PetscViewerFileSetMode(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerBinaryOpen()

161: @*/
162: PetscErrorCode PetscViewerHDF5SetBaseDimension2(PetscViewer viewer,PetscBool flg)
163: {

168:   PetscTryMethod(viewer,"PetscViewerHDF5SetBaseDimension2_C",(PetscViewer,PetscBool),(viewer,flg));
169:   return(0);
170: }

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

176:     Logically Collective on PetscViewer

178:   Input Parameter:
179: .  viewer - the PetscViewer, must be of type HDF5

181:   Output Parameter:
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:   Notes:
185:     Setting this option allegedly makes code that reads the HDF5 in easier since they do not have a "special case" of a bs/dof
186:          of one when the dimension is lower. Others think the option is crazy.

188:   Level: intermediate

190: .seealso: PetscViewerFileSetMode(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerBinaryOpen()

192: @*/
193: PetscErrorCode PetscViewerHDF5GetBaseDimension2(PetscViewer viewer,PetscBool *flg)
194: {
195:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

199:   *flg = hdf5->basedimension2;
200:   return(0);
201: }

203: static PetscErrorCode  PetscViewerHDF5SetSPOutput_HDF5(PetscViewer viewer, PetscBool flg)
204: {
205:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

208:   hdf5->spoutput = flg;
209:   return(0);
210: }

212: /*@
213:      PetscViewerHDF5SetSPOutput - Data is written to disk in single precision even if PETSc is
214:        compiled with double precision PetscReal.

216:     Logically Collective on PetscViewer

218:   Input Parameters:
219: +  viewer - the PetscViewer; if it is not hdf5 then this command is ignored
220: -  flg - if PETSC_TRUE the data will be written to disk with single precision

222:   Options Database:
223: .  -viewer_hdf5_sp_output - turns on (true) or off (false) output in single precision


226:   Notes:
227:     Setting this option does not make any difference if PETSc is compiled with single precision
228:          in the first place. It does not affect reading datasets (HDF5 handle this internally).

230:   Level: intermediate

232: .seealso: PetscViewerFileSetMode(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerBinaryOpen(),
233:           PetscReal

235: @*/
236: PetscErrorCode PetscViewerHDF5SetSPOutput(PetscViewer viewer,PetscBool flg)
237: {

242:   PetscTryMethod(viewer,"PetscViewerHDF5SetSPOutput_C",(PetscViewer,PetscBool),(viewer,flg));
243:   return(0);
244: }

246: /*@
247:      PetscViewerHDF5GetSPOutput - Data is written to disk in single precision even if PETSc is
248:        compiled with double precision PetscReal.

250:     Logically Collective on PetscViewer

252:   Input Parameter:
253: .  viewer - the PetscViewer, must be of type HDF5

255:   Output Parameter:
256: .  flg - if PETSC_TRUE the data will be written to disk with single precision

258:   Notes:
259:     Setting this option does not make any difference if PETSc is compiled with single precision
260:          in the first place. It does not affect reading datasets (HDF5 handle this internally).

262:   Level: intermediate

264: .seealso: PetscViewerFileSetMode(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerBinaryOpen(),
265:           PetscReal

267: @*/
268: PetscErrorCode PetscViewerHDF5GetSPOutput(PetscViewer viewer,PetscBool *flg)
269: {
270:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

274:   *flg = hdf5->spoutput;
275:   return(0);
276: }

278: static PetscErrorCode  PetscViewerHDF5SetCollective_HDF5(PetscViewer viewer, PetscBool flg)
279: {
281:   /* H5FD_MPIO_COLLECTIVE is wrong in hdf5 1.10.2, and is the same as H5FD_MPIO_INDEPENDENT in earlier versions
282:      - see e.g. https://gitlab.cosma.dur.ac.uk/swift/swiftsim/issues/431 */
283: #if H5_VERSION_GE(1,10,3) && defined(H5_HAVE_PARALLEL)
284:   {
285:     PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
286:     PetscStackCallHDF5(H5Pset_dxpl_mpio,(hdf5->dxpl_id, flg ? H5FD_MPIO_COLLECTIVE : H5FD_MPIO_INDEPENDENT));
287:   }
288: #else
289:   if (flg) {
291:     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");
292:   }
293: #endif
294:   return(0);
295: }

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

300:   Logically Collective; flg must contain common value

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

306:   Options Database:
307: . -viewer_hdf5_collective - turns on (true) or off (false) collective transfers

309:   Notes:
310:   Collective mode gives the MPI-IO layer underneath HDF5 a chance to do some additional collective optimizations and hence can perform better.
311:   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.

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

318:   Level: intermediate

320: .seealso: PetscViewerHDF5GetCollective(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerHDF5Open()

322: @*/
323: PetscErrorCode PetscViewerHDF5SetCollective(PetscViewer viewer,PetscBool flg)
324: {

330:   PetscTryMethod(viewer,"PetscViewerHDF5SetCollective_C",(PetscViewer,PetscBool),(viewer,flg));
331:   return(0);
332: }

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

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

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

354:   Not Collective

356:   Input Parameters:
357: . viewer - the HDF5 PetscViewer

359:   Output Parameters:
360: . flg - the flag

362:   Level: intermediate

364:   Notes:
365:   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.
366:   For more details, see PetscViewerHDF5SetCollective().

368: .seealso: PetscViewerHDF5SetCollective(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerHDF5Open()

370: @*/
371: PetscErrorCode PetscViewerHDF5GetCollective(PetscViewer viewer,PetscBool *flg)
372: {


379:   PetscUseMethod(viewer,"PetscViewerHDF5GetCollective_C",(PetscViewer,PetscBool*),(viewer,flg));
380:   return(0);
381: }

383: static PetscErrorCode  PetscViewerFileSetName_HDF5(PetscViewer viewer, const char name[])
384: {
385:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
386:   hid_t             plist_id;
387:   PetscErrorCode    ierr;

390:   if (hdf5->file_id) PetscStackCallHDF5(H5Fclose,(hdf5->file_id));
391:   if (hdf5->filename) {PetscFree(hdf5->filename);}
392:   PetscStrallocpy(name, &hdf5->filename);
393:   /* Set up file access property list with parallel I/O access */
394:   PetscStackCallHDF5Return(plist_id,H5Pcreate,(H5P_FILE_ACCESS));
395: #if defined(H5_HAVE_PARALLEL)
396:   PetscStackCallHDF5(H5Pset_fapl_mpio,(plist_id, PetscObjectComm((PetscObject)viewer), MPI_INFO_NULL));
397: #endif
398:   /* Create or open the file collectively */
399:   switch (hdf5->btype) {
400:   case FILE_MODE_READ:
401:     PetscStackCallHDF5Return(hdf5->file_id,H5Fopen,(name, H5F_ACC_RDONLY, plist_id));
402:     break;
403:   case FILE_MODE_APPEND:
404:   case FILE_MODE_UPDATE:
405:   {
406:     PetscBool flg;
407:     PetscTestFile(hdf5->filename, 'r', &flg);
408:     if (flg) PetscStackCallHDF5Return(hdf5->file_id,H5Fopen,(name, H5F_ACC_RDWR, plist_id));
409:     else     PetscStackCallHDF5Return(hdf5->file_id,H5Fcreate,(name, H5F_ACC_EXCL, H5P_DEFAULT, plist_id));
410:     break;
411:   }
412:   case FILE_MODE_WRITE:
413:     PetscStackCallHDF5Return(hdf5->file_id,H5Fcreate,(name, H5F_ACC_TRUNC, H5P_DEFAULT, plist_id));
414:     break;
415:   case FILE_MODE_UNDEFINED:
416:     SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
417:   default:
418:     SETERRQ1(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP, "Unsupported file mode %s",PetscFileModes[hdf5->btype]);
419:   }
420:   if (hdf5->file_id < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB, "H5Fcreate failed for %s", name);
421:   PetscStackCallHDF5(H5Pclose,(plist_id));
422:   return(0);
423: }

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

430:   *name = vhdf5->filename;
431:   return(0);
432: }

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

442:   return(0);
443: }

445: /*MC
446:    PETSCVIEWERHDF5 - A viewer that writes to an HDF5 file


449: .seealso:  PetscViewerHDF5Open(), PetscViewerStringSPrintf(), PetscViewerSocketOpen(), PetscViewerDrawOpen(), PETSCVIEWERSOCKET,
450:            PetscViewerCreate(), PetscViewerASCIIOpen(), PetscViewerBinaryOpen(), PETSCVIEWERBINARY, PETSCVIEWERDRAW, PETSCVIEWERSTRING,
451:            PetscViewerMatlabOpen(), VecView(), DMView(), PetscViewerMatlabPutArray(), PETSCVIEWERASCII, PETSCVIEWERMATLAB,
452:            PetscViewerFileSetName(), PetscViewerFileSetMode(), PetscViewerFormat, PetscViewerType, PetscViewerSetType()

454:   Level: beginner
455: M*/

457: PETSC_EXTERN PetscErrorCode PetscViewerCreate_HDF5(PetscViewer v)
458: {
459:   PetscViewer_HDF5 *hdf5;
460:   PetscErrorCode   ierr;

463:   PetscNewLog(v,&hdf5);

465:   v->data                = (void*) hdf5;
466:   v->ops->destroy        = PetscViewerDestroy_HDF5;
467:   v->ops->setfromoptions = PetscViewerSetFromOptions_HDF5;
468:   v->ops->setup          = PetscViewerSetUp_HDF5;
469:   v->ops->view           = PetscViewerView_HDF5;
470:   v->ops->flush          = NULL;
471:   hdf5->btype            = FILE_MODE_UNDEFINED;
472:   hdf5->filename         = NULL;
473:   hdf5->timestep         = -1;
474:   hdf5->groups           = NULL;

476:   PetscStackCallHDF5Return(hdf5->dxpl_id,H5Pcreate,(H5P_DATASET_XFER));

478:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileSetName_C",PetscViewerFileSetName_HDF5);
479:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileGetName_C",PetscViewerFileGetName_HDF5);
480:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileSetMode_C",PetscViewerFileSetMode_HDF5);
481:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileGetMode_C",PetscViewerFileGetMode_HDF5);
482:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerHDF5SetBaseDimension2_C",PetscViewerHDF5SetBaseDimension2_HDF5);
483:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerHDF5SetSPOutput_C",PetscViewerHDF5SetSPOutput_HDF5);
484:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerHDF5SetCollective_C",PetscViewerHDF5SetCollective_HDF5);
485:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerHDF5GetCollective_C",PetscViewerHDF5GetCollective_HDF5);
486:   return(0);
487: }

489: /*@C
490:    PetscViewerHDF5Open - Opens a file for HDF5 input/output.

492:    Collective

494:    Input Parameters:
495: +  comm - MPI communicator
496: .  name - name of file
497: -  type - type of file

499:    Output Parameter:
500: .  hdf5v - PetscViewer for HDF5 input/output to use with the specified file

502:   Options Database:
503: +  -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
504: -  -viewer_hdf5_sp_output - forces (if true) the viewer to write data in single precision independent on the precision of PetscReal

506:    Level: beginner

508:    Notes:
509:    Reading is always available, regardless of the mode. Available modes are
510: +  FILE_MODE_READ - open existing HDF5 file for read only access, fail if file does not exist [H5Fopen() with H5F_ACC_RDONLY]
511: .  FILE_MODE_WRITE - if file exists, fully overwrite it, else create new HDF5 file [H5FcreateH5Fcreate() with H5F_ACC_TRUNC]
512: .  FILE_MODE_APPEND - if file exists, keep existing contents [H5Fopen() with H5F_ACC_RDWR], else create new HDF5 file [H5FcreateH5Fcreate() with H5F_ACC_EXCL]
513: -  FILE_MODE_UPDATE - same as FILE_MODE_APPEND

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

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


520: .seealso: PetscViewerASCIIOpen(), PetscViewerPushFormat(), PetscViewerDestroy(), PetscViewerHDF5SetBaseDimension2(),
521:           PetscViewerHDF5SetSPOutput(), PetscViewerHDF5GetBaseDimension2(), VecView(), MatView(), VecLoad(),
522:           MatLoad(), PetscFileMode, PetscViewer, PetscViewerSetType(), PetscViewerFileSetMode(), PetscViewerFileSetName()
523: @*/
524: PetscErrorCode  PetscViewerHDF5Open(MPI_Comm comm, const char name[], PetscFileMode type, PetscViewer *hdf5v)
525: {

529:   PetscViewerCreate(comm, hdf5v);
530:   PetscViewerSetType(*hdf5v, PETSCVIEWERHDF5);
531:   PetscViewerFileSetMode(*hdf5v, type);
532:   PetscViewerFileSetName(*hdf5v, name);
533:   PetscViewerSetFromOptions(*hdf5v);
534:   return(0);
535: }

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

540:   Not collective

542:   Input Parameter:
543: . viewer - the PetscViewer

545:   Output Parameter:
546: . file_id - The file id

548:   Level: intermediate

550: .seealso: PetscViewerHDF5Open()
551: @*/
552: PetscErrorCode  PetscViewerHDF5GetFileId(PetscViewer viewer, hid_t *file_id)
553: {
554:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

558:   if (file_id) *file_id = hdf5->file_id;
559:   return(0);
560: }

562: /*@C
563:   PetscViewerHDF5PushGroup - Set the current HDF5 group for output

565:   Not collective

567:   Input Parameters:
568: + viewer - the PetscViewer
569: - name - The group name

571:   Level: intermediate

573:   Notes:
574:   This is designed to mnemonically resemble the Unix cd command.
575:   + 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.
576:   . NULL, empty string, or any sequence of all slashes (e.g. "///") is interpreted as the root group "/".
577:   - "." means the current group is pushed again.

579:   Example:
580:   Suppose the current group is "/a".
581:   + If name is NULL, empty string, or a sequence of all slashes (e.g. "///"), then the new group will be "/".
582:   . If name is ".", then the new group will be "/a".
583:   . If name is "b", then the new group will be "/a/b".
584:   - If name is "/b", then the new group will be "/b".

586:   Developer Notes:
587:   The root group "/" is internally stored as NULL.

589: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup(),PetscViewerHDF5OpenGroup()
590: @*/
591: PetscErrorCode  PetscViewerHDF5PushGroup(PetscViewer viewer, const char name[])
592: {
593:   PetscViewer_HDF5          *hdf5 = (PetscViewer_HDF5*) viewer->data;
594:   PetscViewerHDF5GroupList  *groupNode;
595:   size_t                    i,len;
596:   char                      buf[PETSC_MAX_PATH_LEN];
597:   const char                *gname;
598:   PetscErrorCode            ierr;

603:   PetscStrlen(name, &len);
604:   gname = NULL;
605:   if (len) {
606:     if (len == 1 && name[0] == '.') {
607:       /* use current name */
608:       gname = (hdf5->groups && hdf5->groups->name) ? hdf5->groups->name : NULL;
609:     } else if (name[0] == '/') {
610:       /* absolute */
611:       for (i=1; i<len; i++) {
612:         if (name[i] != '/') {
613:           gname = name;
614:           break;
615:         }
616:       }
617:     } else {
618:       /* relative */
619:       const char *parent = (hdf5->groups && hdf5->groups->name) ? hdf5->groups->name : "";
620:       PetscSNPrintf(buf, sizeof(buf), "%s/%s", parent, name);
621:       gname = buf;
622:     }
623:   }
624:   PetscNew(&groupNode);
625:   PetscStrallocpy(gname, (char**) &groupNode->name);
626:   groupNode->next = hdf5->groups;
627:   hdf5->groups    = groupNode;
628:   return(0);
629: }

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

634:   Not collective

636:   Input Parameter:
637: . viewer - the PetscViewer

639:   Level: intermediate

641: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5GetGroup(),PetscViewerHDF5OpenGroup()
642: @*/
643: PetscErrorCode  PetscViewerHDF5PopGroup(PetscViewer viewer)
644: {
645:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
646:   PetscViewerHDF5GroupList *groupNode;
647:   PetscErrorCode   ierr;

651:   if (!hdf5->groups) SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "HDF5 group stack is empty, cannot pop");
652:   groupNode    = hdf5->groups;
653:   hdf5->groups = hdf5->groups->next;
654:   PetscFree(groupNode->name);
655:   PetscFree(groupNode);
656:   return(0);
657: }

659: /*@C
660:   PetscViewerHDF5GetGroup - Get the current HDF5 group name (full path), set with PetscViewerHDF5PushGroup()/PetscViewerHDF5PopGroup().
661:   If none has been assigned, returns NULL.

663:   Not collective

665:   Input Parameter:
666: . viewer - the PetscViewer

668:   Output Parameter:
669: . name - The group name

671:   Level: intermediate

673: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5OpenGroup()
674: @*/
675: PetscErrorCode  PetscViewerHDF5GetGroup(PetscViewer viewer, const char *name[])
676: {
677:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;

682:   if (hdf5->groups) *name = hdf5->groups->name;
683:   else *name = NULL;
684:   return(0);
685: }

687: /*@
688:   PetscViewerHDF5OpenGroup - Open the HDF5 group with the name (full path) returned by PetscViewerHDF5GetGroup(),
689:   and return this group's ID and file ID.
690:   If PetscViewerHDF5GetGroup() yields NULL, then group ID is file ID.

692:   Not collective

694:   Input Parameter:
695: . viewer - the PetscViewer

697:   Output Parameter:
698: + fileId - The HDF5 file ID
699: - groupId - The HDF5 group ID

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

704:   Level: intermediate

706: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
707: @*/
708: PetscErrorCode PetscViewerHDF5OpenGroup(PetscViewer viewer, hid_t *fileId, hid_t *groupId)
709: {
710:   hid_t          file_id;
711:   H5O_type_t     type;
712:   const char     *groupName = NULL;
713:   PetscBool      create;

717:   PetscViewerWritable(viewer, &create);
718:   PetscViewerHDF5GetFileId(viewer, &file_id);
719:   PetscViewerHDF5GetGroup(viewer, &groupName);
720:   PetscViewerHDF5Traverse_Internal(viewer, groupName, create, NULL, &type);
721:   if (type != H5O_TYPE_GROUP) SETERRQ1(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Path %s resolves to something which is not a group", groupName);
722:   PetscStackCallHDF5Return(*groupId,H5Gopen2,(file_id, groupName ? groupName : "/", H5P_DEFAULT));
723:   *fileId  = file_id;
724:   return(0);
725: }

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

730:   Not collective

732:   Input Parameter:
733: . viewer - the PetscViewer

735:   Level: intermediate

737: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5SetTimestep(), PetscViewerHDF5GetTimestep()
738: @*/
739: PetscErrorCode PetscViewerHDF5IncrementTimestep(PetscViewer viewer)
740: {
741:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

745:   ++hdf5->timestep;
746:   return(0);
747: }

749: /*@
750:   PetscViewerHDF5SetTimestep - Set the current timestep for the HDF5 output. Fields are stacked in time. A timestep
751:   of -1 disables blocking with timesteps.

753:   Not collective

755:   Input Parameters:
756: + viewer - the PetscViewer
757: - timestep - The timestep number

759:   Level: intermediate

761: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5GetTimestep()
762: @*/
763: PetscErrorCode  PetscViewerHDF5SetTimestep(PetscViewer viewer, PetscInt timestep)
764: {
765:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

769:   hdf5->timestep = timestep;
770:   return(0);
771: }

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

776:   Not collective

778:   Input Parameter:
779: . viewer - the PetscViewer

781:   Output Parameter:
782: . timestep - The timestep number

784:   Level: intermediate

786: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5SetTimestep()
787: @*/
788: PetscErrorCode  PetscViewerHDF5GetTimestep(PetscViewer viewer, PetscInt *timestep)
789: {
790:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

795:   *timestep = hdf5->timestep;
796:   return(0);
797: }

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

802:   Not collective

804:   Input Parameter:
805: . ptype - the PETSc datatype name (for example PETSC_DOUBLE)

807:   Output Parameter:
808: . mtype - the MPI datatype (for example MPI_DOUBLE, ...)

810:   Level: advanced

812: .seealso: PetscDataType, PetscHDF5DataTypeToPetscDataType()
813: @*/
814: PetscErrorCode PetscDataTypeToHDF5DataType(PetscDataType ptype, hid_t *htype)
815: {
817:   if (ptype == PETSC_INT)
818: #if defined(PETSC_USE_64BIT_INDICES)
819:                                        *htype = H5T_NATIVE_LLONG;
820: #else
821:                                        *htype = H5T_NATIVE_INT;
822: #endif
823:   else if (ptype == PETSC_DOUBLE)      *htype = H5T_NATIVE_DOUBLE;
824:   else if (ptype == PETSC_LONG)        *htype = H5T_NATIVE_LONG;
825:   else if (ptype == PETSC_SHORT)       *htype = H5T_NATIVE_SHORT;
826:   else if (ptype == PETSC_ENUM)        *htype = H5T_NATIVE_INT;
827:   else if (ptype == PETSC_BOOL)        *htype = H5T_NATIVE_INT;
828:   else if (ptype == PETSC_FLOAT)       *htype = H5T_NATIVE_FLOAT;
829:   else if (ptype == PETSC_CHAR)        *htype = H5T_NATIVE_CHAR;
830:   else if (ptype == PETSC_BIT_LOGICAL) *htype = H5T_NATIVE_UCHAR;
831:   else if (ptype == PETSC_STRING)      *htype = H5Tcopy(H5T_C_S1);
832:   else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported PETSc datatype");
833:   return(0);
834: }

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

839:   Not collective

841:   Input Parameter:
842: . htype - the HDF5 datatype (for example H5T_NATIVE_DOUBLE, ...)

844:   Output Parameter:
845: . ptype - the PETSc datatype name (for example PETSC_DOUBLE)

847:   Level: advanced

849: .seealso: PetscDataType, PetscHDF5DataTypeToPetscDataType()
850: @*/
851: PetscErrorCode PetscHDF5DataTypeToPetscDataType(hid_t htype, PetscDataType *ptype)
852: {
854: #if defined(PETSC_USE_64BIT_INDICES)
855:   if      (htype == H5T_NATIVE_INT)    *ptype = PETSC_LONG;
856:   else if (htype == H5T_NATIVE_LLONG)  *ptype = PETSC_INT;
857: #else
858:   if      (htype == H5T_NATIVE_INT)    *ptype = PETSC_INT;
859: #endif
860:   else if (htype == H5T_NATIVE_DOUBLE) *ptype = PETSC_DOUBLE;
861:   else if (htype == H5T_NATIVE_LONG)   *ptype = PETSC_LONG;
862:   else if (htype == H5T_NATIVE_SHORT)  *ptype = PETSC_SHORT;
863:   else if (htype == H5T_NATIVE_FLOAT)  *ptype = PETSC_FLOAT;
864:   else if (htype == H5T_NATIVE_CHAR)   *ptype = PETSC_CHAR;
865:   else if (htype == H5T_NATIVE_UCHAR)  *ptype = PETSC_CHAR;
866:   else if (htype == H5T_C_S1)          *ptype = PETSC_STRING;
867:   else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported HDF5 datatype");
868:   return(0);
869: }

871: /*@C
872:  PetscViewerHDF5WriteAttribute - Write an attribute

874:   Input Parameters:
875: + viewer - The HDF5 viewer
876: . parent - The parent dataset/group name
877: . name   - The attribute name
878: . datatype - The attribute type
879: - value    - The attribute value

881:   Level: advanced

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

886: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5WriteObjectAttribute(), PetscViewerHDF5ReadAttribute(), PetscViewerHDF5HasAttribute(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
887: @*/
888: PetscErrorCode PetscViewerHDF5WriteAttribute(PetscViewer viewer, const char parent[], const char name[], PetscDataType datatype, const void *value)
889: {
890:   char           *parentAbsPath;
891:   hid_t          h5, dataspace, obj, attribute, dtype;
892:   PetscBool      has;

901:   PetscViewerHDF5GetAbsolutePath_Internal(viewer, parent, &parentAbsPath);
902:   PetscViewerHDF5Traverse_Internal(viewer, parentAbsPath, PETSC_TRUE, NULL, NULL);
903:   PetscViewerHDF5HasAttribute_Internal(viewer, parentAbsPath, name, &has);
904:   PetscDataTypeToHDF5DataType(datatype, &dtype);
905:   if (datatype == PETSC_STRING) {
906:     size_t len;
907:     PetscStrlen((const char *) value, &len);
908:     PetscStackCallHDF5(H5Tset_size,(dtype, len+1));
909:   }
910:   PetscViewerHDF5GetFileId(viewer, &h5);
911:   PetscStackCallHDF5Return(dataspace,H5Screate,(H5S_SCALAR));
912:   PetscStackCallHDF5Return(obj,H5Oopen,(h5, parentAbsPath, H5P_DEFAULT));
913:   if (has) {
914:     PetscStackCallHDF5Return(attribute,H5Aopen_name,(obj, name));
915:   } else {
916:     PetscStackCallHDF5Return(attribute,H5Acreate2,(obj, name, dtype, dataspace, H5P_DEFAULT, H5P_DEFAULT));
917:   }
918:   PetscStackCallHDF5(H5Awrite,(attribute, dtype, value));
919:   if (datatype == PETSC_STRING) PetscStackCallHDF5(H5Tclose,(dtype));
920:   PetscStackCallHDF5(H5Aclose,(attribute));
921:   PetscStackCallHDF5(H5Oclose,(obj));
922:   PetscStackCallHDF5(H5Sclose,(dataspace));
923:   PetscFree(parentAbsPath);
924:   return(0);
925: }

927: /*@C
928:  PetscViewerHDF5WriteObjectAttribute - Write an attribute to the dataset matching the given PetscObject by name

930:   Input Parameters:
931: + viewer   - The HDF5 viewer
932: . obj      - The object whose name is used to lookup the parent dataset, relative to the current group.
933: . name     - The attribute name
934: . datatype - The attribute type
935: - value    - The attribute value

937:   Notes:
938:   This fails if currentGroup/objectName doesn't resolve to a dataset (the path doesn't exist or is not a dataset).
939:   You might want to check first if it does using PetscViewerHDF5HasObject().

941:   Level: advanced

943: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5WriteAttribute(), PetscViewerHDF5ReadObjectAttribute(), PetscViewerHDF5HasObjectAttribute(), PetscViewerHDF5HasObject(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
944: @*/
945: PetscErrorCode PetscViewerHDF5WriteObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscDataType datatype, const void *value)
946: {

954:   PetscViewerHDF5CheckNamedObject_Internal(viewer, obj);
955:   PetscViewerHDF5WriteAttribute(viewer, obj->name, name, datatype, value);
956:   return(0);
957: }

959: /*@C
960:  PetscViewerHDF5ReadAttribute - Read an attribute

962:   Input Parameters:
963: + viewer - The HDF5 viewer
964: . parent - The parent dataset/group name
965: . name   - The attribute name
966: . datatype - The attribute type
967: - defaultValue - The pointer to the default value

969:   Output Parameter:
970: . value    - The pointer to the read HDF5 attribute value

972:   Notes:
973:   If defaultValue is NULL and the attribute is not found, an error occurs.
974:   If defaultValue is not NULL and the attribute is not found, defaultValue is copied to value.
975:   The pointers defaultValue and value can be the same; for instance
976: $  flg = PETSC_FALSE;
977: $  PetscViewerHDF5ReadAttribute(viewer,name,"attr",PETSC_BOOL,&flg,&flg);
978:   is valid, but make sure the default value is initialized.

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

982:   Level: advanced

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

987: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5ReadObjectAttribute(), PetscViewerHDF5WriteAttribute(), PetscViewerHDF5HasAttribute(), PetscViewerHDF5HasObject(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
988: @*/
989: PetscErrorCode PetscViewerHDF5ReadAttribute(PetscViewer viewer, const char parent[], const char name[], PetscDataType datatype, void *defaultValue, void *value)
990: {
991:   char           *parentAbsPath;
992:   hid_t          h5, obj, attribute, dtype;
993:   PetscBool      has;

1002:   PetscDataTypeToHDF5DataType(datatype, &dtype);
1003:   PetscViewerHDF5GetAbsolutePath_Internal(viewer, parent, &parentAbsPath);
1004:   PetscViewerHDF5Traverse_Internal(viewer, parentAbsPath, PETSC_FALSE, &has, NULL);
1005:   if (has) {PetscViewerHDF5HasAttribute_Internal(viewer, parentAbsPath, name, &has);}
1006:   if (!has) {
1007:     if (defaultValue) {
1008:       if (defaultValue != value) {
1009:         if (datatype == PETSC_STRING) {
1010:           PetscStrallocpy(*(char**)defaultValue, (char**)value);
1011:         } else {
1012:           size_t len;
1013:           PetscStackCallHDF5Return(len,H5Tget_size,(dtype));
1014:           PetscMemcpy(value, defaultValue, len);
1015:         }
1016:       }
1017:       PetscFree(parentAbsPath);
1018:       return(0);
1019:     } else SETERRQ2(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Attribute %s/%s does not exist and default value not provided", parentAbsPath, name);
1020:   }
1021:   PetscViewerHDF5GetFileId(viewer, &h5);
1022:   PetscStackCallHDF5Return(obj,H5Oopen,(h5, parentAbsPath, H5P_DEFAULT));
1023:   PetscStackCallHDF5Return(attribute,H5Aopen_name,(obj, name));
1024:   if (datatype == PETSC_STRING) {
1025:     size_t len;
1026:     hid_t  atype;
1027:     PetscStackCallHDF5Return(atype,H5Aget_type,(attribute));
1028:     PetscStackCallHDF5Return(len,H5Tget_size,(atype));
1029:     PetscMalloc((len+1) * sizeof(char), value);
1030:     PetscStackCallHDF5(H5Tset_size,(dtype, len+1));
1031:     PetscStackCallHDF5(H5Aread,(attribute, dtype, *(char**)value));
1032:   } else {
1033:     PetscStackCallHDF5(H5Aread,(attribute, dtype, value));
1034:   }
1035:   PetscStackCallHDF5(H5Aclose,(attribute));
1036:   /* H5Oclose can be used to close groups, datasets, or committed datatypes */
1037:   PetscStackCallHDF5(H5Oclose,(obj));
1038:   PetscFree(parentAbsPath);
1039:   return(0);
1040: }

1042: /*@C
1043:  PetscViewerHDF5ReadObjectAttribute - Read an attribute from the dataset matching the given PetscObject by name

1045:   Input Parameters:
1046: + viewer   - The HDF5 viewer
1047: . obj      - The object whose name is used to lookup the parent dataset, relative to the current group.
1048: . name     - The attribute name
1049: - datatype - The attribute type

1051:   Output Parameter:
1052: . value    - The attribute value

1054:   Notes:
1055:   This fails if current_group/object_name doesn't resolve to a dataset (the path doesn't exist or is not a dataset).
1056:   You might want to check first if it does using PetscViewerHDF5HasObject().

1058:   Level: advanced

1060: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5ReadAttribute() PetscViewerHDF5WriteObjectAttribute(), PetscViewerHDF5HasObjectAttribute(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
1061: @*/
1062: PetscErrorCode PetscViewerHDF5ReadObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscDataType datatype, void *defaultValue, void *value)
1063: {

1071:   PetscViewerHDF5CheckNamedObject_Internal(viewer, obj);
1072:   PetscViewerHDF5ReadAttribute(viewer, obj->name, name, datatype, defaultValue, value);
1073:   return(0);
1074: }

1076: PETSC_STATIC_INLINE PetscErrorCode PetscViewerHDF5Traverse_Inner_Internal(hid_t h5, const char name[], PetscBool createGroup, PetscBool *exists_)
1077: {
1078:   htri_t exists;
1079:   hid_t group;

1082:   PetscStackCallHDF5Return(exists,H5Lexists,(h5, name, H5P_DEFAULT));
1083:   if (exists) PetscStackCallHDF5Return(exists,H5Oexists_by_name,(h5, name, H5P_DEFAULT));
1084:   if (!exists && createGroup) {
1085:     PetscStackCallHDF5Return(group,H5Gcreate2,(h5, name, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT));
1086:     PetscStackCallHDF5(H5Gclose,(group));
1087:     exists = PETSC_TRUE;
1088:   }
1089:   *exists_ = (PetscBool) exists;
1090:   return(0);
1091: }

1093: static PetscErrorCode PetscViewerHDF5Traverse_Internal(PetscViewer viewer, const char name[], PetscBool createGroup, PetscBool *has, H5O_type_t *otype)
1094: {
1095:   const char     rootGroupName[] = "/";
1096:   hid_t          h5;
1097:   PetscBool      exists=PETSC_FALSE;
1098:   PetscInt       i;
1099:   int            n;
1100:   char           **hierarchy;
1101:   char           buf[PETSC_MAX_PATH_LEN]="";

1107:   else name = rootGroupName;
1108:   if (has) {
1110:     *has = PETSC_FALSE;
1111:   }
1112:   if (otype) {
1114:     *otype = H5O_TYPE_UNKNOWN;
1115:   }
1116:   PetscViewerHDF5GetFileId(viewer, &h5);

1118:   /*
1119:      Unfortunately, H5Oexists_by_name() fails if any object in hierarchy is missing.
1120:      Hence, each of them needs to be tested separately:
1121:      1) whether it's a valid link
1122:      2) whether this link resolves to an object
1123:      See H5Oexists_by_name() documentation.
1124:   */
1125:   PetscStrToArray(name,'/',&n,&hierarchy);
1126:   if (!n) {
1127:     /*  Assume group "/" always exists in accordance with HDF5 >= 1.10.0. See H5Lexists() documentation. */
1128:     if (has)   *has   = PETSC_TRUE;
1129:     if (otype) *otype = H5O_TYPE_GROUP;
1130:     PetscStrToArrayDestroy(n,hierarchy);
1131:     return(0);
1132:   }
1133:   for (i=0; i<n; i++) {
1134:     PetscStrcat(buf,"/");
1135:     PetscStrcat(buf,hierarchy[i]);
1136:     PetscViewerHDF5Traverse_Inner_Internal(h5, buf, createGroup, &exists);
1137:     if (!exists) break;
1138:   }
1139:   PetscStrToArrayDestroy(n,hierarchy);

1141:   /* If the object exists, get its type */
1142:   if (exists && otype) {
1143:     H5O_info_t info;

1145:     /* We could use H5Iget_type() here but that would require opening the object. This way we only need its name. */
1146:     PetscStackCallHDF5(H5Oget_info_by_name,(h5, name, &info, H5P_DEFAULT));
1147:     *otype = info.type;
1148:   }
1149:   if (has) *has = exists;
1150:   return(0);
1151: }

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

1156:   Input Parameters:
1157: + viewer - The HDF5 viewer
1158: - path - The group path

1160:   Output Parameter:
1161: . has - Flag for group existence

1163:   Level: advanced

1165:   Notes:
1166:   If path starts with '/', it is taken as an absolute path overriding currently pushed group, else the dataset/group is relative to the current pushed group. NULL means the current pushed group.
1167:   If the path exists but is not a group, PETSC_FALSE is returned.

1169: .seealso: PetscViewerHDF5HasAttribute(), PetscViewerHDF5HasDataset(), PetscViewerHDF5PushGroup(), PetscViewerHDF5PopGroup(), PetscViewerHDF5GetGroup(), PetscViewerHDF5OpenGroup()
1170: @*/
1171: PetscErrorCode PetscViewerHDF5HasGroup(PetscViewer viewer, const char path[], PetscBool *has)
1172: {
1173:   H5O_type_t     type;
1174:   char           *abspath;

1181:   PetscViewerHDF5GetAbsolutePath_Internal(viewer, path, &abspath);
1182:   PetscViewerHDF5Traverse_Internal(viewer, abspath, PETSC_FALSE, NULL, &type);
1183:   *has = (PetscBool)(type == H5O_TYPE_GROUP);
1184:   PetscFree(abspath);
1185:   return(0);
1186: }

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

1191:   Input Parameters:
1192: + viewer - The HDF5 viewer
1193: - path - The dataset path

1195:   Output Parameter:
1196: . has - Flag whether dataset exists

1198:   Level: advanced

1200:   Notes:
1201:   If parent starts with '/', it is taken as an absolute path overriding currently pushed group, else the dataset/group is relative to the current pushed group.
1202:   If path is NULL or empty, has is set to PETSC_FALSE.
1203:   If the path exists but is not a dataset, has is set to PETSC_FALSE as well.

1205: .seealso: PetscViewerHDF5HasObject(), PetscViewerHDF5HasAttribute(), PetscViewerHDF5HasGroup(), PetscViewerHDF5PushGroup(), PetscViewerHDF5PopGroup(), PetscViewerHDF5GetGroup()
1206: @*/
1207: PetscErrorCode PetscViewerHDF5HasDataset(PetscViewer viewer, const char path[], PetscBool *has)
1208: {
1209:   H5O_type_t     type;
1210:   char           *abspath;

1217:   PetscViewerHDF5GetAbsolutePath_Internal(viewer, path, &abspath);
1218:   PetscViewerHDF5Traverse_Internal(viewer, abspath, PETSC_FALSE, NULL, &type);
1219:   *has = (PetscBool)(type == H5O_TYPE_DATASET);
1220:   PetscFree(abspath);
1221:   return(0);
1222: }

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

1227:   Input Parameters:
1228: + viewer - The HDF5 viewer
1229: - obj    - The named object

1231:   Output Parameter:
1232: . has    - Flag for dataset existence

1234:   Notes:
1235:   If the object is unnamed, an error occurs.
1236:   If the path currentGroup/objectName exists but is not a dataset, has is set to PETSC_FALSE as well.

1238:   Level: advanced

1240: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5HasDataset(), PetscViewerHDF5HasAttribute(), PetscViewerHDF5PushGroup(), PetscViewerHDF5PopGroup(), PetscViewerHDF5GetGroup()
1241: @*/
1242: PetscErrorCode PetscViewerHDF5HasObject(PetscViewer viewer, PetscObject obj, PetscBool *has)
1243: {
1244:   size_t         len;

1251:   PetscStrlen(obj->name, &len);
1252:   if (!len) SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONG, "Object must be named");
1253:   PetscViewerHDF5HasDataset(viewer, obj->name, has);
1254:   return(0);
1255: }

1257: /*@C
1258:  PetscViewerHDF5HasAttribute - Check whether an attribute exists

1260:   Input Parameters:
1261: + viewer - The HDF5 viewer
1262: . parent - The parent dataset/group name
1263: - name   - The attribute name

1265:   Output Parameter:
1266: . has    - Flag for attribute existence

1268:   Level: advanced

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

1273: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5HasObjectAttribute(), PetscViewerHDF5WriteAttribute(), PetscViewerHDF5ReadAttribute(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
1274: @*/
1275: PetscErrorCode PetscViewerHDF5HasAttribute(PetscViewer viewer, const char parent[], const char name[], PetscBool *has)
1276: {
1277:   char           *parentAbsPath;

1285:   PetscViewerHDF5GetAbsolutePath_Internal(viewer, parent, &parentAbsPath);
1286:   PetscViewerHDF5Traverse_Internal(viewer, parentAbsPath, PETSC_FALSE, has, NULL);
1287:   if (*has) {PetscViewerHDF5HasAttribute_Internal(viewer, parentAbsPath, name, has);}
1288:   PetscFree(parentAbsPath);
1289:   return(0);
1290: }

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

1295:   Input Parameters:
1296: + viewer - The HDF5 viewer
1297: . obj    - The object whose name is used to lookup the parent dataset, relative to the current group.
1298: - name   - The attribute name

1300:   Output Parameter:
1301: . has    - Flag for attribute existence

1303:   Notes:
1304:   This fails if current_group/object_name doesn't resolve to a dataset (the path doesn't exist or is not a dataset).
1305:   You might want to check first if it does using PetscViewerHDF5HasObject().

1307:   Level: advanced

1309: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5HasAttribute(), PetscViewerHDF5WriteObjectAttribute(), PetscViewerHDF5ReadObjectAttribute(), PetscViewerHDF5HasObject(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
1310: @*/
1311: PetscErrorCode PetscViewerHDF5HasObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscBool *has)
1312: {

1320:   PetscViewerHDF5CheckNamedObject_Internal(viewer, obj);
1321:   PetscViewerHDF5HasAttribute(viewer, obj->name, name, has);
1322:   return(0);
1323: }

1325: static PetscErrorCode PetscViewerHDF5HasAttribute_Internal(PetscViewer viewer, const char parent[], const char name[], PetscBool *has)
1326: {
1327:   hid_t          h5;
1328:   htri_t         hhas;

1332:   PetscViewerHDF5GetFileId(viewer, &h5);
1333:   PetscStackCallHDF5Return(hhas,H5Aexists_by_name,(h5, parent, name, H5P_DEFAULT));
1334:   *has = hhas ? PETSC_TRUE : PETSC_FALSE;
1335:   return(0);
1336: }

1338: /*
1339:   The variable Petsc_Viewer_HDF5_keyval is used to indicate an MPI attribute that
1340:   is attached to a communicator, in this case the attribute is a PetscViewer.
1341: */
1342: PetscMPIInt Petsc_Viewer_HDF5_keyval = MPI_KEYVAL_INVALID;

1344: /*@C
1345:   PETSC_VIEWER_HDF5_ - Creates an HDF5 PetscViewer shared by all processors in a communicator.

1347:   Collective

1349:   Input Parameter:
1350: . comm - the MPI communicator to share the HDF5 PetscViewer

1352:   Level: intermediate

1354:   Options Database Keys:
1355: . -viewer_hdf5_filename <name>

1357:   Environmental variables:
1358: . PETSC_VIEWER_HDF5_FILENAME

1360:   Notes:
1361:   Unlike almost all other PETSc routines, PETSC_VIEWER_HDF5_ does not return
1362:   an error code.  The HDF5 PetscViewer is usually used in the form
1363: $       XXXView(XXX object, PETSC_VIEWER_HDF5_(comm));

1365: .seealso: PetscViewerHDF5Open(), PetscViewerCreate(), PetscViewerDestroy()
1366: @*/
1367: PetscViewer PETSC_VIEWER_HDF5_(MPI_Comm comm)
1368: {
1370:   PetscBool      flg;
1371:   PetscViewer    viewer;
1372:   char           fname[PETSC_MAX_PATH_LEN];
1373:   MPI_Comm       ncomm;

1376:   PetscCommDuplicate(comm,&ncomm,NULL);if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");PetscFunctionReturn(NULL);}
1377:   if (Petsc_Viewer_HDF5_keyval == MPI_KEYVAL_INVALID) {
1378:     MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN,MPI_COMM_NULL_DELETE_FN,&Petsc_Viewer_HDF5_keyval,NULL);
1379:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");PetscFunctionReturn(NULL);}
1380:   }
1381:   MPI_Comm_get_attr(ncomm,Petsc_Viewer_HDF5_keyval,(void**)&viewer,(int*)&flg);
1382:   if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");PetscFunctionReturn(NULL);}
1383:   if (!flg) { /* PetscViewer not yet created */
1384:     PetscOptionsGetenv(ncomm,"PETSC_VIEWER_HDF5_FILENAME",fname,PETSC_MAX_PATH_LEN,&flg);
1385:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_REPEAT," ");PetscFunctionReturn(NULL);}
1386:     if (!flg) {
1387:       PetscStrcpy(fname,"output.h5");
1388:       if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_REPEAT," ");PetscFunctionReturn(NULL);}
1389:     }
1390:     PetscViewerHDF5Open(ncomm,fname,FILE_MODE_WRITE,&viewer);
1391:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_REPEAT," ");PetscFunctionReturn(NULL);}
1392:     PetscObjectRegisterDestroy((PetscObject)viewer);
1393:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_REPEAT," ");PetscFunctionReturn(NULL);}
1394:     MPI_Comm_set_attr(ncomm,Petsc_Viewer_HDF5_keyval,(void*)viewer);
1395:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");PetscFunctionReturn(NULL);}
1396:   }
1397:   PetscCommDestroy(&ncomm);
1398:   if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_REPEAT," ");PetscFunctionReturn(NULL);}
1399:   PetscFunctionReturn(viewer);
1400: }