From b8a1f4ce54ca58f1020cb8102c9935c1d62cd193 Mon Sep 17 00:00:00 2001 From: Daniel Garcia-Sanchez Date: Thu, 2 Aug 2018 20:16:56 +0200 Subject: [PATCH] Add support for all the scalar and complex types Add FullMatrix support in write_dataset Modify the namespace Place get_hdf5_datatype in the internal namespace --- include/deal.II/base/hdf5.h | 79 +++++--- source/base/hdf5.cc | 347 ++++++++++++++++++++++-------------- 2 files changed, 264 insertions(+), 162 deletions(-) diff --git a/include/deal.II/base/hdf5.h b/include/deal.II/base/hdf5.h index 08e9c981fb..4d34cd6646 100644 --- a/include/deal.II/base/hdf5.h +++ b/include/deal.II/base/hdf5.h @@ -32,7 +32,9 @@ DEAL_II_NAMESPACE_OPEN * This set of classes can be used to store the data in the hdf5 file format. * * The classes can be used to improve the interaction with python scripts. - * + */ +// clang-format off +/** * This python script writes the parameters for a deal.ii simulation: * @code * h5_file = h5py.File('simulation.hdf5','w') @@ -45,8 +47,7 @@ DEAL_II_NAMESPACE_OPEN * * C++ deal.ii simulation with MPI HDF5: * @code - * hdf5::File data_file("simulation.hdf5", MPI_COMM_WORLD, - * hdf5::File::Mode::open); + * hdf5::File data_file("simulation.hdf5", MPI_COMM_WORLD, hdf5::File::Mode::open); * hdf5::Group data = data_file.group("data"); * * // Read some parameters @@ -56,15 +57,13 @@ DEAL_II_NAMESPACE_OPEN * const std::string simulation_type = * data.attr("simulation_type"); * - * // Create some distributed datasets - * hdf5::DataSet frequency_dataset(data.template - * create_dataset("frequency", - * std::vector{nb_frequency_points})); + * // Create some distributed datasets. + * // Note that the keyword template is used between data and create because + * // the function create_dataset returns a dependent template. + * auto frequency_dataset = data.template create_dataset("frequency", std::vector{nb_frequency_points}); * - * hdf5::DataSet> displacement(parameters.data.template - * create_dataset>("displacement", - * std::vector{nb_slice_points, - * nb_frequency_points})); + * auto displacement = data.template create_dataset>("displacement", + * std::vector{nb_slice_points, nb_frequency_points})); * * // ... do some calculations ... * // The displacement data is distributed across the MPI processes. @@ -100,12 +99,26 @@ DEAL_II_NAMESPACE_OPEN * displacement = data['displacement'] # complex128 dtype * print('Degrees of freedom: ' + repr(data.attrs['degrees_of_freedom'])) * @endcode + */ +// clang-format on +/** * * @author Daniel Garcia-Sanchez, 2018 */ -namespace hdf5 +namespace HDF5 { + namespace internal + { + // This function gives the HDF5 datatype corresponding to the C++ type. In + // the case of std::complex types the HDF5 handlers are automatically freed + // using the destructor of std::shared_ptr. + template + std::shared_ptr + get_hdf5_datatype(); + } // namespace internal + + /** * General class for the HDF5 objects. * @@ -204,7 +217,7 @@ namespace hdf5 * User's Guide. */ void - write_data(const dealii::FullMatrix &data) const; + write_data(const FullMatrix &data) const; /** * Writes data to a subset of the dataset. T can be double, int, unsigned @@ -262,9 +275,9 @@ namespace hdf5 * User's Guide. */ void - write_data_hyperslab(const dealii::FullMatrix &data, - const std::vector offset, - const std::vector count) const; + write_data_hyperslab(const FullMatrix & data, + const std::vector offset, + const std::vector count) const; void write_data_none() const; const unsigned int rank; @@ -300,6 +313,25 @@ namespace hdf5 Group create_group(const std::string name); + /** + * Creates a dataset. T can be double, int, unsigned int, bool or + * std::complex. + * + * The keyword template has to be used between create_dataset and the name + * of the object because the function create_dataset returns a dependent + * template. See the example in HDF5. + * + * Datatype conversion takes place at the time of a read or write and is + * automatic. See the "Data + * Transfer: Datatype Conversion and Selection" section in the HDF5 + * User's Guide. + */ + template + DataSet + create_dataset(const std::string name, + const std::vector dimensions) const; + /** * Creates and writes data to a dataset. T can be double, int, unsigned int, * bool or std::complex. @@ -312,12 +344,11 @@ namespace hdf5 */ template void - write_dataset(const std::string name, const T &data) const; - template + write_dataset(const std::string name, const std::vector &data) const; /** - * Creates a dataset. T can be double, int, unsigned int, bool or - * std::complex. + * Creates and writes data to a dataset. T can be double, int, unsigned int, + * bool or std::complex. * * Datatype conversion takes place at the time of a read or write and is * automatic. See the section in the HDF5 * User's Guide. */ - DataSet - create_dataset(const std::string name, - const std::vector dimensions) const; + template + void + write_dataset(const std::string name, const FullMatrix &data) const; }; /** @@ -356,7 +387,7 @@ namespace hdf5 */ File(const std::string name, const Mode mode = Mode::create); }; -} // namespace hdf5 +} // namespace HDF5 DEAL_II_NAMESPACE_CLOSE diff --git a/source/base/hdf5.cc b/source/base/hdf5.cc index 73384b4c78..80c28624d7 100644 --- a/source/base/hdf5.cc +++ b/source/base/hdf5.cc @@ -28,9 +28,95 @@ DEAL_II_NAMESPACE_OPEN -namespace hdf5 +namespace HDF5 { + namespace internal + { + template + std::shared_ptr + get_hdf5_datatype() + { + std::shared_ptr t_type; + if (std::is_same::value) + { + t_type = std::shared_ptr(new hid_t); + *t_type = H5T_NATIVE_FLOAT; + } + else if (std::is_same::value) + { + t_type = std::shared_ptr(new hid_t); + *t_type = H5T_NATIVE_DOUBLE; + } + else if (std::is_same::value) + { + t_type = std::shared_ptr(new hid_t); + *t_type = H5T_NATIVE_LDOUBLE; + } + else if (std::is_same::value) + { + t_type = std::shared_ptr(new hid_t); + *t_type = H5T_NATIVE_INT; + } + else if (std::is_same::value) + { + t_type = std::shared_ptr(new hid_t); + *t_type = H5T_NATIVE_UINT; + } + else if (std::is_same>::value) + { + t_type = std::shared_ptr(new hid_t, [](auto pointer) { + // Relase the HDF5 resource + H5Tclose(*pointer); + delete pointer; + }); + *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex)); + // The C++ standards committee agreed to mandate that the storage + // format used for the std::complex type be binary-compatible with + // the C99 type, i.e. an array T[2] with consecutive real [0] and + // imaginary [1] parts. + H5Tinsert(*t_type, "r", 0, H5T_NATIVE_FLOAT); + H5Tinsert(*t_type, "i", sizeof(float), H5T_NATIVE_FLOAT); + } + else if (std::is_same>::value) + { + t_type = std::shared_ptr(new hid_t, [](auto pointer) { + // Relase the HDF5 resource + H5Tclose(*pointer); + delete pointer; + }); + *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex)); + // The C++ standards committee agreed to mandate that the storage + // format used for the std::complex type be binary-compatible with + // the C99 type, i.e. an array T[2] with consecutive real [0] and + // imaginary [1] parts. + H5Tinsert(*t_type, "r", 0, H5T_NATIVE_DOUBLE); + H5Tinsert(*t_type, "i", sizeof(double), H5T_NATIVE_DOUBLE); + } + else if (std::is_same>::value) + { + t_type = std::shared_ptr(new hid_t, [](auto pointer) { + // Relase the HDF5 resource + H5Tclose(*pointer); + delete pointer; + }); + *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex)); + // The C++ standards committee agreed to mandate that the storage + // format used for the std::complex type be binary-compatible with + // the C99 type, i.e. an array T[2] with consecutive real [0] and + // imaginary [1] parts. + H5Tinsert(*t_type, "r", 0, H5T_NATIVE_LDOUBLE); + H5Tinsert(*t_type, "i", sizeof(long double), H5T_NATIVE_LDOUBLE); + } + else + { + Assert(false, ExcInternalError()); + } + return t_type; + } + } // namespace internal + + HDF5Object::HDF5Object(const std::string name, bool mpi) : name(name) , mpi(mpi) @@ -40,43 +126,10 @@ namespace hdf5 T HDF5Object::attr(const std::string attr_name) const { - std::shared_ptr t_type; + std::shared_ptr t_type = internal::get_hdf5_datatype(); T value; hid_t attr; - if (std::is_same::value) - { - t_type = std::shared_ptr(new hid_t); - *t_type = H5T_NATIVE_DOUBLE; - } - else if (std::is_same::value) - { - t_type = std::shared_ptr(new hid_t); - *t_type = H5T_NATIVE_INT; - } - else if (std::is_same::value) - { - t_type = std::shared_ptr(new hid_t); - *t_type = H5T_NATIVE_UINT; - } - else if (std::is_same>::value) - { - t_type = std::shared_ptr(new hid_t, [](auto pointer) { - // Relase the HDF5 resource - H5Tclose(*pointer); - delete pointer; - }); - *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex)); - // The C++ standards committee agreed to mandate that the storage - // format used for the std::complex type be binary-compatible with the - // C99 type, i.e. an array T[2] with consecutive real [0] and imaginary - // [1] parts. - H5Tinsert(*t_type, "r", 0, H5T_NATIVE_DOUBLE); - H5Tinsert(*t_type, "i", sizeof(double), H5T_NATIVE_DOUBLE); - } - else - { - Assert(false, ExcInternalError()); - } + attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT); H5Aread(attr, *t_type, &value); @@ -151,7 +204,7 @@ namespace hdf5 } template <> - std::vector> + FullMatrix HDF5Object::attr(const std::string attr_name) const { hid_t attr; @@ -162,8 +215,7 @@ namespace hdf5 hsize_t dims[2]; H5Sget_simple_extent_dims(attr_space, dims, NULL); - std::vector> matrix_value(dims[0], - std::vector(dims[1])); + FullMatrix matrix_value(dims[0], dims[1]); double *hdf5_data = (double *)malloc(dims[0] * dims[1] * sizeof(double)); H5Aread(attr, H5T_NATIVE_DOUBLE, hdf5_data); @@ -190,42 +242,8 @@ namespace hdf5 { hid_t attr; hid_t aid; - std::shared_ptr t_type; + std::shared_ptr t_type = internal::get_hdf5_datatype(); - if (std::is_same::value) - { - t_type = std::shared_ptr(new hid_t); - *t_type = H5T_NATIVE_DOUBLE; - } - else if (std::is_same::value) - { - t_type = std::shared_ptr(new hid_t); - *t_type = H5T_NATIVE_INT; - } - else if (std::is_same::value) - { - t_type = std::shared_ptr(new hid_t); - *t_type = H5T_NATIVE_UINT; - } - else if (std::is_same>::value) - { - t_type = std::shared_ptr(new hid_t, [](auto pointer) { - // Relase the HDF5 resource - H5Tclose(*pointer); - delete pointer; - }); - *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex)); - // The C++ standards committee agreed to mandate that the storage - // format used for the std::complex type be binary-compatible with the - // C99 type, i.e. an array T[2] with consecutive real [0] and imaginary - // [1] parts. - H5Tinsert(*t_type, "r", 0, H5T_NATIVE_DOUBLE); - H5Tinsert(*t_type, "i", sizeof(double), H5T_NATIVE_DOUBLE); - } - else - { - Assert(false, ExcInternalError()); - } /* * Create scalar attribute. @@ -256,6 +274,7 @@ namespace hdf5 : HDF5Object(name, mpi) , rank(dimensions.size()) , dimensions(dimensions) + , t_type(internal::get_hdf5_datatype()) { hdf5_reference = std::shared_ptr(new hid_t, [](auto pointer) { // Relase the HDF5 resource @@ -268,40 +287,6 @@ namespace hdf5 delete pointer; }); - if (std::is_same::value) - { - t_type = std::shared_ptr(new hid_t); - *t_type = H5T_NATIVE_DOUBLE; - } - else if (std::is_same::value) - { - t_type = std::shared_ptr(new hid_t); - *t_type = H5T_NATIVE_INT; - } - else if (std::is_same::value) - { - t_type = std::shared_ptr(new hid_t); - *t_type = H5T_NATIVE_UINT; - } - else if (std::is_same>::value) - { - t_type = std::shared_ptr(new hid_t, [](auto pointer) { - // Relase the HDF5 resource - H5Tclose(*pointer); - delete pointer; - }); - *t_type = H5Tcreate(H5T_COMPOUND, sizeof(std::complex)); - // The C++ standards committee agreed to mandate that the storage - // format used for the std::complex type be binary-compatible with the - // C99 type, i.e. an array T[2] with consecutive real [0] and imaginary - // [1] parts. - H5Tinsert(*t_type, "r", 0, H5T_NATIVE_DOUBLE); - H5Tinsert(*t_type, "i", sizeof(double), H5T_NATIVE_DOUBLE); - } - else - { - Assert(false, ExcInternalError()); - } total_size = 1; for (auto &&dimension : dimensions) @@ -368,7 +353,7 @@ namespace hdf5 template void - DataSet::write_data(const dealii::FullMatrix &data) const + DataSet::write_data(const FullMatrix &data) const { AssertDimension(total_size, data.m() * data.n()); hid_t plist; @@ -477,9 +462,9 @@ namespace hdf5 template void - DataSet::write_data_hyperslab(const dealii::FullMatrix &data, - const std::vector offset, - const std::vector count) const + DataSet::write_data_hyperslab(const FullMatrix & data, + const std::vector offset, + const std::vector count) const { AssertDimension(std::accumulate(count.begin(), count.end(), @@ -597,34 +582,33 @@ namespace hdf5 return Group(name, *this, mpi, Mode::create); } - template <> + template + DataSet + Group::create_dataset(const std::string name, + const std::vector dimensions) const + { + return DataSet( + name, *hdf5_reference, dimensions, mpi, DataSet::Mode::create); + } + + template void - Group::write_dataset(const std::string name, - const std::vector &data) const + Group::write_dataset(const std::string name, const std::vector &data) const { std::vector dimensions = {data.size()}; - auto dataset = create_dataset(name, dimensions); + auto dataset = create_dataset(name, dimensions); dataset.write_data(data); } - template <> + template void - Group::write_dataset(const std::string name, - const dealii::FullMatrix &data) const + Group::write_dataset(const std::string name, const FullMatrix &data) const { std::vector dimensions = {data.m(), data.n()}; - auto dataset = create_dataset(name, dimensions); + auto dataset = create_dataset(name, dimensions); dataset.write_data(data); } - template - DataSet - Group::create_dataset(const std::string name, - const std::vector dimensions) const - { - return DataSet( - name, *hdf5_reference, dimensions, mpi, DataSet::Mode::create); - } File::File(const std::string name, const bool mpi, @@ -684,36 +668,84 @@ namespace hdf5 // explicit instantiations of classes + template class DataSet; template class DataSet; + template class DataSet; + template class DataSet>; template class DataSet>; + template class DataSet>; + template class DataSet; + template class DataSet; + // explicit instantiations of functions + template float + HDF5Object::attr(const std::string attr_name) const; template double HDF5Object::attr(const std::string attr_name) const; + template long double + HDF5Object::attr(const std::string attr_name) const; + template std::complex + HDF5Object::attr>(const std::string attr_name) const; + template std::complex + HDF5Object::attr>(const std::string attr_name) const; + template std::complex + HDF5Object::attr>( + const std::string attr_name) const; template int HDF5Object::attr(const std::string attr_name) const; template unsigned int HDF5Object::attr(const std::string attr_name) const; - template std::complex - HDF5Object::attr>(const std::string attr_name) const; // The specialization of HDF5Object::attr has been defined above + template void + HDF5Object::write_attr(const std::string attr_name, float value) const; template void HDF5Object::write_attr(const std::string attr_name, double value) const; template void - HDF5Object::write_attr(const std::string attr_name, int value) const; + HDF5Object::write_attr(const std::string attr_name, + long double value) const; template void - HDF5Object::write_attr(const std::string attr_name, - unsigned int value) const; + HDF5Object::write_attr>(const std::string attr_name, + std::complex value) const; template void HDF5Object::write_attr>( const std::string attr_name, std::complex value) const; + template void + HDF5Object::write_attr>( + const std::string attr_name, + std::complex value) const; + template void + HDF5Object::write_attr(const std::string attr_name, int value) const; + template void + HDF5Object::write_attr(const std::string attr_name, + unsigned int value) const; + + template DataSet + Group::create_dataset(const std::string name, + const std::vector dimensions) const; template DataSet Group::create_dataset(const std::string name, const std::vector dimensions) const; + template DataSet + Group::create_dataset( + const std::string name, + const std::vector dimensions) const; + template DataSet> + Group::create_dataset>( + const std::string name, + const std::vector dimensions) const; + template DataSet> + Group::create_dataset>( + const std::string name, + const std::vector dimensions) const; + template DataSet> + Group::create_dataset>( + const std::string name, + const std::vector dimensions) const; template DataSet Group::create_dataset(const std::string name, const std::vector dimensions) const; @@ -721,12 +753,51 @@ namespace hdf5 Group::create_dataset( const std::string name, const std::vector dimensions) const; - template DataSet> - Group::create_dataset>( - const std::string name, - const std::vector dimensions) const; + template void + Group::write_dataset(const std::string name, + const std::vector &data) const; + template void + Group::write_dataset(const std::string name, + const std::vector &data) const; + template void + Group::write_dataset(const std::string name, + const std::vector &data) const; + template void + Group::write_dataset(const std::string name, + const std::vector> &data) const; + template void + Group::write_dataset(const std::string name, + const std::vector> &data) const; + template void + Group::write_dataset( + const std::string name, + const std::vector> &data) const; + template void + Group::write_dataset(const std::string name, + const std::vector &data) const; + template void + Group::write_dataset(const std::string name, + const std::vector &data) const; -} // namespace hdf5 + template void + Group::write_dataset(const std::string name, + const FullMatrix &data) const; + template void + Group::write_dataset(const std::string name, + const FullMatrix &data) const; + template void + Group::write_dataset(const std::string name, + const FullMatrix &data) const; + template void + Group::write_dataset(const std::string name, + const FullMatrix> &data) const; + template void + Group::write_dataset(const std::string name, + const FullMatrix> &data) const; + template void + Group::write_dataset(const std::string name, + const FullMatrix> &data) const; +} // namespace HDF5 DEAL_II_NAMESPACE_CLOSE -- 2.39.5