DOLFIN-X
DOLFIN-X C++ interface
HDF5File.h
1 // Copyright (C) 2012 Chris N. Richardson
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 "HDF5Interface.h"
10 #include <Eigen/Dense>
11 #include <dolfinx/common/MPI.h>
12 #include <memory>
13 #include <string>
14 #include <utility>
15 #include <vector>
16 
17 namespace dolfinx::io
18 {
19 
21 
22 class HDF5File
23 {
24 
25 public:
28  HDF5File(MPI_Comm comm, const std::string& filename,
29  const std::string& file_mode);
30 
32  ~HDF5File();
33 
35  void close();
36 
38  void flush();
39 
41  bool has_dataset(const std::string& dataset_name) const;
42 
44  void set_mpi_atomicity(bool atomic);
45 
47  bool get_mpi_atomicity() const;
48 
50  hid_t h5_id() const { return _hdf5_file_id; }
51 
53  bool chunking = false;
54 
57  template <typename T>
58  void write_data(const std::string& dataset_name, const std::vector<T>& data,
59  const std::vector<std::int64_t>& global_size,
60  bool use_mpi_io);
61 
64  template <typename T>
65  void write_data(const std::string& dataset_name,
66  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic,
67  Eigen::RowMajor>& data,
68  bool use_mpi_io);
69 
70 private:
71  // HDF5 file descriptor/handle
72  hid_t _hdf5_file_id;
73 
74  // MPI communicator
75  dolfinx::MPI::Comm _mpi_comm;
76 };
77 
78 //---------------------------------------------------------------------------
79 // Needs to go here, because of use in XDMFFile.cpp
80 template <typename T>
81 void HDF5File::write_data(const std::string& dataset_name,
82  const std::vector<T>& data,
83  const std::vector<std::int64_t>& global_size,
84  bool use_mpi_io)
85 {
86  assert(_hdf5_file_id > 0);
87  assert(global_size.size() > 0);
88 
89  // Get number of 'items'
90  std::int64_t num_local_items = 1;
91  for (std::size_t i = 1; i < global_size.size(); ++i)
92  num_local_items *= global_size[i];
93  num_local_items = data.size() / num_local_items;
94 
95  // Compute offset
96  const std::int64_t offset
97  = MPI::global_offset(_mpi_comm.comm(), num_local_items, true);
98  std::array<std::int64_t, 2> range = {{offset, offset + num_local_items}};
99 
100  // Write data to HDF5 file. Ensure dataset starts with '/'.
101  std::string dset_name(dataset_name);
102  if (dset_name[0] != '/')
103  dset_name = "/" + dataset_name;
104 
105  HDF5Interface::write_dataset(_hdf5_file_id, dset_name, data.data(), range,
106  global_size, use_mpi_io, chunking);
107 }
108 //-----------------------------------------------------------------------------
109 template <typename T>
110 void HDF5File::write_data(const std::string& dataset_name,
111  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic,
112  Eigen::RowMajor>& data,
113  bool use_mpi_io)
114 {
115  assert(_hdf5_file_id > 0);
116 
117  // Compute offset
118  const std::int64_t offset
119  = MPI::global_offset(_mpi_comm.comm(), data.rows(), true);
120  std::array<std::int64_t, 2> range = {{offset, offset + data.rows()}};
121 
122  // Write data to HDF5 file. Ensure dataset starts with '/'.
123  std::string dset_name(dataset_name);
124  if (dset_name[0] != '/')
125  dset_name = "/" + dataset_name;
126 
127  std::int64_t global_rows = 0;
128  const std::int64_t local_rows = data.rows();
129  MPI_Allreduce(&local_rows, &global_rows, 1, MPI_INT64_T, MPI_SUM,
130  _mpi_comm.comm());
131 
132  std::vector<std::int64_t> global_size = {global_rows, data.cols()};
133  if (data.cols() == 1)
134  global_size = {global_rows};
135 
136  HDF5Interface::write_dataset(_hdf5_file_id, dset_name, data.data(), range,
137  global_size, use_mpi_io, chunking);
138 }
139 //---------------------------------------------------------------------------
140 } // namespace dolfinx::io
dolfinx::MPI::global_offset
static std::size_t global_offset(MPI_Comm comm, std::size_t range, bool exclusive)
Find global offset (index) (wrapper for MPI_(Ex)Scan with MPI_SUM as reduction op)
Definition: MPI.cpp:90
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::HDF5File::has_dataset
bool has_dataset(const std::string &dataset_name) const
Check if dataset exists in HDF5 file.
Definition: HDF5File.cpp:76
dolfinx::io
Support for file IO.
Definition: cells.h:14
dolfinx::io::HDF5File::write_data
void write_data(const std::string &dataset_name, const std::vector< T > &data, const std::vector< std::int64_t > &global_size, bool use_mpi_io)
Write contiguous data to HDF5 data set. Data is flattened into a 1D array, e.g. [x0,...
Definition: HDF5File.h:81
dolfinx::io::HDF5File::set_mpi_atomicity
void set_mpi_atomicity(bool atomic)
Set the MPI atomicity.
Definition: HDF5File.cpp:82
dolfinx::io::HDF5File
Interface to HDF5 files.
Definition: HDF5File.h:22
dolfinx::io::HDF5File::~HDF5File
~HDF5File()
Destructor.
Definition: HDF5File.cpp:60
dolfinx::MPI::Comm
A duplicate MPI communicator and manage lifetime of the communicator.
Definition: MPI.h:34
dolfinx::io::HDF5File::HDF5File
HDF5File(MPI_Comm comm, const std::string &filename, const std::string &file_mode)
Constructor. file_mode should be "a" (append), "w" (write) or "r" (read).
Definition: HDF5File.cpp:20
dolfinx::io::HDF5File::chunking
bool chunking
Chunking parameter - partition data into fixed size blocks for efficiency.
Definition: HDF5File.h:53
dolfinx::io::HDF5File::get_mpi_atomicity
bool get_mpi_atomicity() const
Get the MPI atomicity.
Definition: HDF5File.cpp:88
dolfinx::io::HDF5File::flush
void flush()
Flush buffered I/O to disk.
Definition: HDF5File.cpp:70
dolfinx::io::HDF5File::h5_id
hid_t h5_id() const
Get the file ID.
Definition: HDF5File.h:50
dolfinx::io::HDF5File::close
void close()
Close file.
Definition: HDF5File.cpp:62