Actual source code: hdf5v.c
1: #include <private/viewerimpl.h> /*I "petscsys.h" I*/
2: #include <hdf5.h>
4: typedef struct GroupList {
5: const char *name;
6: struct GroupList *next;
7: } GroupList;
9: typedef struct {
10: char *filename;
11: PetscFileMode btype;
12: hid_t file_id;
13: PetscInt timestep;
14: GroupList *groups;
15: } PetscViewer_HDF5;
19: PetscErrorCode PetscViewerDestroy_HDF5(PetscViewer viewer)
20: {
21: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
22: PetscErrorCode ierr;
25: PetscFree(hdf5->filename);
26: if (hdf5->file_id) {
27: H5Fclose(hdf5->file_id);
28: }
29: if (hdf5->groups) {
30: while(hdf5->groups) {
31: GroupList *tmp = hdf5->groups->next;
33: PetscFree(hdf5->groups->name);
34: PetscFree(hdf5->groups);
35: hdf5->groups = tmp;
36: }
37: }
38: PetscFree(hdf5);
39: return(0);
40: }
45: PetscErrorCode PetscViewerFileSetMode_HDF5(PetscViewer viewer, PetscFileMode type)
46: {
47: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
51: hdf5->btype = type;
52: return(0);
53: }
59: PetscErrorCode PetscViewerFileSetName_HDF5(PetscViewer viewer, const char name[])
60: {
61: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
62: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
63: MPI_Info info = MPI_INFO_NULL;
64: #endif
65: hid_t plist_id;
66: herr_t herr;
67: PetscErrorCode ierr;
70: PetscStrallocpy(name, &hdf5->filename);
71: /* Set up file access property list with parallel I/O access */
72: plist_id = H5Pcreate(H5P_FILE_ACCESS);
73: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
74: herr = H5Pset_fapl_mpio(plist_id, ((PetscObject)viewer)->comm, info);CHKERRQ(herr);
75: #endif
76: /* Create or open the file collectively */
77: switch (hdf5->btype) {
78: case FILE_MODE_READ:
79: hdf5->file_id = H5Fopen(name, H5F_ACC_RDONLY, plist_id);
80: break;
81: case FILE_MODE_WRITE:
82: hdf5->file_id = H5Fcreate(name, H5F_ACC_TRUNC, H5P_DEFAULT, plist_id);
83: break;
84: default:
85: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
86: }
87: if (hdf5->file_id < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB, "H5Fcreate failed for %s", name);
88: viewer->format = PETSC_VIEWER_NOFORMAT;
89: H5Pclose(plist_id);
90: return(0);
91: }
97: PetscErrorCode PetscViewerCreate_HDF5(PetscViewer v)
98: {
99: PetscViewer_HDF5 *hdf5;
100: PetscErrorCode ierr;
101:
103: PetscNewLog(v, PetscViewer_HDF5, &hdf5);
105: v->data = (void *) hdf5;
106: v->ops->destroy = PetscViewerDestroy_HDF5;
107: v->ops->flush = 0;
108: v->iformat = 0;
109: hdf5->btype = (PetscFileMode) -1;
110: hdf5->filename = 0;
111: hdf5->timestep = -1;
112: hdf5->groups = PETSC_NULL;
114: PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscViewerFileSetName_C","PetscViewerFileSetName_HDF5",
115: PetscViewerFileSetName_HDF5);
116: PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscViewerFileSetMode_C","PetscViewerFileSetMode_HDF5",
117: PetscViewerFileSetMode_HDF5);
118: return(0);
119: }
124: /*@C
125: PetscViewerHDF5Open - Opens a file for HDF5 input/output.
127: Collective on MPI_Comm
129: Input Parameters:
130: + comm - MPI communicator
131: . name - name of file
132: - type - type of file
133: $ FILE_MODE_WRITE - create new file for binary output
134: $ FILE_MODE_READ - open existing file for binary input
135: $ FILE_MODE_APPEND - open existing file for binary output
137: Output Parameter:
138: . hdf5v - PetscViewer for HDF5 input/output to use with the specified file
140: Level: beginner
142: Note:
143: This PetscViewer should be destroyed with PetscViewerDestroy().
145: Concepts: HDF5 files
146: Concepts: PetscViewerHDF5^creating
148: .seealso: PetscViewerASCIIOpen(), PetscViewerSetFormat(), PetscViewerDestroy(),
149: VecView(), MatView(), VecLoad(), MatLoad(),
150: PetscFileMode, PetscViewer
151: @*/
152: PetscErrorCode PetscViewerHDF5Open(MPI_Comm comm, const char name[], PetscFileMode type, PetscViewer *hdf5v)
153: {
155:
157: PetscViewerCreate(comm, hdf5v);
158: PetscViewerSetType(*hdf5v, PETSCVIEWERHDF5);
159: PetscViewerFileSetMode(*hdf5v, type);
160: PetscViewerFileSetName(*hdf5v, name);
161: return(0);
162: }
166: /*@C
167: PetscViewerHDF5GetFileId - Retrieve the file id, this file ID then can be used in direct HDF5 calls
169: Not collective
171: Input Parameter:
172: . viewer - the PetscViewer
174: Output Parameter:
175: . file_id - The file id
177: Level: intermediate
179: .seealso: PetscViewerHDF5Open()
180: @*/
181: PetscErrorCode PetscViewerHDF5GetFileId(PetscViewer viewer, hid_t *file_id)
182: {
183: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
184:
187: if (file_id) *file_id = hdf5->file_id;
188: return(0);
189: }
193: /*@C
194: PetscViewerHDF5PushGroup - Set the current HDF5 group for output
196: Not collective
198: Input Parameters:
199: + viewer - the PetscViewer
200: - name - The group name
202: Level: intermediate
204: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
205: @*/
206: PetscErrorCode PetscViewerHDF5PushGroup(PetscViewer viewer, const char *name)
207: {
208: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
209: GroupList *groupNode;
210: PetscErrorCode ierr;
211:
215: PetscMalloc(sizeof(GroupList), &groupNode);
216: PetscStrallocpy(name, (char **) &groupNode->name);
217: groupNode->next = hdf5->groups;
218: hdf5->groups = groupNode;
219: return(0);
220: }
224: /*@C
225: PetscViewerHDF5PopGroup - Return the current HDF5 group for output to the previous value
227: Not collective
229: Input Parameter:
230: . viewer - the PetscViewer
232: Level: intermediate
234: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5GetGroup()
235: @*/
236: PetscErrorCode PetscViewerHDF5PopGroup(PetscViewer viewer)
237: {
238: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
239: GroupList *groupNode;
240: PetscErrorCode ierr;
241:
244: if (!hdf5->groups) SETERRQ(((PetscObject) viewer)->comm, PETSC_ERR_ARG_WRONGSTATE, "HDF5 group stack is empty, cannot pop");
245: groupNode = hdf5->groups;
246: hdf5->groups = hdf5->groups->next;
247: PetscFree(groupNode->name);
248: PetscFree(groupNode);
249: return(0);
250: }
254: /*@C
255: PetscViewerHDF5GetGroup - Get the current HDF5 group for output. If none has been assigned, returns PETSC_NULL.
257: Not collective
259: Input Parameter:
260: . viewer - the PetscViewer
262: Output Parameter:
263: . name - The group name
265: Level: intermediate
267: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup()
268: @*/
269: PetscErrorCode PetscViewerHDF5GetGroup(PetscViewer viewer, const char **name)
270: {
271: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
272:
276: if (hdf5->groups) {
277: *name = hdf5->groups->name;
278: } else {
279: *name = PETSC_NULL;
280: }
281: return(0);
282: }
286: /*@C
287: PetscViewerHDF5IncrementTimestep - Increments the current timestep for the HDF5 output. Fields are stacked in time.
289: Not collective
291: Input Parameter:
292: . viewer - the PetscViewer
294: Level: intermediate
296: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5SetTimestep(), PetscViewerHDF5GetTimestep()
297: @*/
298: PetscErrorCode PetscViewerHDF5IncrementTimestep(PetscViewer viewer)
299: {
300: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
301:
304: ++hdf5->timestep;
305: return(0);
306: }
310: /*@C
311: PetscViewerHDF5SetTimestep - Set the current timestep for the HDF5 output. Fields are stacked in time. A timestep
312: of -1 disables blocking with timesteps.
314: Not collective
316: Input Parameters:
317: + viewer - the PetscViewer
318: - timestep - The timestep number
320: Level: intermediate
322: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5GetTimestep()
323: @*/
324: PetscErrorCode PetscViewerHDF5SetTimestep(PetscViewer viewer, PetscInt timestep)
325: {
326: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
327:
330: hdf5->timestep = timestep;
331: return(0);
332: }
336: /*@C
337: PetscViewerHDF5GetTimestep - Get the current timestep for the HDF5 output. Fields are stacked in time.
339: Not collective
341: Input Parameter:
342: . viewer - the PetscViewer
344: Output Parameter:
345: . timestep - The timestep number
347: Level: intermediate
349: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5SetTimestep()
350: @*/
351: PetscErrorCode PetscViewerHDF5GetTimestep(PetscViewer viewer, PetscInt *timestep)
352: {
353: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
354:
358: *timestep = hdf5->timestep;
359: return(0);
360: }
362: #if defined(oldhdf4stuff)
365: PetscErrorCode PetscViewerHDF5WriteSDS(PetscViewer viewer, float *xf, int d, int *dims,int bs)
366: {
367: int i;
368: PetscViewer_HDF5 *vhdf5 = (PetscViewer_HDF5*)viewer->data;
369: int32 sds_id,zero32[3],dims32[3];
373: for (i = 0; i < d; i++) {
374: zero32[i] = 0;
375: dims32[i] = dims[i];
376: }
377: sds_id = SDcreate(vhdf5->sd_id, "Vec", DFNT_FLOAT32, d, dims32);
378: if (sds_id < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"SDcreate failed");
379: SDwritedata(sds_id, zero32, 0, dims32, xf);
380: SDendaccess(sds_id);
381: return(0);
382: }
384: #endif