Actual source code: hdf5v.c
1: #include "../src/sys/viewer/viewerimpl.h" /*I "petsc.h" I*/
2: #include <hdf5.h>
4: typedef struct {
5: char *filename;
6: PetscFileMode btype;
7: hid_t file_id;
8: } PetscViewer_HDF5;
12: PetscErrorCode PetscViewerDestroy_HDF5(PetscViewer viewer)
13: {
14: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
15: PetscErrorCode ierr;
18: PetscFree(hdf5->filename);
19: if (hdf5->file_id) {
20: H5Fclose(hdf5->file_id);
21: }
22: PetscFree(hdf5);
23: return(0);
24: }
29: PetscErrorCode PetscViewerFileSetMode_HDF5(PetscViewer viewer, PetscFileMode type)
30: {
31: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
35: hdf5->btype = type;
36: return(0);
37: }
43: PetscErrorCode PetscViewerFileSetName_HDF5(PetscViewer viewer, const char name[])
44: {
45: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
46: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
47: MPI_Info info = MPI_INFO_NULL;
48: #endif
49: hid_t plist_id;
50: PetscErrorCode ierr;
53: PetscStrallocpy(name, &hdf5->filename);
54: /* Set up file access property list with parallel I/O access */
55: plist_id = H5Pcreate(H5P_FILE_ACCESS);
56: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
57: H5Pset_fapl_mpio(plist_id, ((PetscObject)viewer)->comm, info);
58: #endif
59: /* Create or open the file collectively */
60: switch (hdf5->btype) {
61: case FILE_MODE_READ:
62: hdf5->file_id = H5Fopen(name, H5F_ACC_RDONLY, plist_id);
63: break;
64: case FILE_MODE_WRITE:
65: hdf5->file_id = H5Fcreate(name, H5F_ACC_TRUNC, H5P_DEFAULT, plist_id);
66: break;
67: default:
68: SETERRQ(PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
69: }
70: if (hdf5->file_id < 0) SETERRQ1(PETSC_ERR_LIB, "H5Fcreate failed for %s", name);
71: viewer->format = PETSC_VIEWER_NOFORMAT;
72: H5Pclose(plist_id);
73: return(0);
74: }
80: PetscErrorCode PetscViewerCreate_HDF5(PetscViewer v)
81: {
82: PetscViewer_HDF5 *hdf5;
83: PetscErrorCode ierr;
84:
86: PetscNewLog(v, PetscViewer_HDF5, &hdf5);
87: v->data = (void *) hdf5;
88: v->ops->destroy = PetscViewerDestroy_HDF5;
89: v->ops->flush = 0;
90: v->iformat = 0;
91: hdf5->btype = (PetscFileMode) -1;
92: hdf5->filename = 0;
94: PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscViewerFileSetName_C","PetscViewerFileSetName_HDF5",
95: PetscViewerFileSetName_HDF5);
96: PetscObjectComposeFunctionDynamic((PetscObject)v,"PetscViewerFileSetMode_C","PetscViewerFileSetMode_HDF5",
97: PetscViewerFileSetMode_HDF5);
98: return(0);
99: }
104: /*@C
105: PetscViewerHDF5Open - Opens a file for HDF5 input/output.
107: Collective on MPI_Comm
109: Input Parameters:
110: + comm - MPI communicator
111: . name - name of file
112: - type - type of file
113: $ FILE_MODE_WRITE - create new file for binary output
114: $ FILE_MODE_READ - open existing file for binary input
115: $ FILE_MODE_APPEND - open existing file for binary output
117: Output Parameter:
118: . hd5v - PetscViewer for HDF5 input/output to use with the specified file
120: Level: beginner
122: Note:
123: This PetscViewer should be destroyed with PetscViewerDestroy().
125: Concepts: HDF5 files
126: Concepts: PetscViewerHDF5^creating
128: .seealso: PetscViewerASCIIOpen(), PetscViewerSetFormat(), PetscViewerDestroy(),
129: VecView(), MatView(), VecLoad(), MatLoad(),
130: PetscFileMode, PetscViewer
131: @*/
132: PetscErrorCode PetscViewerHDF5Open(MPI_Comm comm, const char name[], PetscFileMode type, PetscViewer *hdf5v)
133: {
135:
137: PetscViewerCreate(comm, hdf5v);
138: PetscViewerSetType(*hdf5v, PETSC_VIEWER_HDF5);
139: PetscViewerFileSetMode(*hdf5v, type);
140: PetscViewerFileSetName(*hdf5v, name);
141: return(0);
142: }
146: /*@C
147: PetscViewerHDF5GetFileId - Retrieve the file id, this file ID then can be used in direct HDF5 calls
149: Not collective
151: Input Parameter:
152: . viewer - the PetscViewer
154: Output Parameter:
155: . file_id - The file id
157: Level: intermediate
159: .seealso: PetscViewerHDF5Open()
160: @*/
161: PetscErrorCode PetscViewerHDF5GetFileId(PetscViewer viewer, hid_t *file_id)
162: {
163: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data;
164:
167: if (file_id) *file_id = hdf5->file_id;
168: return(0);
169: }