DOLFIN
DOLFIN C++ interface
HDF5File.h
1// Copyright (C) 2012 Chris N. Richardson
2//
3// This file is part of DOLFIN.
4//
5// DOLFIN is free software: you can redistribute it and/or modify
6// it under the terms of the GNU Lesser General Public License as published by
7// the Free Software Foundation, either version 3 of the License, or
8// (at your option) any later version.
9//
10// DOLFIN is distributed in the hope that it will be useful,
11// but WITHOUT ANY WARRANTY; without even the implied warranty of
12// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13// GNU Lesser General Public License for more details.
14//
15// You should have received a copy of the GNU Lesser General Public License
16// along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
17//
18// Modified by Garth N. Wells, 2012
19
20#ifndef __DOLFIN_HDF5FILE_H
21#define __DOLFIN_HDF5FILE_H
22
23#ifdef HAS_HDF5
24
25#include <string>
26#include <utility>
27#include <vector>
28
29#include <dolfin/common/MPI.h>
30#include <dolfin/common/Variable.h>
31#include <dolfin/geometry/Point.h>
32#include "HDF5Attribute.h"
33#include "HDF5Interface.h"
34
35namespace dolfin
36{
37
38 class CellType;
39 class Function;
40 class GenericVector;
41 class LocalMeshData;
42 class Mesh;
43 template<typename T> class MeshFunction;
44 template<typename T> class MeshValueCollection;
45 class HDF5Attribute;
46
47 class HDF5File : public Variable
48 {
49
50 public:
51
54 HDF5File(MPI_Comm comm, const std::string filename,
55 const std::string file_mode);
56
58 ~HDF5File();
59
61 void close();
62
64 void flush();
65
67 void write(const std::vector<Point>& points, const std::string name);
68
70 void write(const std::vector<double>& values, const std::string name);
71
73 void write(const GenericVector& x, const std::string name);
74
77 void read(GenericVector& x, const std::string dataset_name,
78 const bool use_partition_from_file) const;
79
81 void write(const Mesh& mesh, const std::string name);
82
85 void write(const Mesh& mesh, const std::size_t cell_dim,
86 const std::string name);
87
89 void write(const Function& u, const std::string name);
90
92 void write(const Function& u, const std::string name, double timestamp);
93
101 void read(Function& u, const std::string name);
102
107 void read(Mesh& mesh, const std::string data_path,
108 bool use_partition_from_file) const;
109
119 void read(Mesh& input_mesh,
120 const std::string topology_path,
121 const std::string geometry_path,
122 const int gdim , const CellType& cell_type,
123 const std::int64_t expected_num_global_cells,
124 const std::int64_t expected_num_global_points,
125 bool use_partition_from_file) const;
126
128 void write(const MeshFunction<std::size_t>& meshfunction,
129 const std::string name);
130
132 void write(const MeshFunction<int>& meshfunction, const std::string name);
133
135 void write(const MeshFunction<double>& meshfunction,
136 const std::string name);
137
139 void write(const MeshFunction<bool>& meshfunction, const std::string name);
140
142 void read(MeshFunction<std::size_t>& meshfunction,
143 const std::string name) const;
144
146 void read(MeshFunction<int>& meshfunction, const std::string name) const;
147
149 void read(MeshFunction<double>& meshfunction,
150 const std::string name) const;
151
153 void read(MeshFunction<bool>& meshfunction, const std::string name) const;
154
156 void write(const MeshValueCollection<std::size_t>& mesh_values,
157 const std::string name);
158
160 void write(const MeshValueCollection<double>& mesh_values,
161 const std::string name);
162
164 void write(const MeshValueCollection<bool>& mesh_values,
165 const std::string name);
166
168 void read(MeshValueCollection<std::size_t>& mesh_values,
169 const std::string name) const;
170
172 void read(MeshValueCollection<double>& mesh_values,
173 const std::string name) const;
174
176 void read(MeshValueCollection<bool>& mesh_values,
177 const std::string name) const;
178
180 bool has_dataset(const std::string dataset_name) const;
181
182 // Get/set attributes of an existing dataset
183 HDF5Attribute attributes(const std::string dataset_name);
184
186 void set_mpi_atomicity(bool atomic);
187
189 bool get_mpi_atomicity() const;
190
191 hid_t h5_id() const
192 { return _hdf5_file_id; }
193
194 private:
195
196 // Friend
197 friend class XDMFFile;
198 friend class TimeSeries;
199
200 // Write a MeshFunction to file
201 template <typename T>
202 void write_mesh_function(const MeshFunction<T>& meshfunction,
203 const std::string name);
204
205 // Read a MeshFunction from file
206 template <typename T>
207 void read_mesh_function(MeshFunction<T>& meshfunction,
208 const std::string name) const;
209
210 // Write a MeshValueCollection to file (old format)
211 template <typename T>
212 void write_mesh_value_collection_old(
213 const MeshValueCollection<T>& mesh_values,
214 const std::string name);
215
216 // Write a MeshValueCollection to file (new version using vertex
217 // indices)
218 template <typename T>
219 void write_mesh_value_collection(const MeshValueCollection<T>& mesh_values,
220 const std::string name);
221
222 // Read a MeshValueCollection from file
223 template <typename T>
224 void read_mesh_value_collection(MeshValueCollection<T>& mesh_values,
225 const std::string name) const;
226
227 // Read a MeshValueCollection (old format)
228 template <typename T>
229 void read_mesh_value_collection_old(MeshValueCollection<T>& mesh_values,
230 const std::string name) const;
231
232 // Write contiguous data to HDF5 data set. Data is flattened into
233 // a 1D array, e.g. [x0, y0, z0, x1, y1, z1] for a vector in 3D
234 template <typename T>
235 void write_data(const std::string dataset_name,
236 const std::vector<T>& data,
237 const std::vector<std::int64_t> global_size,
238 bool use_mpi_io);
239
240 // HDF5 file descriptor/handle
241 hid_t _hdf5_file_id;
242
243 // MPI communicator
244 dolfin::MPI::Comm _mpi_comm;
245 };
246
247 //---------------------------------------------------------------------------
248 // Needs to go here, because of use in XDMFFile.cpp
249 template <typename T>
250 void HDF5File::write_data(const std::string dataset_name,
251 const std::vector<T>& data,
252 const std::vector<std::int64_t> global_size,
253 bool use_mpi_io)
254 {
255 dolfin_assert(_hdf5_file_id > 0);
256 dolfin_assert(global_size.size() > 0);
257
258 // Get number of 'items'
259 std::size_t num_local_items = 1;
260 for (std::size_t i = 1; i < global_size.size(); ++i)
261 num_local_items *= global_size[i];
262 num_local_items = data.size()/num_local_items;
263
264 // Compute offset
265 const std::size_t offset = MPI::global_offset(_mpi_comm.comm(), num_local_items,
266 true);
267 std::pair<std::size_t, std::size_t> range(offset,
268 offset + num_local_items);
269
270 // Write data to HDF5 file
271 const bool chunking = parameters["chunking"];
272 // Ensure dataset starts with '/'
273 std::string dset_name(dataset_name);
274 if (dset_name[0] != '/')
275 dset_name = "/" + dataset_name;
276
277 HDF5Interface::write_dataset(_hdf5_file_id, dset_name, data,
278 range, global_size, use_mpi_io, chunking);
279 }
280 //---------------------------------------------------------------------------
281
282}
283
284#endif
285#endif
Definition: CellType.h:47
Definition: Function.h:66
This class defines a common interface for vectors.
Definition: GenericVector.h:48
Definition: HDF5Attribute.h:39
Definition: HDF5File.h:48
HDF5File(MPI_Comm comm, const std::string filename, const std::string file_mode)
Definition: HDF5File.cpp:58
void read(GenericVector &x, const std::string dataset_name, const bool use_partition_from_file) const
Definition: HDF5File.cpp:185
void flush()
Flush buffered I/O to disk.
Definition: HDF5File.cpp:117
void close()
Close file.
Definition: HDF5File.cpp:109
~HDF5File()
Destructor.
Definition: HDF5File.cpp:104
bool get_mpi_atomicity() const
Get the MPI atomicity.
Definition: HDF5File.cpp:1974
void set_mpi_atomicity(bool atomic)
Set the MPI atomicity.
Definition: HDF5File.cpp:1968
bool has_dataset(const std::string dataset_name) const
Check if dataset exists in HDF5 file.
Definition: HDF5File.cpp:1949
void write(const std::vector< Point > &points, const std::string name)
Write points to file.
Definition: HDF5File.cpp:123
static void write_dataset(const hid_t file_handle, const std::string dataset_path, const std::vector< T > &data, const std::pair< std::int64_t, std::int64_t > range, const std::vector< std::int64_t > global_size, bool use_mpio, bool use_chunking)
Definition: MPI.h:77
MPI_Comm comm() const
Return the underlying MPI_Comm object.
Definition: MPI.cpp:117
static std::size_t global_offset(MPI_Comm comm, std::size_t range, bool exclusive)
Definition: MPI.cpp:185
Definition: MeshValueCollection.h:51
Definition: Mesh.h:84
Definition: TimeSeries.h:47
Common base class for DOLFIN variables.
Definition: Variable.h:36
Parameters parameters
Parameters.
Definition: Variable.h:74
std::string name() const
Return name.
Definition: Variable.cpp:71
Read and write Mesh, Function, MeshFunction and other objects in XDMF.
Definition: XDMFFile.h:78
Definition: adapt.h:30