16 #include <dolfinx/common/MPI.h>
18 #include <dolfinx/common/log.h>
37 static hid_t
open_file(MPI_Comm mpi_comm,
const std::string filename,
38 const std::string mode,
const bool use_mpi_io);
41 static void close_file(
const hid_t hdf5_file_handle);
45 static void flush_file(
const hid_t hdf5_file_handle);
59 const std::string dataset_path,
const T* data,
60 const std::array<std::int64_t, 2> range,
61 const std::vector<std::int64_t> global_size,
62 bool use_mpio,
bool use_chunking);
69 static std::vector<T>
read_dataset(
const hid_t file_handle,
70 const std::string dataset_path,
71 const std::array<std::int64_t, 2> range);
74 static bool has_group(
const hid_t hdf5_file_handle,
75 const std::string group_name);
78 static bool has_dataset(
const hid_t hdf5_file_handle,
79 const std::string dataset_path);
82 static void add_group(
const hid_t hdf5_file_handle,
83 const std::string dataset_path);
87 const std::string dataset_path);
91 const std::string group_name);
94 static std::vector<std::int64_t>
96 const std::string dataset_path);
99 static std::vector<std::string>
dataset_list(
const hid_t hdf5_file_handle,
100 const std::string group_name);
104 const std::string dataset_path,
105 const std::string attribute_name);
108 template <
typename T>
110 const std::string dataset_path,
111 const std::string attribute_name);
114 template <
typename T>
116 add_attribute(
const hid_t hdf5_file_handle,
const std::string dataset_path,
117 const std::string attribute_name,
const T& attribute_value);
121 const std::string dataset_path,
122 const std::string attribute_name);
126 const std::string dataset_path,
127 const std::string attribute_name);
130 static const std::vector<std::string>
131 list_attributes(
const hid_t hdf5_file_handle,
const std::string dataset_path);
149 static herr_t attribute_iteration_function(hid_t loc_id,
const char* name,
150 const H5A_info_t* info,
void* str);
152 template <
typename T>
153 static void add_attribute_value(
const hid_t dset_id,
154 const std::string attribute_name,
155 const T& attribute_value);
157 template <
typename T>
158 static void add_attribute_value(
const hid_t dset_id,
159 const std::string attribute_name,
160 const std::vector<T>& attribute_value);
162 template <
typename T>
163 static void get_attribute_value(
const hid_t attr_type,
const hid_t attr_id,
166 template <
typename T>
167 static void get_attribute_value(
const hid_t attr_type,
const hid_t attr_id,
168 std::vector<T>& attribute_value);
171 template <
typename T>
172 static hid_t hdf5_type()
174 throw std::runtime_error(
"Cannot get HDF5 primitive data type. "
175 "No specialised function for this data type");
182 inline hid_t HDF5Interface::hdf5_type<float>()
184 return H5T_NATIVE_FLOAT;
188 inline hid_t HDF5Interface::hdf5_type<double>()
190 return H5T_NATIVE_DOUBLE;
194 inline hid_t HDF5Interface::hdf5_type<int>()
196 return H5T_NATIVE_INT;
200 inline hid_t HDF5Interface::hdf5_type<std::int64_t>()
202 return H5T_NATIVE_INT64;
206 inline hid_t HDF5Interface::hdf5_type<std::size_t>()
208 if (
sizeof(std::size_t) ==
sizeof(
unsigned long))
209 return H5T_NATIVE_ULONG;
210 else if (
sizeof(std::size_t) ==
sizeof(
unsigned int))
211 return H5T_NATIVE_UINT;
212 throw std::runtime_error(
"Cannot determine size of std::size_t. "
213 "std::size_t is not the same size as long or int");
217 template <
typename T>
219 const hid_t file_handle,
const std::string dataset_path,
const T* data,
220 const std::array<std::int64_t, 2> range,
221 const std::vector<int64_t> global_size,
bool use_mpi_io,
bool use_chunking)
224 const std::size_t rank = global_size.size();
229 throw std::runtime_error(
"Cannot write dataset to HDF5 file"
230 "Only rank 1 and rank 2 dataset are supported");
234 const hid_t h5type = hdf5_type<T>();
237 std::vector<hsize_t> count(global_size.begin(), global_size.end());
238 count[0] = range[1] - range[0];
241 std::vector<hsize_t> offset(rank, 0);
242 offset[0] = range[0];
245 const std::vector<hsize_t> dimsf(global_size.begin(), global_size.end());
251 const hid_t filespace0 = H5Screate_simple(rank, dimsf.data(),
nullptr);
252 assert(filespace0 != HDF5_FAIL);
255 hid_t chunking_properties;
259 hsize_t chunk_size = dimsf[0] / 2;
260 if (chunk_size > 1048576)
261 chunk_size = 1048576;
262 if (chunk_size < 1024)
265 hsize_t chunk_dims[2] = {chunk_size, dimsf[1]};
266 chunking_properties = H5Pcreate(H5P_DATASET_CREATE);
267 H5Pset_chunk(chunking_properties, rank, chunk_dims);
270 chunking_properties = H5P_DEFAULT;
273 const std::string group_name(dataset_path, 0, dataset_path.rfind(
'/'));
278 = H5Dcreate2(file_handle, dataset_path.c_str(), h5type, filespace0,
279 H5P_DEFAULT, chunking_properties, H5P_DEFAULT);
280 assert(dset_id != HDF5_FAIL);
283 status = H5Sclose(filespace0);
284 assert(status != HDF5_FAIL);
287 const hid_t memspace = H5Screate_simple(rank, count.data(),
nullptr);
288 assert(memspace != HDF5_FAIL);
291 const hid_t filespace1 = H5Dget_space(dset_id);
292 status = H5Sselect_hyperslab(filespace1, H5S_SELECT_SET, offset.data(),
293 nullptr, count.data(),
nullptr);
294 assert(status != HDF5_FAIL);
297 const hid_t plist_id = H5Pcreate(H5P_DATASET_XFER);
300 #ifdef H5_HAVE_PARALLEL
301 status = H5Pset_dxpl_mpio(plist_id, H5FD_MPIO_COLLECTIVE);
302 assert(status != HDF5_FAIL);
304 throw std::runtime_error(
"HDF5 library has not been configured with MPI");
309 status = H5Dwrite(dset_id, h5type, memspace, filespace1, plist_id, data);
310 assert(status != HDF5_FAIL);
315 status = H5Pclose(chunking_properties);
316 assert(status != HDF5_FAIL);
320 status = H5Dclose(dset_id);
321 assert(status != HDF5_FAIL);
324 status = H5Sclose(filespace1);
325 assert(status != HDF5_FAIL);
328 status = H5Sclose(memspace);
329 assert(status != HDF5_FAIL);
332 status = H5Pclose(plist_id);
333 assert(status != HDF5_FAIL);
336 template <
typename T>
337 inline std::vector<T>
339 const std::string dataset_path,
340 const std::array<std::int64_t, 2> range)
344 = H5Dopen2(file_handle, dataset_path.c_str(), H5P_DEFAULT);
345 assert(dset_id != HDF5_FAIL);
348 const hid_t dataspace = H5Dget_space(dset_id);
349 assert(dataspace != HDF5_FAIL);
352 const int rank = H5Sget_simple_extent_ndims(dataspace);
356 LOG(WARNING) <<
"HDF5Interface::read_dataset untested for rank > 2.";
359 std::vector<hsize_t> shape(rank);
362 const int ndims = H5Sget_simple_extent_dims(dataspace, shape.data(),
nullptr);
363 assert(ndims == rank);
366 std::vector<hsize_t> offset(rank, 0);
367 std::vector<hsize_t> count = shape;
368 if (range[0] != -1 and range[1] != -1)
370 offset[0] = range[0];
371 count[0] = range[1] - range[0];
378 herr_t status = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET, offset.data(),
379 nullptr, count.data(),
nullptr);
380 assert(status != HDF5_FAIL);
383 const hid_t memspace = H5Screate_simple(rank, count.data(),
nullptr);
384 assert(memspace != HDF5_FAIL);
387 std::size_t data_size = 1;
388 for (std::size_t i = 0; i < count.size(); ++i)
389 data_size *= count[i];
390 std::vector<T> data(data_size);
393 const hid_t h5type = hdf5_type<T>();
395 = H5Dread(dset_id, h5type, memspace, dataspace, H5P_DEFAULT, data.data());
396 assert(status != HDF5_FAIL);
399 status = H5Sclose(dataspace);
400 assert(status != HDF5_FAIL);
403 status = H5Sclose(memspace);
404 assert(status != HDF5_FAIL);
407 status = H5Dclose(dset_id);
408 assert(status != HDF5_FAIL);
413 template <
typename T>
415 const std::string dataset_path,
416 const std::string attribute_name)
422 = H5Oopen(hdf5_file_handle, dataset_path.c_str(), H5P_DEFAULT);
423 assert(dset_id != HDF5_FAIL);
426 const hid_t attr_id = H5Aopen(dset_id, attribute_name.c_str(), H5P_DEFAULT);
427 assert(attr_id != HDF5_FAIL);
428 const hid_t attr_type = H5Aget_type(attr_id);
429 assert(attr_type != HDF5_FAIL);
433 get_attribute_value(attr_type, attr_id, attribute_value);
436 status = H5Tclose(attr_type);
437 assert(status != HDF5_FAIL);
440 status = H5Aclose(attr_id);
441 assert(status != HDF5_FAIL);
444 status = H5Oclose(dset_id);
445 assert(status != HDF5_FAIL);
447 return attribute_value;
450 template <
typename T>
452 const std::string dataset_path,
453 const std::string attribute_name,
454 const T& attribute_value)
458 hid_t dset_id = H5Oopen(hdf5_file_handle, dataset_path.c_str(), H5P_DEFAULT);
459 assert(dset_id != HDF5_FAIL);
462 htri_t has_attr = H5Aexists(dset_id, attribute_name.c_str());
463 assert(has_attr != HDF5_FAIL);
466 herr_t status = H5Adelete(dset_id, attribute_name.c_str());
467 assert(status != HDF5_FAIL);
471 add_attribute_value(dset_id, attribute_name, attribute_value);
474 herr_t status = H5Oclose(dset_id);
475 assert(status != HDF5_FAIL);
484 template <
typename T>
485 inline void HDF5Interface::add_attribute_value(
const hid_t dset_id,
486 const std::string attribute_name,
487 const T& attribute_value)
490 hid_t dataspace_id = H5Screate(H5S_SCALAR);
491 assert(dataspace_id != HDF5_FAIL);
493 const hid_t h5type = hdf5_type<T>();
496 hid_t attribute_id = H5Acreate2(dset_id, attribute_name.c_str(), h5type,
497 dataspace_id, H5P_DEFAULT, H5P_DEFAULT);
498 assert(attribute_id != HDF5_FAIL);
501 herr_t status = H5Awrite(attribute_id, h5type, &attribute_value);
502 assert(status != HDF5_FAIL);
505 status = H5Sclose(dataspace_id);
506 assert(status != HDF5_FAIL);
509 status = H5Aclose(attribute_id);
510 assert(status != HDF5_FAIL);
513 template <
typename T>
515 HDF5Interface::add_attribute_value(
const hid_t dset_id,
516 const std::string attribute_name,
517 const std::vector<T>& attribute_value)
520 const hid_t h5type = hdf5_type<T>();
523 const hsize_t dimsf = attribute_value.size();
524 const hid_t dataspace_id = H5Screate_simple(1, &dimsf,
nullptr);
525 assert(dataspace_id != HDF5_FAIL);
528 const hid_t attribute_id = H5Acreate2(dset_id, attribute_name.c_str(), h5type,
529 dataspace_id, H5P_DEFAULT, H5P_DEFAULT);
530 assert(attribute_id != HDF5_FAIL);
533 herr_t status = H5Awrite(attribute_id, h5type, attribute_value.data());
534 assert(status != HDF5_FAIL);
537 status = H5Sclose(dataspace_id);
538 assert(status != HDF5_FAIL);
541 status = H5Aclose(attribute_id);
542 assert(status != HDF5_FAIL);
547 HDF5Interface::add_attribute_value(
const hid_t dset_id,
548 const std::string attribute_name,
549 const std::string& attribute_value)
552 const hid_t dataspace_id = H5Screate(H5S_SCALAR);
553 assert(dataspace_id != HDF5_FAIL);
556 const hid_t datatype_id = H5Tcopy(H5T_C_S1);
557 herr_t status = H5Tset_size(datatype_id, attribute_value.size());
558 assert(status != HDF5_FAIL);
561 const hid_t attribute_id
562 = H5Acreate2(dset_id, attribute_name.c_str(), datatype_id, dataspace_id,
563 H5P_DEFAULT, H5P_DEFAULT);
564 assert(attribute_id != HDF5_FAIL);
567 status = H5Awrite(attribute_id, datatype_id, attribute_value.c_str());
568 assert(status != HDF5_FAIL);
571 status = H5Sclose(dataspace_id);
572 assert(status != HDF5_FAIL);
575 status = H5Tclose(datatype_id);
576 assert(status != HDF5_FAIL);
579 status = H5Aclose(attribute_id);
580 assert(status != HDF5_FAIL);
583 template <
typename T>
584 inline void HDF5Interface::get_attribute_value(
const hid_t attr_type,
588 const hid_t h5type = hdf5_type<T>();
591 assert(H5Tget_class(attr_type) == H5Tget_class(h5type));
594 herr_t status = H5Aread(attr_id, h5type, &attribute_value);
595 assert(status != HDF5_FAIL);
598 template <
typename T>
599 inline void HDF5Interface::get_attribute_value(
const hid_t attr_type,
601 std::vector<T>& attribute_value)
603 const hid_t h5type = hdf5_type<T>();
606 assert(H5Tget_class(attr_type) == H5Tget_class(h5type));
609 const hid_t dataspace = H5Aget_space(attr_id);
610 assert(dataspace != HDF5_FAIL);
612 hsize_t cur_size[10];
613 hsize_t max_size[10];
614 const int ndims = H5Sget_simple_extent_dims(dataspace, cur_size, max_size);
617 attribute_value.resize(cur_size[0]);
620 herr_t status = H5Aread(attr_id, h5type, attribute_value.data());
621 assert(status != HDF5_FAIL);
624 status = H5Sclose(dataspace);
625 assert(status != HDF5_FAIL);
629 inline void HDF5Interface::get_attribute_value(
const hid_t attr_type,
631 std::string& attribute_value)
634 assert(H5Tget_class(attr_type) == H5T_STRING);
637 const hid_t memtype = H5Tcopy(H5T_C_S1);
638 const int string_length = H5Tget_size(attr_type) + 1;
639 herr_t status = H5Tset_size(memtype, string_length);
640 assert(status != HDF5_FAIL);
645 std::vector<char> attribute_data(string_length);
646 status = H5Aread(attr_id, memtype, attribute_data.data());
647 assert(status != HDF5_FAIL);
649 attribute_value.assign(attribute_data.data());
652 status = H5Tclose(memtype);
653 assert(status != HDF5_FAIL);