Actual source code: hdf5v.c

petsc-3.11.4 2019-09-28
Report Typos and Errors
  1: #include <petsc/private/viewerimpl.h>
  2: #include <petscviewerhdf5.h>    /*I   "petscviewerhdf5.h"   I*/
  3: #if (H5_VERS_MAJOR * 10000 + H5_VERS_MINOR * 100 + H5_VERS_RELEASE < 10800)
  4: #error "PETSc needs HDF5 version >= 1.8.0"
  5: #endif

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

 10: typedef struct GroupList {
 11:   const char       *name;
 12:   struct GroupList *next;
 13: } GroupList;

 15: typedef struct {
 16:   char          *filename;
 17:   PetscFileMode btype;
 18:   hid_t         file_id;
 19:   PetscInt      timestep;
 20:   GroupList     *groups;
 21:   PetscBool     basedimension2;  /* save vectors and DMDA vectors with a dimension of at least 2 even if the bs/dof is 1 */
 22:   PetscBool     spoutput;  /* write data in single precision even if PETSc is compiled with double precision PetscReal */
 23:   char          *mataij_iname;
 24:   char          *mataij_jname;
 25:   char          *mataij_aname;
 26:   char          *mataij_cname;
 27:   PetscBool     mataij_names_set;
 28: } PetscViewer_HDF5;

 30: static PetscErrorCode PetscViewerHDF5GetAbsolutePath_Internal(PetscViewer viewer, const char objname[], char **fullpath)
 31: {
 32:   const char *group;
 33:   char buf[PETSC_MAX_PATH_LEN]="";

 37:   PetscViewerHDF5GetGroup(viewer, &group);
 38:   PetscStrcat(buf, group);
 39:   PetscStrcat(buf, "/");
 40:   PetscStrcat(buf, objname);
 41:   PetscStrallocpy(buf, fullpath);
 42:   return(0);
 43: }

 45: static PetscErrorCode PetscViewerHDF5CheckNamedObject_Internal(PetscViewer viewer, PetscObject obj)
 46: {
 47:   PetscBool has;
 48:   const char *group;

 52:   if (!obj->name) SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONG, "Object must be named");
 53:   PetscViewerHDF5HasObject(viewer, obj, &has);
 54:   if (!has) {
 55:     PetscViewerHDF5GetGroup(viewer, &group);
 56:     SETERRQ2(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Object (dataset) %s not stored in group %s", obj->name, group);
 57:   }
 58:   return(0);
 59: }

 61: static PetscErrorCode PetscViewerSetFromOptions_HDF5(PetscOptionItems *PetscOptionsObject,PetscViewer v)
 62: {
 63:   PetscErrorCode   ierr;
 64:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*)v->data;

 67:   PetscOptionsHead(PetscOptionsObject,"HDF5 PetscViewer Options");
 68:   PetscOptionsBool("-viewer_hdf5_base_dimension2","1d Vectors get 2 dimensions in HDF5","PetscViewerHDF5SetBaseDimension2",hdf5->basedimension2,&hdf5->basedimension2,NULL);
 69:   PetscOptionsBool("-viewer_hdf5_sp_output","Force data to be written in single precision","PetscViewerHDF5SetSPOutput",hdf5->spoutput,&hdf5->spoutput,NULL);
 70:   PetscOptionsTail();
 71:   return(0);
 72: }

 74: static PetscErrorCode PetscViewerFileClose_HDF5(PetscViewer viewer)
 75: {
 76:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*)viewer->data;
 77:   PetscErrorCode   ierr;

 80:   PetscFree(hdf5->filename);
 81:   if (hdf5->file_id) PetscStackCallHDF5(H5Fclose,(hdf5->file_id));
 82:   return(0);
 83: }

 85: static PetscErrorCode PetscViewerDestroy_HDF5(PetscViewer viewer)
 86: {
 87:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
 88:   PetscErrorCode   ierr;

 91:   PetscViewerFileClose_HDF5(viewer);
 92:   while (hdf5->groups) {
 93:     GroupList *tmp = hdf5->groups->next;

 95:     PetscFree(hdf5->groups->name);
 96:     PetscFree(hdf5->groups);
 97:     hdf5->groups = tmp;
 98:   }
 99:   PetscFree(hdf5->mataij_iname);
100:   PetscFree(hdf5->mataij_jname);
101:   PetscFree(hdf5->mataij_aname);
102:   PetscFree(hdf5->mataij_cname);
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:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerHDF5SetAIJNames_C",NULL);
110:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerHDF5GetAIJNames_C",NULL);
111:   return(0);
112: }

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

119:   hdf5->btype = type;
120:   return(0);
121: }

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

128:   *type = hdf5->btype;
129:   return(0);
130: }

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

137:   hdf5->basedimension2 = flg;
138:   return(0);
139: }

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

145:     Logically Collective on PetscViewer

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

151:   Options Database:
152: .  -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


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

159:   Level: intermediate

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

163: @*/
164: PetscErrorCode PetscViewerHDF5SetBaseDimension2(PetscViewer viewer,PetscBool flg)
165: {

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

174: /*@
175:      PetscViewerHDF5GetBaseDimension2 - 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 on PetscViewer

180:   Input Parameter:
181: .  viewer - the PetscViewer, must be of type HDF5

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

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

190:   Level: intermediate

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

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

201:   *flg = hdf5->basedimension2;
202:   return(0);
203: }

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

210:   hdf5->spoutput = flg;
211:   return(0);
212: }

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

218:     Logically Collective on PetscViewer

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

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


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

232:   Level: intermediate

234: .seealso: PetscViewerFileSetMode(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerBinaryOpen(),
235:           PetscReal

237: @*/
238: PetscErrorCode PetscViewerHDF5SetSPOutput(PetscViewer viewer,PetscBool flg)
239: {

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

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

252:     Logically Collective on PetscViewer

254:   Input Parameter:
255: .  viewer - the PetscViewer, must be of type HDF5

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

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

264:   Level: intermediate

266: .seealso: PetscViewerFileSetMode(), PetscViewerCreate(), PetscViewerSetType(), PetscViewerBinaryOpen(),
267:           PetscReal

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

276:   *flg = hdf5->spoutput;
277:   return(0);
278: }

280: static PetscErrorCode  PetscViewerFileSetName_HDF5(PetscViewer viewer, const char name[])
281: {
282:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
283: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
284:   MPI_Info          info = MPI_INFO_NULL;
285: #endif
286:   hid_t             plist_id;
287:   PetscErrorCode    ierr;

290:   if (hdf5->file_id) PetscStackCallHDF5(H5Fclose,(hdf5->file_id));
291:   if (hdf5->filename) {PetscFree(hdf5->filename);}
292:   PetscStrallocpy(name, &hdf5->filename);
293:   /* Set up file access property list with parallel I/O access */
294:   PetscStackCallHDF5Return(plist_id,H5Pcreate,(H5P_FILE_ACCESS));
295: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
296:   PetscStackCallHDF5(H5Pset_fapl_mpio,(plist_id, PetscObjectComm((PetscObject)viewer), info));
297: #endif
298:   /* Create or open the file collectively */
299:   switch (hdf5->btype) {
300:   case FILE_MODE_READ:
301:     PetscStackCallHDF5Return(hdf5->file_id,H5Fopen,(name, H5F_ACC_RDONLY, plist_id));
302:     break;
303:   case FILE_MODE_APPEND:
304:     PetscStackCallHDF5Return(hdf5->file_id,H5Fopen,(name, H5F_ACC_RDWR, plist_id));
305:     break;
306:   case FILE_MODE_WRITE:
307:     PetscStackCallHDF5Return(hdf5->file_id,H5Fcreate,(name, H5F_ACC_TRUNC, H5P_DEFAULT, plist_id));
308:     break;
309:   default:
310:     SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
311:   }
312:   if (hdf5->file_id < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB, "H5Fcreate failed for %s", name);
313:   PetscStackCallHDF5(H5Pclose,(plist_id));
314:   return(0);
315: }

317: static PetscErrorCode PetscViewerFileGetName_HDF5(PetscViewer viewer,const char **name)
318: {
319:   PetscViewer_HDF5 *vhdf5 = (PetscViewer_HDF5*)viewer->data;

322:   *name = vhdf5->filename;
323:   return(0);
324: }

326: static PetscErrorCode  PetscViewerHDF5SetAIJNames_HDF5(PetscViewer viewer, const char iname[], const char jname[], const char aname[], const char cname[])
327: {
328:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

332:   PetscFree(hdf5->mataij_iname);
333:   PetscFree(hdf5->mataij_jname);
334:   PetscFree(hdf5->mataij_aname);
335:   PetscFree(hdf5->mataij_cname);
336:   PetscStrallocpy(iname,&hdf5->mataij_iname);
337:   PetscStrallocpy(jname,&hdf5->mataij_jname);
338:   PetscStrallocpy(aname,&hdf5->mataij_aname);
339:   PetscStrallocpy(cname,&hdf5->mataij_cname);
340:   hdf5->mataij_names_set = PETSC_TRUE;
341:   return(0);
342: }

344: /*@C
345:   PetscViewerHDF5SetAIJNames - Set the names of the datasets representing the three AIJ (CRS) arrays and the name of the attribute storing the number of columns within the HDF5 file.

347:   Collective on PetscViewer

349:   Input Parameters:
350: +  viewer - the PetscViewer; either ASCII or binary
351: .  iname - name of dataset i representing row pointers; that is i[0] = 0, i[row] = i[row-1] + number of elements in that row of the matrix
352: .  jname - name of dataset j representing column indices
353: .  aname - name of dataset a representing matrix values
354: -  cname - name of attribute stoting column count

356:   Level: advanced

358:   Notes:
359:   Current defaults are (iname, jname, aname, cname) = ("jc", "ir", "data", "MATLAB_sparse") so that MAT files can be readily loaded.

361: .seealso: MatLoad(), PetscViewerCreate(), PetscViewerSetType(), PETSCVIEWERHDF5, PetscViewerHDF5GetAIJNames()
362: @*/
363: /* TODO Once corresponding MatView is implemented, add this:
364:   Current defaults are (iname, jname, aname, cname) = ("i", "j", "a", "ncols").
365:   For PetscViewerFormat PETSC_VIEWER_HDF5_MAT they are ("jc", "ir", "data", "MATLAB_sparse") so that MAT files can be loaded.
366: */
367: PetscErrorCode  PetscViewerHDF5SetAIJNames(PetscViewer viewer, const char iname[], const char jname[], const char aname[], const char cname[])
368: {

377:   PetscTryMethod(viewer,"PetscViewerHDF5SetAIJNames_C",(PetscViewer,const char[],const char[],const char[],const char[]),(viewer,iname,jname,aname,cname));
378:   return(0);
379: }

381: static PetscErrorCode  PetscViewerHDF5GetAIJNames_HDF5(PetscViewer viewer, const char *iname[], const char *jname[], const char *aname[], const char *cname[])
382: {
383:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

386:   *iname = hdf5->mataij_iname;
387:   *jname = hdf5->mataij_jname;
388:   *aname = hdf5->mataij_aname;
389:   *cname = hdf5->mataij_cname;
390:   return(0);
391: }

393: /*@C
394:   PetscViewerHDF5GetAIJNames - Get the names of the datasets representing the three AIJ (CRS) arrays and the name of the attribute storing the number of columns within the HDF5 file.

396:   Collective on PetscViewer

398:   Input Parameters:
399: .  viewer - the PetscViewer; either ASCII or binary

401:   Output Parameters:
402: +  iname - name of dataset i representing row pointers; that is i[0] = 0, i[row] = i[row-1] + number of elements in that row of the matrix
403: .  jname - name of dataset j representing column indices
404: .  aname - name of dataset a representing matrix values
405: -  cname - name of attribute stoting column count

407:   Level: advanced

409:   Notes:
410:   Current defaults are (iname, jname, aname, cname) = ("jc", "ir", "data", "MATLAB_sparse") so that MAT files can be readily loaded.

412: .seealso: MatLoad(), PetscViewerCreate(), PetscViewerSetType(), PETSCVIEWERHDF5, PetscViewerHDF5SetAIJNames()
413: @*/
414: /* TODO Once corresponding MatView is implemented, add this:
415:   Current defaults are (iname, jname, aname, cname) = ("i", "j", "a", "ncols").
416:   For PetscViewerFormat PETSC_VIEWER_HDF5_MAT they are ("jc", "ir", "data", "MATLAB_sparse") so that MAT files can be loaded.
417: */
418: PetscErrorCode  PetscViewerHDF5GetAIJNames(PetscViewer viewer, const char *iname[], const char *jname[], const char *aname[], const char *cname[])
419: {

428:   PetscUseMethod(viewer,"PetscViewerHDF5GetAIJNames_C",(PetscViewer,const char*[],const char*[],const char*[],const char*[]),(viewer,iname,jname,aname,cname));
429:   return(0);
430: }

432: static PetscErrorCode PetscViewerSetUp_HDF5(PetscViewer viewer)
433: {
434:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
435:   PetscErrorCode   ierr;

438:   if (!hdf5->mataij_names_set) {
439: /* TODO Once corresponding MatView is implemented, uncomment */
440: #if 0
441:     if (viewer->format == PETSC_VIEWER_HDF5_MAT) {
442: #endif
443:       PetscViewerHDF5SetAIJNames_HDF5(viewer,"jc","ir","data","MATLAB_sparse");
444: #if 0
445:     } else {
446:       PetscViewerHDF5SetAIJNames_HDF5(viewer,"i","j","a","ncols");
447:     }
448: #endif
449:   }
450:   return(0);
451: }

453: /*MC
454:    PETSCVIEWERHDF5 - A viewer that writes to an HDF5 file


457: .seealso:  PetscViewerHDF5Open(), PetscViewerStringSPrintf(), PetscViewerSocketOpen(), PetscViewerDrawOpen(), PETSCVIEWERSOCKET,
458:            PetscViewerCreate(), PetscViewerASCIIOpen(), PetscViewerBinaryOpen(), PETSCVIEWERBINARY, PETSCVIEWERDRAW, PETSCVIEWERSTRING,
459:            PetscViewerMatlabOpen(), VecView(), DMView(), PetscViewerMatlabPutArray(), PETSCVIEWERASCII, PETSCVIEWERMATLAB,
460:            PetscViewerFileSetName(), PetscViewerFileSetMode(), PetscViewerFormat, PetscViewerType, PetscViewerSetType()

462:   Level: beginner
463: M*/

465: PETSC_EXTERN PetscErrorCode PetscViewerCreate_HDF5(PetscViewer v)
466: {
467:   PetscViewer_HDF5 *hdf5;
468:   PetscErrorCode   ierr;

471:   PetscNewLog(v,&hdf5);

473:   v->data                = (void*) hdf5;
474:   v->ops->destroy        = PetscViewerDestroy_HDF5;
475:   v->ops->setfromoptions = PetscViewerSetFromOptions_HDF5;
476:   v->ops->setup          = PetscViewerSetUp_HDF5;
477:   v->ops->flush          = 0;
478:   hdf5->btype            = (PetscFileMode) -1;
479:   hdf5->filename         = 0;
480:   hdf5->timestep         = -1;
481:   hdf5->groups           = NULL;

483:   hdf5->mataij_iname     = NULL;
484:   hdf5->mataij_jname     = NULL;
485:   hdf5->mataij_aname     = NULL;
486:   hdf5->mataij_cname     = NULL;
487:   hdf5->mataij_names_set = PETSC_FALSE;

489:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileSetName_C",PetscViewerFileSetName_HDF5);
490:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileGetName_C",PetscViewerFileGetName_HDF5);
491:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileSetMode_C",PetscViewerFileSetMode_HDF5);
492:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileGetMode_C",PetscViewerFileGetMode_HDF5);
493:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerHDF5SetBaseDimension2_C",PetscViewerHDF5SetBaseDimension2_HDF5);
494:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerHDF5SetSPOutput_C",PetscViewerHDF5SetSPOutput_HDF5);
495:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerHDF5SetAIJNames_C",PetscViewerHDF5SetAIJNames_HDF5);
496:   PetscObjectComposeFunction((PetscObject)v,"PetscViewerHDF5GetAIJNames_C",PetscViewerHDF5GetAIJNames_HDF5);
497:   return(0);
498: }

500: /*@C
501:    PetscViewerHDF5Open - Opens a file for HDF5 input/output.

503:    Collective on MPI_Comm

505:    Input Parameters:
506: +  comm - MPI communicator
507: .  name - name of file
508: -  type - type of file
509: $    FILE_MODE_WRITE - create new file for binary output
510: $    FILE_MODE_READ - open existing file for binary input
511: $    FILE_MODE_APPEND - open existing file for binary output

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

516:   Options Database:
517: .  -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
518: .  -viewer_hdf5_sp_output - forces (if true) the viewer to write data in single precision independent on the precision of PetscReal

520:    Level: beginner

522:    Note:
523:    This PetscViewer should be destroyed with PetscViewerDestroy().

525:    Concepts: HDF5 files
526:    Concepts: PetscViewerHDF5^creating

528: .seealso: PetscViewerASCIIOpen(), PetscViewerPushFormat(), PetscViewerDestroy(), PetscViewerHDF5SetBaseDimension2(),
529:           PetscViewerHDF5SetSPOutput(), PetscViewerHDF5GetBaseDimension2(), VecView(), MatView(), VecLoad(),
530:           MatLoad(), PetscFileMode, PetscViewer, PetscViewerSetType(), PetscViewerFileSetMode(), PetscViewerFileSetName()
531: @*/
532: PetscErrorCode  PetscViewerHDF5Open(MPI_Comm comm, const char name[], PetscFileMode type, PetscViewer *hdf5v)
533: {

537:   PetscViewerCreate(comm, hdf5v);
538:   PetscViewerSetType(*hdf5v, PETSCVIEWERHDF5);
539:   PetscViewerFileSetMode(*hdf5v, type);
540:   PetscViewerFileSetName(*hdf5v, name);
541:   return(0);
542: }

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

547:   Not collective

549:   Input Parameter:
550: . viewer - the PetscViewer

552:   Output Parameter:
553: . file_id - The file id

555:   Level: intermediate

557: .seealso: PetscViewerHDF5Open()
558: @*/
559: PetscErrorCode  PetscViewerHDF5GetFileId(PetscViewer viewer, hid_t *file_id)
560: {
561:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

565:   if (file_id) *file_id = hdf5->file_id;
566:   return(0);
567: }

569: /*@C
570:   PetscViewerHDF5PushGroup - Set the current HDF5 group for output

572:   Not collective

574:   Input Parameters:
575: + viewer - the PetscViewer
576: - name - The group name

578:   Level: intermediate

580: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup(),PetscViewerHDF5OpenGroup()
581: @*/
582: PetscErrorCode  PetscViewerHDF5PushGroup(PetscViewer viewer, const char *name)
583: {
584:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
585:   GroupList        *groupNode;
586:   PetscErrorCode   ierr;

591:   PetscNew(&groupNode);
592:   PetscStrallocpy(name, (char**) &groupNode->name);

594:   groupNode->next = hdf5->groups;
595:   hdf5->groups    = groupNode;
596:   return(0);
597: }

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

602:   Not collective

604:   Input Parameter:
605: . viewer - the PetscViewer

607:   Level: intermediate

609: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5GetGroup(),PetscViewerHDF5OpenGroup()
610: @*/
611: PetscErrorCode  PetscViewerHDF5PopGroup(PetscViewer viewer)
612: {
613:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;
614:   GroupList        *groupNode;
615:   PetscErrorCode   ierr;

619:   if (!hdf5->groups) SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "HDF5 group stack is empty, cannot pop");
620:   groupNode    = hdf5->groups;
621:   hdf5->groups = hdf5->groups->next;
622:   PetscFree(groupNode->name);
623:   PetscFree(groupNode);
624:   return(0);
625: }

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

631:   Not collective

633:   Input Parameter:
634: . viewer - the PetscViewer

636:   Output Parameter:
637: . name - The group name

639:   Level: intermediate

641: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5OpenGroup()
642: @*/
643: PetscErrorCode  PetscViewerHDF5GetGroup(PetscViewer viewer, const char **name)
644: {
645:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;

650:   if (hdf5->groups) *name = hdf5->groups->name;
651:   else *name = NULL;
652:   return(0);
653: }

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

660:   Not collective

662:   Input Parameter:
663: . viewer - the PetscViewer

665:   Output Parameter:
666: + fileId - The HDF5 file ID
667: - groupId - The HDF5 group ID

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

672:   Level: intermediate

674: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
675: @*/
676: PetscErrorCode PetscViewerHDF5OpenGroup(PetscViewer viewer, hid_t *fileId, hid_t *groupId)
677: {
678:   hid_t          file_id;
679:   H5O_type_t     type;
680:   const char     *groupName = NULL;
681:   PetscBool      create;

685:   PetscViewerWritable(viewer, &create);
686:   PetscViewerHDF5GetFileId(viewer, &file_id);
687:   PetscViewerHDF5GetGroup(viewer, &groupName);
688:   PetscViewerHDF5Traverse_Internal(viewer, groupName, create, NULL, &type);
689:   if (type != H5O_TYPE_GROUP) SETERRQ1(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Path %s resolves to something which is not a group", groupName);
690:   PetscStackCallHDF5Return(*groupId,H5Gopen2,(file_id, groupName ? groupName : "/", H5P_DEFAULT));
691:   *fileId  = file_id;
692:   return(0);
693: }

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

698:   Not collective

700:   Input Parameter:
701: . viewer - the PetscViewer

703:   Level: intermediate

705: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5SetTimestep(), PetscViewerHDF5GetTimestep()
706: @*/
707: PetscErrorCode PetscViewerHDF5IncrementTimestep(PetscViewer viewer)
708: {
709:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

713:   ++hdf5->timestep;
714:   return(0);
715: }

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

721:   Not collective

723:   Input Parameters:
724: + viewer - the PetscViewer
725: - timestep - The timestep number

727:   Level: intermediate

729: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5GetTimestep()
730: @*/
731: PetscErrorCode  PetscViewerHDF5SetTimestep(PetscViewer viewer, PetscInt timestep)
732: {
733:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

737:   hdf5->timestep = timestep;
738:   return(0);
739: }

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

744:   Not collective

746:   Input Parameter:
747: . viewer - the PetscViewer

749:   Output Parameter:
750: . timestep - The timestep number

752:   Level: intermediate

754: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5SetTimestep()
755: @*/
756: PetscErrorCode  PetscViewerHDF5GetTimestep(PetscViewer viewer, PetscInt *timestep)
757: {
758:   PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data;

763:   *timestep = hdf5->timestep;
764:   return(0);
765: }

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

770:   Not collective

772:   Input Parameter:
773: . ptype - the PETSc datatype name (for example PETSC_DOUBLE)

775:   Output Parameter:
776: . mtype - the MPI datatype (for example MPI_DOUBLE, ...)

778:   Level: advanced

780: .seealso: PetscDataType, PetscHDF5DataTypeToPetscDataType()
781: @*/
782: PetscErrorCode PetscDataTypeToHDF5DataType(PetscDataType ptype, hid_t *htype)
783: {
785:   if (ptype == PETSC_INT)
786: #if defined(PETSC_USE_64BIT_INDICES)
787:                                        *htype = H5T_NATIVE_LLONG;
788: #else
789:                                        *htype = H5T_NATIVE_INT;
790: #endif
791:   else if (ptype == PETSC_DOUBLE)      *htype = H5T_NATIVE_DOUBLE;
792:   else if (ptype == PETSC_LONG)        *htype = H5T_NATIVE_LONG;
793:   else if (ptype == PETSC_SHORT)       *htype = H5T_NATIVE_SHORT;
794:   else if (ptype == PETSC_ENUM)        *htype = H5T_NATIVE_INT;
795:   else if (ptype == PETSC_BOOL)        *htype = H5T_NATIVE_INT;
796:   else if (ptype == PETSC_FLOAT)       *htype = H5T_NATIVE_FLOAT;
797:   else if (ptype == PETSC_CHAR)        *htype = H5T_NATIVE_CHAR;
798:   else if (ptype == PETSC_BIT_LOGICAL) *htype = H5T_NATIVE_UCHAR;
799:   else if (ptype == PETSC_STRING)      *htype = H5Tcopy(H5T_C_S1);
800:   else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported PETSc datatype");
801:   return(0);
802: }

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

807:   Not collective

809:   Input Parameter:
810: . htype - the HDF5 datatype (for example H5T_NATIVE_DOUBLE, ...)

812:   Output Parameter:
813: . ptype - the PETSc datatype name (for example PETSC_DOUBLE)

815:   Level: advanced

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

839: /*@C
840:  PetscViewerHDF5WriteAttribute - Write an attribute

842:   Input Parameters:
843: + viewer - The HDF5 viewer
844: . dataset - The parent dataset name, relative to the current group. NULL means a group-wise attribute.
845: . name   - The attribute name
846: . datatype - The attribute type
847: - value    - The attribute value

849:   Level: advanced

851: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5WriteObjectAttribute(), PetscViewerHDF5ReadAttribute(), PetscViewerHDF5HasAttribute(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
852: @*/
853: PetscErrorCode PetscViewerHDF5WriteAttribute(PetscViewer viewer, const char dataset[], const char name[], PetscDataType datatype, const void *value)
854: {
855:   char           *parent;
856:   hid_t          h5, dataspace, obj, attribute, dtype;
857:   PetscBool      has;

865:   PetscViewerHDF5GetAbsolutePath_Internal(viewer, dataset, &parent);
866:   PetscViewerHDF5Traverse_Internal(viewer, parent, PETSC_TRUE, NULL, NULL);
867:   PetscViewerHDF5HasAttribute_Internal(viewer, parent, name, &has);
868:   PetscDataTypeToHDF5DataType(datatype, &dtype);
869:   if (datatype == PETSC_STRING) {
870:     size_t len;
871:     PetscStrlen((const char *) value, &len);
872:     PetscStackCallHDF5(H5Tset_size,(dtype, len+1));
873:   }
874:   PetscViewerHDF5GetFileId(viewer, &h5);
875:   PetscStackCallHDF5Return(dataspace,H5Screate,(H5S_SCALAR));
876:   PetscStackCallHDF5Return(obj,H5Oopen,(h5, parent, H5P_DEFAULT));
877:   if (has) {
878:     PetscStackCallHDF5Return(attribute,H5Aopen_name,(obj, name));
879:   } else {
880:     PetscStackCallHDF5Return(attribute,H5Acreate2,(obj, name, dtype, dataspace, H5P_DEFAULT, H5P_DEFAULT));
881:   }
882:   PetscStackCallHDF5(H5Awrite,(attribute, dtype, value));
883:   if (datatype == PETSC_STRING) PetscStackCallHDF5(H5Tclose,(dtype));
884:   PetscStackCallHDF5(H5Aclose,(attribute));
885:   PetscStackCallHDF5(H5Oclose,(obj));
886:   PetscStackCallHDF5(H5Sclose,(dataspace));
887:   PetscFree(parent);
888:   return(0);
889: }

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

894:   Input Parameters:
895: + viewer   - The HDF5 viewer
896: . obj      - The object whose name is used to lookup the parent dataset, relative to the current group.
897: . name     - The attribute name
898: . datatype - The attribute type
899: - value    - The attribute value

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

905:   Level: advanced

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

918:   PetscViewerHDF5CheckNamedObject_Internal(viewer, obj);
919:   PetscViewerHDF5WriteAttribute(viewer, obj->name, name, datatype, value);
920:   return(0);
921: }

923: /*@C
924:  PetscViewerHDF5ReadAttribute - Read an attribute

926:   Input Parameters:
927: + viewer - The HDF5 viewer
928: . dataset - The parent dataset name, relative to the current group. NULL means a group-wise attribute.
929: . name   - The attribute name
930: - datatype - The attribute type

932:   Output Parameter:
933: . value    - The attribute value

935:   Level: advanced

937: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5ReadObjectAttribute(), PetscViewerHDF5WriteAttribute(), PetscViewerHDF5HasAttribute(), PetscViewerHDF5HasObject(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
938: @*/
939: PetscErrorCode PetscViewerHDF5ReadAttribute(PetscViewer viewer, const char dataset[], const char name[], PetscDataType datatype, void *value)
940: {
941:   char           *parent;
942:   hid_t          h5, obj, attribute, atype, dtype;
943:   PetscBool      has;

951:   PetscViewerHDF5GetAbsolutePath_Internal(viewer, dataset, &parent);
952:   PetscViewerHDF5Traverse_Internal(viewer, parent, PETSC_FALSE, &has, NULL);
953:   if (has) {PetscViewerHDF5HasAttribute_Internal(viewer, parent, name, &has);}
954:   if (!has) SETERRQ2(PetscObjectComm((PetscObject)viewer), PETSC_ERR_FILE_UNEXPECTED, "Attribute %s/%s does not exist", parent, name);
955:   PetscDataTypeToHDF5DataType(datatype, &dtype);
956:   PetscViewerHDF5GetFileId(viewer, &h5);
957:   PetscStackCallHDF5Return(obj,H5Oopen,(h5, parent, H5P_DEFAULT));
958:   PetscStackCallHDF5Return(attribute,H5Aopen_name,(obj, name));
959:   if (datatype == PETSC_STRING) {
960:     size_t len;
961:     PetscStackCallHDF5Return(atype,H5Aget_type,(attribute));
962:     PetscStackCall("H5Tget_size",len = H5Tget_size(atype));
963:     PetscStackCallHDF5(H5Tclose,(atype));
964:     PetscMalloc((len+1) * sizeof(char *), &value);
965:   }
966:   PetscStackCallHDF5(H5Aread,(attribute, dtype, value));
967:   PetscStackCallHDF5(H5Aclose,(attribute));
968:   /* H5Oclose can be used to close groups, datasets, or committed datatypes */
969:   PetscStackCallHDF5(H5Oclose,(obj));
970:   PetscFree(parent);
971:   return(0);
972: }

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

977:   Input Parameters:
978: + viewer   - The HDF5 viewer
979: . obj      - The object whose name is used to lookup the parent dataset, relative to the current group.
980: . name     - The attribute name
981: - datatype - The attribute type

983:   Output Parameter:
984: . value    - The attribute value

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

990:   Level: advanced

992: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5ReadAttribute() PetscViewerHDF5WriteObjectAttribute(), PetscViewerHDF5HasObjectAttribute(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
993: @*/
994: PetscErrorCode PetscViewerHDF5ReadObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscDataType datatype, void *value)
995: {

1003:   PetscViewerHDF5CheckNamedObject_Internal(viewer, obj);
1004:   PetscViewerHDF5ReadAttribute(viewer, obj->name, name, datatype, value);
1005:   return(0);
1006: }

1008: PETSC_STATIC_INLINE PetscErrorCode PetscViewerHDF5Traverse_Inner_Internal(hid_t h5, const char name[], PetscBool createGroup, PetscBool *exists_)
1009: {
1010:   htri_t exists;
1011:   hid_t group;

1014:   PetscStackCallHDF5Return(exists,H5Lexists,(h5, name, H5P_DEFAULT));
1015:   if (exists) PetscStackCallHDF5Return(exists,H5Oexists_by_name,(h5, name, H5P_DEFAULT));
1016:   if (!exists && createGroup) {
1017:     PetscStackCallHDF5Return(group,H5Gcreate2,(h5, name, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT));
1018:     PetscStackCallHDF5(H5Gclose,(group));
1019:     exists = PETSC_TRUE;
1020:   }
1021:   *exists_ = (PetscBool) exists;
1022:   return(0);
1023: }

1025: static PetscErrorCode PetscViewerHDF5Traverse_Internal(PetscViewer viewer, const char name[], PetscBool createGroup, PetscBool *has, H5O_type_t *otype)
1026: {
1027:   const char     rootGroupName[] = "/";
1028:   hid_t          h5;
1029:   PetscBool      exists=PETSC_FALSE;
1030:   PetscInt       i;
1031:   int            n;
1032:   char           **hierarchy;
1033:   char           buf[PETSC_MAX_PATH_LEN]="";

1039:   else name = rootGroupName;
1040:   if (has) {
1042:     *has = PETSC_FALSE;
1043:   }
1044:   if (otype) {
1046:     *otype = H5O_TYPE_UNKNOWN;
1047:   }
1048:   PetscViewerHDF5GetFileId(viewer, &h5);

1050:   /*
1051:      Unfortunately, H5Oexists_by_name() fails if any object in hierarchy is missing.
1052:      Hence, each of them needs to be tested separately:
1053:      1) whether it's a valid link
1054:      2) whether this link resolves to an object
1055:      See H5Oexists_by_name() documentation.
1056:   */
1057:   PetscStrToArray(name,'/',&n,&hierarchy);
1058:   if (!n) {
1059:     /*  Assume group "/" always exists in accordance with HDF5 >= 1.10.0. See H5Lexists() documentation. */
1060:     if (has)   *has   = PETSC_TRUE;
1061:     if (otype) *otype = H5O_TYPE_GROUP;
1062:     PetscStrToArrayDestroy(n,hierarchy);
1063:     return(0);
1064:   }
1065:   for (i=0; i<n; i++) {
1066:     PetscStrcat(buf,"/");
1067:     PetscStrcat(buf,hierarchy[i]);
1068:     PetscViewerHDF5Traverse_Inner_Internal(h5, buf, createGroup, &exists);
1069:     if (!exists) break;
1070:   }
1071:   PetscStrToArrayDestroy(n,hierarchy);

1073:   /* If the object exists, get its type */
1074:   if (exists && otype) {
1075:     H5O_info_t info;

1077:     /* We could use H5Iget_type() here but that would require opening the object. This way we only need its name. */
1078:     PetscStackCallHDF5(H5Oget_info_by_name,(h5, name, &info, H5P_DEFAULT));
1079:     *otype = info.type;
1080:   }
1081:   if (has) *has = exists;
1082:   return(0);
1083: }

1085: /*@
1086:  PetscViewerHDF5HasGroup - Check whether the current (pushed) group exists in the HDF5 file

1088:   Input Parameters:
1089: . viewer - The HDF5 viewer

1091:   Output Parameter:
1092: . has    - Flag for group existence

1094:   Notes:
1095:   If the path exists but is not a group, this returns PETSC_FALSE as well.

1097:   Level: advanced

1099: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5PushGroup(), PetscViewerHDF5PopGroup(), PetscViewerHDF5OpenGroup()
1100: @*/
1101: PetscErrorCode PetscViewerHDF5HasGroup(PetscViewer viewer, PetscBool *has)
1102: {
1103:   H5O_type_t type;
1104:   const char *name;

1110:   PetscViewerHDF5GetGroup(viewer, &name);
1111:   PetscViewerHDF5Traverse_Internal(viewer, name, PETSC_FALSE, has, &type);
1112:   *has = (type == H5O_TYPE_GROUP) ? PETSC_TRUE : PETSC_FALSE;
1113:   return(0);
1114: }

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

1119:   Input Parameters:
1120: + viewer - The HDF5 viewer
1121: - obj    - The named object

1123:   Output Parameter:
1124: . has    - Flag for dataset existence; PETSC_FALSE for unnamed object

1126:   Notes:
1127:   If the path exists but is not a dataset, this returns PETSC_FALSE as well.

1129:   Level: advanced

1131: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5HasAttribute(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
1132: @*/
1133: PetscErrorCode PetscViewerHDF5HasObject(PetscViewer viewer, PetscObject obj, PetscBool *has)
1134: {
1135:   H5O_type_t type;
1136:   char *path;

1143:   *has = PETSC_FALSE;
1144:   if (!obj->name) return(0);
1145:   PetscViewerHDF5GetAbsolutePath_Internal(viewer, obj->name, &path);
1146:   PetscViewerHDF5Traverse_Internal(viewer, path, PETSC_FALSE, has, &type);
1147:   *has = (type == H5O_TYPE_DATASET) ? PETSC_TRUE : PETSC_FALSE;
1148:   PetscFree(path);
1149:   return(0);
1150: }

1152: /*@C
1153:  PetscViewerHDF5HasAttribute - Check whether an attribute exists

1155:   Input Parameters:
1156: + viewer - The HDF5 viewer
1157: . dataset - The parent dataset name, relative to the current group. NULL means a group-wise attribute.
1158: - name   - The attribute name

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

1163:   Level: advanced

1165: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5HasObjectAttribute(), PetscViewerHDF5WriteAttribute(), PetscViewerHDF5ReadAttribute(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
1166: @*/
1167: PetscErrorCode PetscViewerHDF5HasAttribute(PetscViewer viewer, const char dataset[], const char name[], PetscBool *has)
1168: {
1169:   char           *parent;

1177:   PetscViewerHDF5GetAbsolutePath_Internal(viewer, dataset, &parent);
1178:   PetscViewerHDF5Traverse_Internal(viewer, parent, PETSC_FALSE, has, NULL);
1179:   if (*has) {PetscViewerHDF5HasAttribute_Internal(viewer, parent, name, has);}
1180:   PetscFree(parent);
1181:   return(0);
1182: }

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

1187:   Input Parameters:
1188: + viewer - The HDF5 viewer
1189: . obj    - The object whose name is used to lookup the parent dataset, relative to the current group.
1190: - name   - The attribute name

1192:   Output Parameter:
1193: . has    - Flag for attribute existence

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

1199:   Level: advanced

1201: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5HasAttribute(), PetscViewerHDF5WriteObjectAttribute(), PetscViewerHDF5ReadObjectAttribute(), PetscViewerHDF5HasObject(), PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
1202: @*/
1203: PetscErrorCode PetscViewerHDF5HasObjectAttribute(PetscViewer viewer, PetscObject obj, const char name[], PetscBool *has)
1204: {

1212:   PetscViewerHDF5CheckNamedObject_Internal(viewer, obj);
1213:   PetscViewerHDF5HasAttribute(viewer, obj->name, name, has);
1214:   return(0);
1215: }

1217: static PetscErrorCode PetscViewerHDF5HasAttribute_Internal(PetscViewer viewer, const char parent[], const char name[], PetscBool *has)
1218: {
1219:   hid_t          h5;
1220:   htri_t         hhas;

1224:   PetscViewerHDF5GetFileId(viewer, &h5);
1225:   PetscStackCallHDF5Return(hhas,H5Aexists_by_name,(h5, parent, name, H5P_DEFAULT));
1226:   *has = hhas ? PETSC_TRUE : PETSC_FALSE;
1227:   return(0);
1228: }

1230: /*
1231:   The variable Petsc_Viewer_HDF5_keyval is used to indicate an MPI attribute that
1232:   is attached to a communicator, in this case the attribute is a PetscViewer.
1233: */
1234: PetscMPIInt Petsc_Viewer_HDF5_keyval = MPI_KEYVAL_INVALID;

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

1239:   Collective on MPI_Comm

1241:   Input Parameter:
1242: . comm - the MPI communicator to share the HDF5 PetscViewer

1244:   Level: intermediate

1246:   Options Database Keys:
1247: . -viewer_hdf5_filename <name>

1249:   Environmental variables:
1250: . PETSC_VIEWER_HDF5_FILENAME

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

1257: .seealso: PetscViewerHDF5Open(), PetscViewerCreate(), PetscViewerDestroy()
1258: @*/
1259: PetscViewer PETSC_VIEWER_HDF5_(MPI_Comm comm)
1260: {
1262:   PetscBool      flg;
1263:   PetscViewer    viewer;
1264:   char           fname[PETSC_MAX_PATH_LEN];
1265:   MPI_Comm       ncomm;

1268:   PetscCommDuplicate(comm,&ncomm,NULL);if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
1269:   if (Petsc_Viewer_HDF5_keyval == MPI_KEYVAL_INVALID) {
1270:     MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN,MPI_COMM_NULL_DELETE_FN,&Petsc_Viewer_HDF5_keyval,0);
1271:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
1272:   }
1273:   MPI_Comm_get_attr(ncomm,Petsc_Viewer_HDF5_keyval,(void**)&viewer,(int*)&flg);
1274:   if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
1275:   if (!flg) { /* PetscViewer not yet created */
1276:     PetscOptionsGetenv(ncomm,"PETSC_VIEWER_HDF5_FILENAME",fname,PETSC_MAX_PATH_LEN,&flg);
1277:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
1278:     if (!flg) {
1279:       PetscStrcpy(fname,"output.h5");
1280:       if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
1281:     }
1282:     PetscViewerHDF5Open(ncomm,fname,FILE_MODE_WRITE,&viewer);
1283:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
1284:     PetscObjectRegisterDestroy((PetscObject)viewer);
1285:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
1286:     MPI_Comm_set_attr(ncomm,Petsc_Viewer_HDF5_keyval,(void*)viewer);
1287:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
1288:   }
1289:   PetscCommDestroy(&ncomm);
1290:   if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_HDF5_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
1291:   PetscFunctionReturn(viewer);
1292: }