Actual source code: hdf5v.c
petsc-3.3-p7 2013-05-11
1: #include <petsc-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: static PetscErrorCode PetscViewerFileClose_HDF5(PetscViewer viewer)
20: {
21: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*)viewer->data;
25: PetscFree(hdf5->filename);
26: if (hdf5->file_id) {
27: H5Fclose(hdf5->file_id);
28: }
29: return(0);
30: }
34: PetscErrorCode PetscViewerDestroy_HDF5(PetscViewer viewer)
35: {
36: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
37: PetscErrorCode ierr;
40: PetscViewerFileClose_HDF5(viewer);
41: if (hdf5->groups) {
42: while(hdf5->groups) {
43: GroupList *tmp = hdf5->groups->next;
45: PetscFree(hdf5->groups->name);
46: PetscFree(hdf5->groups);
47: hdf5->groups = tmp;
48: }
49: }
50: PetscFree(hdf5);
51: return(0);
52: }
54: EXTERN_C_BEGIN
57: PetscErrorCode PetscViewerFileSetMode_HDF5(PetscViewer viewer, PetscFileMode type)
58: {
59: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
63: hdf5->btype = type;
64: return(0);
65: }
66: EXTERN_C_END
68: EXTERN_C_BEGIN
71: PetscErrorCode PetscViewerFileSetName_HDF5(PetscViewer viewer, const char name[])
72: {
73: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
74: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
75: MPI_Info info = MPI_INFO_NULL;
76: #endif
77: hid_t plist_id;
78: herr_t herr;
79: PetscErrorCode ierr;
82: PetscStrallocpy(name, &hdf5->filename);
83: /* Set up file access property list with parallel I/O access */
84: plist_id = H5Pcreate(H5P_FILE_ACCESS);
85: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
86: herr = H5Pset_fapl_mpio(plist_id, ((PetscObject)viewer)->comm, info);CHKERRQ(herr);
87: #endif
88: /* Create or open the file collectively */
89: switch (hdf5->btype) {
90: case FILE_MODE_READ:
91: hdf5->file_id = H5Fopen(name, H5F_ACC_RDONLY, plist_id);
92: break;
93: case FILE_MODE_APPEND:
94: hdf5->file_id = H5Fopen(name, H5F_ACC_RDWR, plist_id);
95: break;
96: case FILE_MODE_WRITE:
97: hdf5->file_id = H5Fcreate(name, H5F_ACC_TRUNC, H5P_DEFAULT, plist_id);
98: break;
99: default:
100: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
101: }
102: if (hdf5->file_id < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB, "H5Fcreate failed for %s", name);
103: viewer->format = PETSC_VIEWER_NOFORMAT;
104: H5Pclose(plist_id);
105: return(0);
106: }
107: EXTERN_C_END
109: EXTERN_C_BEGIN
112: PetscErrorCode PetscViewerCreate_HDF5(PetscViewer v)
113: {
114: PetscViewer_HDF5 *hdf5;
115: PetscErrorCode ierr;
116:
118: PetscNewLog(v, PetscViewer_HDF5, &hdf5);
120: v->data = (void *) hdf5;
121: v->ops->destroy = PetscViewerDestroy_HDF5;
122: v->ops->flush = 0;
123: v->iformat = 0;
124: hdf5->btype = (PetscFileMode) -1;
125: hdf5->filename = 0;
126: hdf5->timestep = -1;
127: hdf5->groups = PETSC_NULL;
129: PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscViewerFileSetName_C","PetscViewerFileSetName_HDF5",
130: PetscViewerFileSetName_HDF5);
131: PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscViewerFileSetMode_C","PetscViewerFileSetMode_HDF5",
132: PetscViewerFileSetMode_HDF5);
133: return(0);
134: }
135: EXTERN_C_END
139: /*@C
140: PetscViewerHDF5Open - Opens a file for HDF5 input/output.
142: Collective on MPI_Comm
144: Input Parameters:
145: + comm - MPI communicator
146: . name - name of file
147: - type - type of file
148: $ FILE_MODE_WRITE - create new file for binary output
149: $ FILE_MODE_READ - open existing file for binary input
150: $ FILE_MODE_APPEND - open existing file for binary output
152: Output Parameter:
153: . hdf5v - PetscViewer for HDF5 input/output to use with the specified file
155: Level: beginner
157: Note:
158: This PetscViewer should be destroyed with PetscViewerDestroy().
160: Concepts: HDF5 files
161: Concepts: PetscViewerHDF5^creating
163: .seealso: PetscViewerASCIIOpen(), PetscViewerSetFormat(), PetscViewerDestroy(),
164: VecView(), MatView(), VecLoad(), MatLoad(),
165: PetscFileMode, PetscViewer
166: @*/
167: PetscErrorCode PetscViewerHDF5Open(MPI_Comm comm, const char name[], PetscFileMode type, PetscViewer *hdf5v)
168: {
170:
172: PetscViewerCreate(comm, hdf5v);
173: PetscViewerSetType(*hdf5v, PETSCVIEWERHDF5);
174: PetscViewerFileSetMode(*hdf5v, type);
175: PetscViewerFileSetName(*hdf5v, name);
176: return(0);
177: }
181: /*@C
182: PetscViewerHDF5GetFileId - Retrieve the file id, this file ID then can be used in direct HDF5 calls
184: Not collective
186: Input Parameter:
187: . viewer - the PetscViewer
189: Output Parameter:
190: . file_id - The file id
192: Level: intermediate
194: .seealso: PetscViewerHDF5Open()
195: @*/
196: PetscErrorCode PetscViewerHDF5GetFileId(PetscViewer viewer, hid_t *file_id)
197: {
198: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
199:
202: if (file_id) *file_id = hdf5->file_id;
203: return(0);
204: }
208: /*@C
209: PetscViewerHDF5PushGroup - Set the current HDF5 group for output
211: Not collective
213: Input Parameters:
214: + viewer - the PetscViewer
215: - name - The group name
217: Level: intermediate
219: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup()
220: @*/
221: PetscErrorCode PetscViewerHDF5PushGroup(PetscViewer viewer, const char *name)
222: {
223: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
224: GroupList *groupNode;
225: PetscErrorCode ierr;
226:
230: PetscMalloc(sizeof(GroupList), &groupNode);
231: PetscStrallocpy(name, (char **) &groupNode->name);
232: groupNode->next = hdf5->groups;
233: hdf5->groups = groupNode;
234: return(0);
235: }
239: /*@C
240: PetscViewerHDF5PopGroup - Return the current HDF5 group for output to the previous value
242: Not collective
244: Input Parameter:
245: . viewer - the PetscViewer
247: Level: intermediate
249: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5GetGroup()
250: @*/
251: PetscErrorCode PetscViewerHDF5PopGroup(PetscViewer viewer)
252: {
253: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
254: GroupList *groupNode;
255: PetscErrorCode ierr;
256:
259: if (!hdf5->groups) SETERRQ(((PetscObject) viewer)->comm, PETSC_ERR_ARG_WRONGSTATE, "HDF5 group stack is empty, cannot pop");
260: groupNode = hdf5->groups;
261: hdf5->groups = hdf5->groups->next;
262: PetscFree(groupNode->name);
263: PetscFree(groupNode);
264: return(0);
265: }
269: /*@C
270: PetscViewerHDF5GetGroup - Get the current HDF5 group for output. If none has been assigned, returns PETSC_NULL.
272: Not collective
274: Input Parameter:
275: . viewer - the PetscViewer
277: Output Parameter:
278: . name - The group name
280: Level: intermediate
282: .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup()
283: @*/
284: PetscErrorCode PetscViewerHDF5GetGroup(PetscViewer viewer, const char **name)
285: {
286: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
287:
291: if (hdf5->groups) {
292: *name = hdf5->groups->name;
293: } else {
294: *name = PETSC_NULL;
295: }
296: return(0);
297: }
301: /*@C
302: PetscViewerHDF5IncrementTimestep - Increments the current timestep for the HDF5 output. Fields are stacked in time.
304: Not collective
306: Input Parameter:
307: . viewer - the PetscViewer
309: Level: intermediate
311: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5SetTimestep(), PetscViewerHDF5GetTimestep()
312: @*/
313: PetscErrorCode PetscViewerHDF5IncrementTimestep(PetscViewer viewer)
314: {
315: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
316:
319: ++hdf5->timestep;
320: return(0);
321: }
325: /*@C
326: PetscViewerHDF5SetTimestep - Set the current timestep for the HDF5 output. Fields are stacked in time. A timestep
327: of -1 disables blocking with timesteps.
329: Not collective
331: Input Parameters:
332: + viewer - the PetscViewer
333: - timestep - The timestep number
335: Level: intermediate
337: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5GetTimestep()
338: @*/
339: PetscErrorCode PetscViewerHDF5SetTimestep(PetscViewer viewer, PetscInt timestep)
340: {
341: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
342:
345: hdf5->timestep = timestep;
346: return(0);
347: }
351: /*@C
352: PetscViewerHDF5GetTimestep - Get the current timestep for the HDF5 output. Fields are stacked in time.
354: Not collective
356: Input Parameter:
357: . viewer - the PetscViewer
359: Output Parameter:
360: . timestep - The timestep number
362: Level: intermediate
364: .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5SetTimestep()
365: @*/
366: PetscErrorCode PetscViewerHDF5GetTimestep(PetscViewer viewer, PetscInt *timestep)
367: {
368: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
369:
373: *timestep = hdf5->timestep;
374: return(0);
375: }
377: #if defined(oldhdf4stuff)
380: PetscErrorCode PetscViewerHDF5WriteSDS(PetscViewer viewer, float *xf, int d, int *dims,int bs)
381: {
382: int i;
383: PetscViewer_HDF5 *vhdf5 = (PetscViewer_HDF5*)viewer->data;
384: int32 sds_id,zero32[3],dims32[3];
388: for (i = 0; i < d; i++) {
389: zero32[i] = 0;
390: dims32[i] = dims[i];
391: }
392: sds_id = SDcreate(vhdf5->sd_id, "Vec", DFNT_FLOAT32, d, dims32);
393: if (sds_id < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"SDcreate failed");
394: SDwritedata(sds_id, zero32, 0, dims32, xf);
395: SDendaccess(sds_id);
396: return(0);
397: }
399: #endif