DOLFIN-X
DOLFIN-X C++ interface
HDF5Interface.h
1 // Copyright (C) 2012 Chris N. Richardson and Garth N. Wells
2 //
3 // This file is part of DOLFINX (https://www.fenicsproject.org)
4 //
5 // SPDX-License-Identifier: LGPL-3.0-or-later
6 
7 #pragma once
8 
9 #include <array>
10 #include <cstdint>
11 #include <string>
12 #include <vector>
13 
14 // Note: dolfin/common/MPI.h is included before hdf5.h to avoid the
15 // MPICH_IGNORE_CXX_SEEK issue
16 #include <dolfinx/common/MPI.h>
17 
18 #include <dolfinx/common/log.h>
19 #include <hdf5.h>
20 
21 namespace dolfinx
22 {
23 
24 namespace io
25 {
26 class HDF5File;
27 
31 
33 {
34 #define HDF5_FAIL -1
35 public:
37  static hid_t open_file(MPI_Comm mpi_comm, const std::string filename,
38  const std::string mode, const bool use_mpi_io);
39 
41  static void close_file(const hid_t hdf5_file_handle);
42 
45  static void flush_file(const hid_t hdf5_file_handle);
46 
48  static std::string get_filename(hid_t hdf5_file_handle);
49 
57  template <typename T>
58  static void write_dataset(const hid_t 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);
63 
68  template <typename T>
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);
72 
74  static bool has_group(const hid_t hdf5_file_handle,
75  const std::string group_name);
76 
78  static bool has_dataset(const hid_t hdf5_file_handle,
79  const std::string dataset_path);
80 
82  static void add_group(const hid_t hdf5_file_handle,
83  const std::string dataset_path);
84 
86  static int dataset_rank(const hid_t hdf5_file_handle,
87  const std::string dataset_path);
88 
90  static int num_datasets_in_group(const hid_t hdf5_file_handle,
91  const std::string group_name);
92 
94  static std::vector<std::int64_t>
95  get_dataset_shape(const hid_t hdf5_file_handle,
96  const std::string dataset_path);
97 
99  static std::vector<std::string> dataset_list(const hid_t hdf5_file_handle,
100  const std::string group_name);
101 
103  static const std::string get_attribute_type(const hid_t hdf5_file_handle,
104  const std::string dataset_path,
105  const std::string attribute_name);
106 
108  template <typename T>
109  static T get_attribute(const hid_t hdf5_file_handle,
110  const std::string dataset_path,
111  const std::string attribute_name);
112 
114  template <typename T>
115  static void
116  add_attribute(const hid_t hdf5_file_handle, const std::string dataset_path,
117  const std::string attribute_name, const T& attribute_value);
118 
120  static void delete_attribute(const hid_t hdf5_file_handle,
121  const std::string dataset_path,
122  const std::string attribute_name);
123 
125  static bool has_attribute(const hid_t hdf5_file_handle,
126  const std::string dataset_path,
127  const std::string attribute_name);
128 
130  static const std::vector<std::string>
131  list_attributes(const hid_t hdf5_file_handle, const std::string dataset_path);
132 
139  static void set_mpi_atomicity(const hid_t hdf5_file_handle,
140  const bool atomic);
141 
146  static bool get_mpi_atomicity(const hid_t hdf5_file_handle);
147 
148 private:
149  static herr_t attribute_iteration_function(hid_t loc_id, const char* name,
150  const H5A_info_t* info, void* str);
151 
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);
156 
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);
161 
162  template <typename T>
163  static void get_attribute_value(const hid_t attr_type, const hid_t attr_id,
164  T& attribute_value);
165 
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);
169 
170  // Return HDF5 data type
171  template <typename T>
172  static hid_t hdf5_type()
173  {
174  throw std::runtime_error("Cannot get HDF5 primitive data type. "
175  "No specialised function for this data type");
176  return 0;
177  }
178 };
180 //---------------------------------------------------------------------------
181 template <>
182 inline hid_t HDF5Interface::hdf5_type<float>()
183 {
184  return H5T_NATIVE_FLOAT;
185 }
186 //---------------------------------------------------------------------------
187 template <>
188 inline hid_t HDF5Interface::hdf5_type<double>()
189 {
190  return H5T_NATIVE_DOUBLE;
191 }
192 //---------------------------------------------------------------------------
193 template <>
194 inline hid_t HDF5Interface::hdf5_type<int>()
195 {
196  return H5T_NATIVE_INT;
197 }
198 //---------------------------------------------------------------------------
199 template <>
200 inline hid_t HDF5Interface::hdf5_type<std::int64_t>()
201 {
202  return H5T_NATIVE_INT64;
203 }
204 //---------------------------------------------------------------------------
205 template <>
206 inline hid_t HDF5Interface::hdf5_type<std::size_t>()
207 {
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");
214  return 0;
215 }
216 //---------------------------------------------------------------------------
217 template <typename T>
218 inline void HDF5Interface::write_dataset(
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)
222 {
223  // Data rank
224  const std::size_t rank = global_size.size();
225  assert(rank != 0);
226 
227  if (rank > 2)
228  {
229  throw std::runtime_error("Cannot write dataset to HDF5 file"
230  "Only rank 1 and rank 2 dataset are supported");
231  }
232 
233  // Get HDF5 data type
234  const hid_t h5type = hdf5_type<T>();
235 
236  // Hyperslab selection parameters
237  std::vector<hsize_t> count(global_size.begin(), global_size.end());
238  count[0] = range[1] - range[0];
239 
240  // Data offsets
241  std::vector<hsize_t> offset(rank, 0);
242  offset[0] = range[0];
243 
244  // Dataset dimensions
245  const std::vector<hsize_t> dimsf(global_size.begin(), global_size.end());
246 
247  // Generic status report
248  herr_t status;
249 
250  // Create a global data space
251  const hid_t filespace0 = H5Screate_simple(rank, dimsf.data(), nullptr);
252  assert(filespace0 != HDF5_FAIL);
253 
254  // Set chunking parameters
255  hid_t chunking_properties;
256  if (use_chunking)
257  {
258  // Set chunk size and limit to 1kB min/1MB max
259  hsize_t chunk_size = dimsf[0] / 2;
260  if (chunk_size > 1048576)
261  chunk_size = 1048576;
262  if (chunk_size < 1024)
263  chunk_size = 1024;
264 
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);
268  }
269  else
270  chunking_properties = H5P_DEFAULT;
271 
272  // Check that group exists and recursively create if required
273  const std::string group_name(dataset_path, 0, dataset_path.rfind('/'));
274  add_group(file_handle, group_name);
275 
276  // Create global dataset (using dataset_path)
277  const hid_t dset_id
278  = H5Dcreate2(file_handle, dataset_path.c_str(), h5type, filespace0,
279  H5P_DEFAULT, chunking_properties, H5P_DEFAULT);
280  assert(dset_id != HDF5_FAIL);
281 
282  // Close global data space
283  status = H5Sclose(filespace0);
284  assert(status != HDF5_FAIL);
285 
286  // Create a local data space
287  const hid_t memspace = H5Screate_simple(rank, count.data(), nullptr);
288  assert(memspace != HDF5_FAIL);
289 
290  // Create a file dataspace within the global space - a hyperslab
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);
295 
296  // Set parallel access
297  const hid_t plist_id = H5Pcreate(H5P_DATASET_XFER);
298  if (use_mpi_io)
299  {
300 #ifdef H5_HAVE_PARALLEL
301  status = H5Pset_dxpl_mpio(plist_id, H5FD_MPIO_COLLECTIVE);
302  assert(status != HDF5_FAIL);
303 #else
304  throw std::runtime_error("HDF5 library has not been configured with MPI");
305 #endif
306  }
307 
308  // Write local dataset into selected hyperslab
309  status = H5Dwrite(dset_id, h5type, memspace, filespace1, plist_id, data);
310  assert(status != HDF5_FAIL);
311 
312  if (use_chunking)
313  {
314  // Close chunking properties
315  status = H5Pclose(chunking_properties);
316  assert(status != HDF5_FAIL);
317  }
318 
319  // Close dataset collectively
320  status = H5Dclose(dset_id);
321  assert(status != HDF5_FAIL);
322 
323  // Close hyperslab
324  status = H5Sclose(filespace1);
325  assert(status != HDF5_FAIL);
326 
327  // Close local dataset
328  status = H5Sclose(memspace);
329  assert(status != HDF5_FAIL);
330 
331  // Release file-access template
332  status = H5Pclose(plist_id);
333  assert(status != HDF5_FAIL);
334 }
335 //---------------------------------------------------------------------------
336 template <typename T>
337 inline std::vector<T>
338 HDF5Interface::read_dataset(const hid_t file_handle,
339  const std::string dataset_path,
340  const std::array<std::int64_t, 2> range)
341 {
342  // Open the dataset
343  const hid_t dset_id
344  = H5Dopen2(file_handle, dataset_path.c_str(), H5P_DEFAULT);
345  assert(dset_id != HDF5_FAIL);
346 
347  // Open dataspace
348  const hid_t dataspace = H5Dget_space(dset_id);
349  assert(dataspace != HDF5_FAIL);
350 
351  // Get rank of data set
352  const int rank = H5Sget_simple_extent_ndims(dataspace);
353  assert(rank >= 0);
354 
355  if (rank > 2)
356  LOG(WARNING) << "HDF5Interface::read_dataset untested for rank > 2.";
357 
358  // Allocate data for shape
359  std::vector<hsize_t> shape(rank);
360 
361  // Get size in each dimension
362  const int ndims = H5Sget_simple_extent_dims(dataspace, shape.data(), nullptr);
363  assert(ndims == rank);
364 
365  // Hyperslab selection
366  std::vector<hsize_t> offset(rank, 0);
367  std::vector<hsize_t> count = shape;
368  if (range[0] != -1 and range[1] != -1)
369  {
370  offset[0] = range[0];
371  count[0] = range[1] - range[0];
372  }
373  else
374  offset[0] = 0;
375 
376  // Select a block in the dataset beginning at offset[], with
377  // size=count[]
378  herr_t status = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET, offset.data(),
379  nullptr, count.data(), nullptr);
380  assert(status != HDF5_FAIL);
381 
382  // Create a memory dataspace
383  const hid_t memspace = H5Screate_simple(rank, count.data(), nullptr);
384  assert(memspace != HDF5_FAIL);
385 
386  // Create local data to read into
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);
391 
392  // Read data on each process
393  const hid_t h5type = hdf5_type<T>();
394  status
395  = H5Dread(dset_id, h5type, memspace, dataspace, H5P_DEFAULT, data.data());
396  assert(status != HDF5_FAIL);
397 
398  // Close dataspace
399  status = H5Sclose(dataspace);
400  assert(status != HDF5_FAIL);
401 
402  // Close memspace
403  status = H5Sclose(memspace);
404  assert(status != HDF5_FAIL);
405 
406  // Close dataset
407  status = H5Dclose(dset_id);
408  assert(status != HDF5_FAIL);
409 
410  return data;
411 }
412 //---------------------------------------------------------------------------
413 template <typename T>
414 inline T HDF5Interface::get_attribute(hid_t hdf5_file_handle,
415  const std::string dataset_path,
416  const std::string attribute_name)
417 {
418  herr_t status;
419 
420  // Open dataset or group by name
421  const hid_t dset_id
422  = H5Oopen(hdf5_file_handle, dataset_path.c_str(), H5P_DEFAULT);
423  assert(dset_id != HDF5_FAIL);
424 
425  // Open attribute by name and get its type
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);
430 
431  // Specific code for each type of data template
432  T attribute_value;
433  get_attribute_value(attr_type, attr_id, attribute_value);
434 
435  // Close attribute type
436  status = H5Tclose(attr_type);
437  assert(status != HDF5_FAIL);
438 
439  // Close attribute
440  status = H5Aclose(attr_id);
441  assert(status != HDF5_FAIL);
442 
443  // Close dataset or group
444  status = H5Oclose(dset_id);
445  assert(status != HDF5_FAIL);
446 
447  return attribute_value;
448 }
449 //--------------------------------------------------------------------------
450 template <typename T>
451 inline void HDF5Interface::add_attribute(const hid_t hdf5_file_handle,
452  const std::string dataset_path,
453  const std::string attribute_name,
454  const T& attribute_value)
455 {
456 
457  // Open named dataset or group
458  hid_t dset_id = H5Oopen(hdf5_file_handle, dataset_path.c_str(), H5P_DEFAULT);
459  assert(dset_id != HDF5_FAIL);
460 
461  // Check if attribute already exists and delete if so
462  htri_t has_attr = H5Aexists(dset_id, attribute_name.c_str());
463  assert(has_attr != HDF5_FAIL);
464  if (has_attr > 0)
465  {
466  herr_t status = H5Adelete(dset_id, attribute_name.c_str());
467  assert(status != HDF5_FAIL);
468  }
469 
470  // Add attribute of appropriate type
471  add_attribute_value(dset_id, attribute_name, attribute_value);
472 
473  // Close dataset or group
474  herr_t status = H5Oclose(dset_id);
475  assert(status != HDF5_FAIL);
476 }
477 //---------------------------------------------------------------------------
478 // Specialised member functions (must be inlined to avoid link errors)
479 //---------------------------------------------------------------------------
480 
481 // Template for simple types (e.g. size_t, double, int etc.) and
482 // vectors of these
483 // Specialization below for string
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)
488 {
489  // Create a scalar dataspace
490  hid_t dataspace_id = H5Screate(H5S_SCALAR);
491  assert(dataspace_id != HDF5_FAIL);
492 
493  const hid_t h5type = hdf5_type<T>();
494 
495  // Create attribute of type std::size_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);
499 
500  // Write attribute to dataset
501  herr_t status = H5Awrite(attribute_id, h5type, &attribute_value);
502  assert(status != HDF5_FAIL);
503 
504  // Close dataspace
505  status = H5Sclose(dataspace_id);
506  assert(status != HDF5_FAIL);
507 
508  // Close attribute
509  status = H5Aclose(attribute_id);
510  assert(status != HDF5_FAIL);
511 }
512 //---------------------------------------------------------------------------
513 template <typename T>
514 inline void
515 HDF5Interface::add_attribute_value(const hid_t dset_id,
516  const std::string attribute_name,
517  const std::vector<T>& attribute_value)
518 {
519 
520  const hid_t h5type = hdf5_type<T>();
521 
522  // Create a vector dataspace
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);
526 
527  // Create an attribute of type size_t in the dataspace
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);
531 
532  // Write attribute to dataset
533  herr_t status = H5Awrite(attribute_id, h5type, attribute_value.data());
534  assert(status != HDF5_FAIL);
535 
536  // Close dataspace
537  status = H5Sclose(dataspace_id);
538  assert(status != HDF5_FAIL);
539 
540  // Close attribute
541  status = H5Aclose(attribute_id);
542  assert(status != HDF5_FAIL);
543 }
544 //---------------------------------------------------------------------------
545 template <>
546 inline void
547 HDF5Interface::add_attribute_value(const hid_t dset_id,
548  const std::string attribute_name,
549  const std::string& attribute_value)
550 {
551  // Create a scalar dataspace
552  const hid_t dataspace_id = H5Screate(H5S_SCALAR);
553  assert(dataspace_id != HDF5_FAIL);
554 
555  // Copy basic string type from HDF5 types and set string length
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);
559 
560  // Create attribute in the dataspace with the given string
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);
565 
566  // Write attribute to dataset
567  status = H5Awrite(attribute_id, datatype_id, attribute_value.c_str());
568  assert(status != HDF5_FAIL);
569 
570  // Close dataspace
571  status = H5Sclose(dataspace_id);
572  assert(status != HDF5_FAIL);
573 
574  // Close string type
575  status = H5Tclose(datatype_id);
576  assert(status != HDF5_FAIL);
577 
578  // Close attribute
579  status = H5Aclose(attribute_id);
580  assert(status != HDF5_FAIL);
581 }
582 //--------------------------------------------------------------------------
583 template <typename T>
584 inline void HDF5Interface::get_attribute_value(const hid_t attr_type,
585  const hid_t attr_id,
586  T& attribute_value)
587 {
588  const hid_t h5type = hdf5_type<T>();
589 
590  // FIXME: more complete check of type
591  assert(H5Tget_class(attr_type) == H5Tget_class(h5type));
592 
593  // Read value
594  herr_t status = H5Aread(attr_id, h5type, &attribute_value);
595  assert(status != HDF5_FAIL);
596 }
597 //---------------------------------------------------------------------------
598 template <typename T>
599 inline void HDF5Interface::get_attribute_value(const hid_t attr_type,
600  const hid_t attr_id,
601  std::vector<T>& attribute_value)
602 {
603  const hid_t h5type = hdf5_type<T>();
604 
605  // FIXME: more complete check of type
606  assert(H5Tget_class(attr_type) == H5Tget_class(h5type));
607 
608  // get dimensions of attribute array, check it is one-dimensional
609  const hid_t dataspace = H5Aget_space(attr_id);
610  assert(dataspace != HDF5_FAIL);
611 
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);
615  assert(ndims == 1);
616 
617  attribute_value.resize(cur_size[0]);
618 
619  // Read value to vector
620  herr_t status = H5Aread(attr_id, h5type, attribute_value.data());
621  assert(status != HDF5_FAIL);
622 
623  // Close dataspace
624  status = H5Sclose(dataspace);
625  assert(status != HDF5_FAIL);
626 }
627 //---------------------------------------------------------------------------
628 template <>
629 inline void HDF5Interface::get_attribute_value(const hid_t attr_type,
630  const hid_t attr_id,
631  std::string& attribute_value)
632 {
633  // Check this attribute is a string
634  assert(H5Tget_class(attr_type) == H5T_STRING);
635 
636  // Copy string type from HDF5 types and set length accordingly
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);
641 
642  // FIXME: messy
643  // Copy string value into temporary vector std::vector::data can
644  // be copied into (std::string::data cannot)
645  std::vector<char> attribute_data(string_length);
646  status = H5Aread(attr_id, memtype, attribute_data.data());
647  assert(status != HDF5_FAIL);
648 
649  attribute_value.assign(attribute_data.data());
650 
651  // Close memory type
652  status = H5Tclose(memtype);
653  assert(status != HDF5_FAIL);
654 }
655 //---------------------------------------------------------------------------
657 } // namespace io
658 } // namespace dolfinx
dolfinx::io::HDF5Interface::dataset_list
static std::vector< std::string > dataset_list(const hid_t hdf5_file_handle, const std::string group_name)
Return list all datasets in named group of file.
Definition: HDF5Interface.cpp:418
dolfinx::io::HDF5Interface
This class wraps HDF5 function calls. HDF5 function calls should only appear in a member function of ...
Definition: HDF5Interface.h:32
dolfinx::io::HDF5Interface::add_group
static void add_group(const hid_t hdf5_file_handle, const std::string dataset_path)
Add group to HDF5 file.
Definition: HDF5Interface.cpp:300
dolfinx::io::HDF5Interface::dataset_rank
static int dataset_rank(const hid_t hdf5_file_handle, const std::string dataset_path)
Get dataset rank.
Definition: HDF5Interface.cpp:335
dolfinx::io::HDF5Interface::write_dataset
static void write_dataset(const hid_t file_handle, const std::string dataset_path, const T *data, const std::array< std::int64_t, 2 > range, const std::vector< std::int64_t > global_size, bool use_mpio, bool use_chunking)
Write data to existing HDF file as defined by range blocks on each process data: data to be written,...
dolfinx::io::HDF5Interface::has_dataset
static bool has_dataset(const hid_t hdf5_file_handle, const std::string dataset_path)
Check for existence of dataset in HDF5 file.
Definition: HDF5Interface.cpp:282
dolfinx::io::HDF5Interface::set_mpi_atomicity
static void set_mpi_atomicity(const hid_t hdf5_file_handle, const bool atomic)
Set MPI atomicity. See https://support.hdfgroup.org/HDF5/doc/RM/RM_H5F.html#File-SetMpiAtomicity and ...
Definition: HDF5Interface.cpp:450
dolfinx::io::HDF5Interface::has_group
static bool has_group(const hid_t hdf5_file_handle, const std::string group_name)
Check for existence of group in HDF5 file.
Definition: HDF5Interface.cpp:249
dolfinx::io::HDF5Interface::has_attribute
static bool has_attribute(const hid_t hdf5_file_handle, const std::string dataset_path, const std::string attribute_name)
Check if an attribute exists on a dataset or group.
Definition: HDF5Interface.cpp:227
dolfinx::io::HDF5Interface::get_filename
static std::string get_filename(hid_t hdf5_file_handle)
Get filename.
Definition: HDF5Interface.cpp:88
dolfinx::io::HDF5Interface::get_dataset_shape
static std::vector< std::int64_t > get_dataset_shape(const hid_t hdf5_file_handle, const std::string dataset_path)
Get dataset shape (size of each dimension)
Definition: HDF5Interface.cpp:364
dolfinx::io::HDF5Interface::get_mpi_atomicity
static bool get_mpi_atomicity(const hid_t hdf5_file_handle)
Get MPI atomicity. See https://support.hdfgroup.org/HDF5/doc/RM/RM_H5F.html#File-GetMpiAtomicity and ...
Definition: HDF5Interface.cpp:459
dolfinx::io::HDF5Interface::open_file
static hid_t open_file(MPI_Comm mpi_comm, const std::string filename, const std::string mode, const bool use_mpi_io)
Open HDF5 and return file descriptor.
Definition: HDF5Interface.cpp:18
dolfinx::io::HDF5Interface::close_file
static void close_file(const hid_t hdf5_file_handle)
Close HDF5 file.
Definition: HDF5Interface.cpp:76
dolfinx::io::HDF5Interface::get_attribute
static T get_attribute(const hid_t hdf5_file_handle, const std::string dataset_path, const std::string attribute_name)
Get a named attribute of a dataset of known type.
dolfinx::io::HDF5Interface::delete_attribute
static void delete_attribute(const hid_t hdf5_file_handle, const std::string dataset_path, const std::string attribute_name)
Delete an attribute from a dataset or group.
Definition: HDF5Interface.cpp:173
dolfinx::io::HDF5Interface::list_attributes
static const std::vector< std::string > list_attributes(const hid_t hdf5_file_handle, const std::string dataset_path)
List attributes of dataset or group.
Definition: HDF5Interface.cpp:202
dolfinx::io::HDF5Interface::read_dataset
static std::vector< T > read_dataset(const hid_t file_handle, const std::string dataset_path, const std::array< std::int64_t, 2 > range)
Read data from a HDF5 dataset "dataset_path" as defined by range blocks on each process range: the lo...
dolfinx::io::HDF5Interface::flush_file
static void flush_file(const hid_t hdf5_file_handle)
Flush data to file to improve data integrity after interruption.
Definition: HDF5Interface.cpp:82
dolfinx::io::HDF5Interface::add_attribute
static void add_attribute(const hid_t hdf5_file_handle, const std::string dataset_path, const std::string attribute_name, const T &attribute_value)
Add attribute to dataset or group.
dolfinx::io::HDF5Interface::num_datasets_in_group
static int num_datasets_in_group(const hid_t hdf5_file_handle, const std::string group_name)
Return number of data sets in a group.
Definition: HDF5Interface.cpp:398
dolfinx::io::HDF5Interface::get_attribute_type
static const std::string get_attribute_type(const hid_t hdf5_file_handle, const std::string dataset_path, const std::string attribute_name)
Get type of attribute.
Definition: HDF5Interface.cpp:106