From 73adc1f5879cab53fa8fc5953f1be83aa44e4c02 Mon Sep 17 00:00:00 2001 From: Daniel Garcia-Sanchez Date: Thu, 2 Aug 2018 13:09:54 +0200 Subject: [PATCH] Add HDF5 classes --- include/deal.II/base/hdf5.h | 363 ++++++++++++++++++ source/base/CMakeLists.txt | 1 + source/base/hdf5.cc | 732 ++++++++++++++++++++++++++++++++++++ 3 files changed, 1096 insertions(+) create mode 100644 include/deal.II/base/hdf5.h create mode 100644 source/base/hdf5.cc diff --git a/include/deal.II/base/hdf5.h b/include/deal.II/base/hdf5.h new file mode 100644 index 0000000000..08e9c981fb --- /dev/null +++ b/include/deal.II/base/hdf5.h @@ -0,0 +1,363 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 1998 - 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +#ifndef dealii_hdf5_h +#define dealii_hdf5_h + +#include + +#include + +#include + +#include + +DEAL_II_NAMESPACE_OPEN + +/** + * Namespace containing the HDF5 wrappers. + * + * 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. + * + * This python script writes the parameters for a deal.ii simulation: + * @code + * h5_file = h5py.File('simulation.hdf5','w') + * data = h5_file.create_group('data') + * data.attrs['nb_frequency_points'] = 50 # int + * data.attrs['nb_slice_points'] = 10 # int + * data.attrs['rho'] = 2300.5 # double + * data.attrs['simulation_type'] = 'elastic_equation' # utf8 string + * @endcode + * + * C++ deal.ii simulation with MPI HDF5: + * @code + * hdf5::File data_file("simulation.hdf5", MPI_COMM_WORLD, + * hdf5::File::Mode::open); + * hdf5::Group data = data_file.group("data"); + * + * // Read some parameters + * const int nb_frequency_points = data.attr("nb_frequency_points"); + * const int nb_slice_points = data.attr("nb_cut_points"); + * const double rho = data.attr("rho"); + * 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})); + * + * hdf5::DataSet> displacement(parameters.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. + * // The displacement_data vector is local. The coordinates vector corresponds + * // to the coordinates of the data in the hdf5 dataset. + * if (coordinates.size() > 0) { + * displacement.write_data_selection(displacement_data, coordinates); + * } else { + * // Even if this MPI process has no displacement data to write, + * // it has to contribute to the collective call + * displacement.write_data_none(); + * } + * + * // ... do some calculations ... + * // The whole frequency dataset can be written by one single MPI process, + * // but the rest of the processes have to contribute to the collective call + * if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0) { + * frequency_dataset.write_data(frequency); + * } else { + * frequency_dataset.write_data_none(); + * } + * + * // Write the simulation metadata + * data.write_attr("active_cells", triangulation.n_active_cells()); + * data.write_attr("degrees_of_freedom", dof_handler.n_dofs()); + * @endcode + * + * Read the simulation data with python: + * @code + * h5_file = h5py.File('simulation.hdf5','r+') + * data = h5_file['data'] + * frequency = data['frequency'] # float64 dtype + * displacement = data['displacement'] # complex128 dtype + * print('Degrees of freedom: ' + repr(data.attrs['degrees_of_freedom'])) + * @endcode + * + * @author Daniel Garcia-Sanchez, 2018 + */ +namespace hdf5 + +{ + /** + * General class for the HDF5 objects. + * + * @author Daniel Garcia-Sanchez, 2018 + */ + class HDF5Object + { + protected: + HDF5Object(const std::string name, bool mpi); + + public: + enum class Mode + { + create, + open + }; + + /** + * Reads an attribute. T can be double, int, unsigned int, bool, + * std::complex or std::string. Note that the encoding of + * std::string is UTF8 in order to be compatible with python3. + * + * 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 + T + attr(const std::string attr_name) const; + + /** + * Writes an attribute. T can be double, int, unsigned int, bool, + * std::complex or std::string. Note that the encoding of + * std::string is UTF8 in order to be compatible with python3. + * + * 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 + void + write_attr(const std::string attr_name, T value) const; + const std::string name; + + protected: + std::shared_ptr hdf5_reference; + // If true use mpi hdf5. If false use serial hdf5. + const bool mpi; + }; + + /** + * This class implements a HDF5 DataSet. + * + * @author Daniel Garcia-Sanchez, 2018 + */ + template + class DataSet : public HDF5Object + { + friend class Group; + + protected: + DataSet(const std::string name, + const hid_t & parent_group_id, + std::vector dimensions, + bool mpi, + const Mode mode); + + public: + ~DataSet(); + + /** + * Writes data in the 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 "Data + * Transfer: Datatype Conversion and Selection" section in the HDF5 + * User's Guide. + */ + void + write_data(const std::vector &data) const; + + /** + * Writes data in the 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 "Data + * Transfer: Datatype Conversion and Selection" section in the HDF5 + * User's Guide. + */ + void + write_data(const dealii::FullMatrix &data) const; + + /** + * Writes data to a subset of the dataset. T can be double, int, unsigned + * int, bool or std::complex. + * + * The selected elements can be scattered and take any shape in the dataset. + * For examples, in the 4D case, we will be selecting three points and a 4D + * dataspace has rank 4, so the selection will be described in a 3-by-4 + * array. To select the points (1,1,1,1), (14,6,12,18), and (8,22,30,22), + * the point selection array would be as follows: + * + * 0 0 0 0 + * + * 13 5 11 17 + * + * 7 21 29 21 + * + * 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. + */ + void + write_data_selection(const std::vector & data, + const std::vector coordinates) const; + + /** + * Writes data to a subset of the dataset. T can be double, int, unsigned + * int, bool or std::complex. + * + * The selected elements form a hyperslab in the dataset. + * + * 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. + */ + void + write_data_hyperslab(const std::vector & data, + const std::vector offset, + const std::vector count) const; + + /** + * Writes data to a subset of the dataset. T can be double, int, unsigned + * int, bool or std::complex. + * + * The selected elements form a hyperslab in the dataset. + * + * 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. + */ + void + write_data_hyperslab(const dealii::FullMatrix &data, + const std::vector offset, + const std::vector count) const; + void + write_data_none() const; + const unsigned int rank; + const std::vector dimensions; + + private: + std::shared_ptr dataspace; + std::shared_ptr t_type; + unsigned int total_size; + }; + + /** + * This class implements a HDF5 Group + * + * @author Daniel Garcia-Sanchez, 2018 + */ + class Group : public HDF5Object + { + protected: + Group(const std::string name, + const Group & parent_group, + const bool mpi, + const Mode mode); + Group(const std::string name, const bool mpi); + + public: + Group + group(const std::string name); + + /** + * Creates a group. + */ + Group + create_group(const std::string name); + + /** + * 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 "Data + * Transfer: Datatype Conversion and Selection" section in the HDF5 + * User's Guide. + */ + template + void + write_dataset(const std::string name, const T &data) const; + template + + /** + * Creates 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 "Data + * Transfer: Datatype Conversion and Selection" section in the HDF5 + * User's Guide. + */ + DataSet + create_dataset(const std::string name, + const std::vector dimensions) const; + }; + + /** + * This class implements a HDF5 Group + * + * @author Daniel Garcia-Sanchez, 2018 + */ + class File : public Group + { + private: + File(const std::string name, + const bool mpi, + const MPI_Comm mpi_communicator, + const Mode mode); + + public: + /** + * Creates or opens a hdf5 file. + */ + File(const std::string name, + const MPI_Comm mpi_communicator, + const Mode mode = Mode::create); + + /** + * Creates or opens a hdf5 file in parallel. + */ + File(const std::string name, const Mode mode = Mode::create); + }; +} // namespace hdf5 + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/base/CMakeLists.txt b/source/base/CMakeLists.txt index 1379d6b6e2..34b1e95033 100644 --- a/source/base/CMakeLists.txt +++ b/source/base/CMakeLists.txt @@ -43,6 +43,7 @@ SET(_unity_include_src index_set.cc job_identifier.cc logstream.cc + hdf5.cc mpi.cc multithread_info.cc named_selection.cc diff --git a/source/base/hdf5.cc b/source/base/hdf5.cc new file mode 100644 index 0000000000..73384b4c78 --- /dev/null +++ b/source/base/hdf5.cc @@ -0,0 +1,732 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 - 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +#include + +#include + +#include + +#include + +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + +namespace hdf5 + +{ + HDF5Object::HDF5Object(const std::string name, bool mpi) + : name(name) + , mpi(mpi) + {} + + template + T + HDF5Object::attr(const std::string attr_name) const + { + std::shared_ptr t_type; + 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); + H5Aclose(attr); + return value; + } + + template <> + bool + HDF5Object::attr(const std::string attr_name) const + { + // The enum field generated by h5py can be casted to int + int int_value; + hid_t attr; + attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT); + H5Aread(attr, H5T_NATIVE_INT, &int_value); + H5Aclose(attr); + // The int can be casted to a bool + bool bool_value = (bool)int_value; + return bool_value; + } + + template <> + std::string + HDF5Object::attr(const std::string attr_name) const + { + // Reads a UTF8 python string + // + // code from + // https://support.hdfgroup.org/ftp/HDF5/examples/misc-examples/vlstratt.c + // + // In the case of a variable length string the user does not have to reserve + // memory for string_out. The call HAread will reserve the memory and the + // user has to free the memory. + // + // Todo: + // - In debug mode, check the type of the attribute. If it is not a variable + // length string raise an exception + // - Use valgrind to check that there are no memory leaks + // - Use H5Dvlen_reclaim instead of free + // - Use collective + + char * string_out; + hid_t attr; + hid_t type; + herr_t ret; + + /* Create a datatype to refer to. */ + type = H5Tcopy(H5T_C_S1); + Assert(type >= 0, ExcInternalError()); + + // Python strings are encoded in UTF8 + ret = H5Tset_cset(type, H5T_CSET_UTF8); + Assert(type >= 0, ExcInternalError()); + + ret = H5Tset_size(type, H5T_VARIABLE); + Assert(ret >= 0, ExcInternalError()); + + attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT); + Assert(attr >= 0, ExcInternalError()); + + ret = H5Aread(attr, type, &string_out); + Assert(ret >= 0, ExcInternalError()); + + std::string string_value(string_out); + // The memory of the variable length string has to be freed. + // H5Dvlen_reclaim could be also used + free(string_out); + H5Tclose(type); + H5Aclose(attr); + return string_value; + } + + template <> + std::vector> + HDF5Object::attr(const std::string attr_name) const + { + hid_t attr; + hid_t attr_space; + attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT); + attr_space = H5Aget_space(attr); + AssertDimension(H5Sget_simple_extent_ndims(attr_space), 2); + hsize_t dims[2]; + H5Sget_simple_extent_dims(attr_space, dims, NULL); + + std::vector> matrix_value(dims[0], + std::vector(dims[1])); + + double *hdf5_data = (double *)malloc(dims[0] * dims[1] * sizeof(double)); + H5Aread(attr, H5T_NATIVE_DOUBLE, hdf5_data); + for (unsigned int dim0_idx = 0; dim0_idx < dims[0]; ++dim0_idx) + { + for (unsigned int dim1_idx = 0; dim1_idx < dims[1]; ++dim1_idx) + { + matrix_value[dim0_idx][dim1_idx] = + hdf5_data[dim0_idx * dims[1] + dim1_idx]; + } + } + + H5Sclose(attr_space); + H5Aclose(attr); + + free(hdf5_data); + + return matrix_value; + } + + template + void + HDF5Object::write_attr(const std::string attr_name, T value) const + { + hid_t attr; + hid_t aid; + std::shared_ptr t_type; + + 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. + */ + aid = H5Screate(H5S_SCALAR); + attr = H5Acreate2(*hdf5_reference, + attr_name.data(), + *t_type, + aid, + H5P_DEFAULT, + H5P_DEFAULT); + + /* + * Write scalar attribute. + */ + H5Awrite(attr, *t_type, &value); + + H5Sclose(aid); + H5Aclose(attr); + } + + template + DataSet::DataSet(const std::string name, + const hid_t & parent_group_id, + std::vector dimensions, + const bool mpi, + const Mode mode) + : HDF5Object(name, mpi) + , rank(dimensions.size()) + , dimensions(dimensions) + { + hdf5_reference = std::shared_ptr(new hid_t, [](auto pointer) { + // Relase the HDF5 resource + H5Dclose(*pointer); + delete pointer; + }); + dataspace = std::shared_ptr(new hid_t, [](auto pointer) { + // Relase the HDF5 resource + H5Sclose(*pointer); + 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) + { + total_size *= dimension; + } + + switch (mode) + { + case (Mode::create): + *dataspace = H5Screate_simple(rank, dimensions.data(), NULL); + + *hdf5_reference = H5Dcreate2(parent_group_id, + name.data(), + *t_type, + *dataspace, + H5P_DEFAULT, + H5P_DEFAULT, + H5P_DEFAULT); + break; + case (Mode::open): + // Write the code for open + *hdf5_reference = H5Gopen2(parent_group_id, name.data(), H5P_DEFAULT); + break; + default: + Assert(false, ExcInternalError()); + break; + } + } + + template + DataSet::~DataSet() + { + // Close first the dataset + hdf5_reference.reset(); + // Then close the dataspace + dataspace.reset(); + } + + template + void + DataSet::write_data(const std::vector &data) const + { + AssertDimension(total_size, data.size()); + hid_t plist; + + if (mpi) + { + plist = H5Pcreate(H5P_DATASET_XFER); + H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE); + } + else + { + plist = H5P_DEFAULT; + } + + H5Dwrite(*hdf5_reference, *t_type, H5S_ALL, H5S_ALL, plist, data.data()); + + if (mpi) + { + H5Pclose(plist); + } + } + + template + void + DataSet::write_data(const dealii::FullMatrix &data) const + { + AssertDimension(total_size, data.m() * data.n()); + hid_t plist; + + if (mpi) + { + plist = H5Pcreate(H5P_DATASET_XFER); + H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE); + } + else + { + plist = H5P_DEFAULT; + } + + // The iterator of FullMatrix has to be converted to a pointer to the raw + // data + H5Dwrite(*hdf5_reference, *t_type, H5S_ALL, H5S_ALL, plist, &*data.begin()); + + if (mpi) + { + H5Pclose(plist); + } + } + + template + void + DataSet::write_data_selection(const std::vector & data, + const std::vector coordinates) const + { + AssertDimension(coordinates.size(), data.size() * rank); + hid_t memory_dataspace; + hid_t plist; + std::vector data_dimensions = {data.size()}; + + memory_dataspace = H5Screate_simple(1, data_dimensions.data(), NULL); + H5Sselect_elements(*dataspace, + H5S_SELECT_SET, + data.size(), + coordinates.data()); + + if (mpi) + { + plist = H5Pcreate(H5P_DATASET_XFER); + H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE); + } + else + { + plist = H5P_DEFAULT; + } + + H5Dwrite(*hdf5_reference, + *t_type, + memory_dataspace, + *dataspace, + plist, + data.data()); + + if (mpi) + { + H5Pclose(plist); + } + H5Sclose(memory_dataspace); + } + + template + void + DataSet::write_data_hyperslab(const std::vector & data, + const std::vector offset, + const std::vector count) const + { + AssertDimension(std::accumulate(count.begin(), + count.end(), + 1, + std::multiplies()), + data.size()); + hid_t memory_dataspace; + hid_t plist; + std::vector data_dimensions = {data.size()}; + + memory_dataspace = H5Screate_simple(1, data_dimensions.data(), NULL); + H5Sselect_hyperslab( + *dataspace, H5S_SELECT_SET, offset.data(), NULL, count.data(), NULL); + + if (mpi) + { + plist = H5Pcreate(H5P_DATASET_XFER); + H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE); + } + else + { + plist = H5P_DEFAULT; + } + H5Dwrite(*hdf5_reference, + *t_type, + memory_dataspace, + *dataspace, + plist, + &*data.begin()); + + if (mpi) + { + H5Pclose(plist); + } + H5Sclose(memory_dataspace); + } + + template + void + DataSet::write_data_hyperslab(const dealii::FullMatrix &data, + const std::vector offset, + const std::vector count) const + { + AssertDimension(std::accumulate(count.begin(), + count.end(), + 1, + std::multiplies()), + data.m() * data.n()); + hid_t memory_dataspace; + hid_t plist; + std::vector data_dimensions = {data.m(), data.n()}; + + memory_dataspace = H5Screate_simple(2, data_dimensions.data(), NULL); + H5Sselect_hyperslab( + *dataspace, H5S_SELECT_SET, offset.data(), NULL, count.data(), NULL); + + if (mpi) + { + plist = H5Pcreate(H5P_DATASET_XFER); + H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE); + } + else + { + plist = H5P_DEFAULT; + } + H5Dwrite(*hdf5_reference, + *t_type, + memory_dataspace, + *dataspace, + plist, + &*data.begin()); + + if (mpi) + { + H5Pclose(plist); + } + H5Sclose(memory_dataspace); + } + + template + void + DataSet::write_data_none() const + { + hid_t memory_dataspace; + hid_t plist; + std::vector data_dimensions = {0}; + + memory_dataspace = H5Screate_simple(1, data_dimensions.data(), NULL); + H5Sselect_none(*dataspace); + + if (mpi) + { + plist = H5Pcreate(H5P_DATASET_XFER); + H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE); + } + else + { + plist = H5P_DEFAULT; + } + + // The pointer of data can safely be NULL, see the discussion at the HDF5 + // forum: + // https://forum.hdfgroup.org/t/parallel-i-o-does-not-support-filters-yet/884/17 + H5Dwrite( + *hdf5_reference, *t_type, memory_dataspace, *dataspace, plist, NULL); + + if (mpi) + { + H5Pclose(plist); + } + H5Sclose(memory_dataspace); + } + + Group::Group(const std::string name, + const Group & parentGroup, + const bool mpi, + const Mode mode) + : HDF5Object(name, mpi) + { + hdf5_reference = std::shared_ptr(new hid_t, [](auto pointer) { + // Relase the HDF5 resource + H5Gclose(*pointer); + delete pointer; + }); + switch (mode) + { + case (Mode::create): + *hdf5_reference = H5Gcreate2(*(parentGroup.hdf5_reference), + name.data(), + H5P_DEFAULT, + H5P_DEFAULT, + H5P_DEFAULT); + break; + case (Mode::open): + *hdf5_reference = + H5Gopen2(*(parentGroup.hdf5_reference), name.data(), H5P_DEFAULT); + break; + default: + Assert(false, ExcInternalError()); + break; + } + } + + Group::Group(const std::string name, const bool mpi) + : HDF5Object(name, mpi) + {} + + Group + Group::group(const std::string name) + { + return Group(name, *this, mpi, Mode::open); + } + + Group + Group::create_group(const std::string name) + { + return Group(name, *this, mpi, Mode::create); + } + + template <> + void + Group::write_dataset(const std::string name, + const std::vector &data) const + { + std::vector dimensions = {data.size()}; + auto dataset = create_dataset(name, dimensions); + dataset.write_data(data); + } + + template <> + void + Group::write_dataset(const std::string name, + const dealii::FullMatrix &data) const + { + std::vector dimensions = {data.m(), data.n()}; + 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, + const MPI_Comm mpi_communicator, + const Mode mode) + : Group(name, mpi) + { + hdf5_reference = std::shared_ptr(new hid_t, [](auto pointer) { + // Relase the HDF5 resource + H5Fclose(*pointer); + delete pointer; + }); + + hid_t plist; + const MPI_Info info = MPI_INFO_NULL; + + if (mpi) + { + plist = H5Pcreate(H5P_FILE_ACCESS); + H5Pset_fapl_mpio(plist, mpi_communicator, info); + } + else + { + plist = H5P_DEFAULT; + } + + switch (mode) + { + case (Mode::create): + *hdf5_reference = + H5Fcreate(name.data(), H5F_ACC_TRUNC, H5P_DEFAULT, plist); + break; + case (Mode::open): + *hdf5_reference = H5Fopen(name.data(), H5F_ACC_RDWR, plist); + break; + default: + Assert(false, ExcInternalError()); + break; + } + + if (mpi) + { + // Relase the HDF5 resource + H5Pclose(plist); + } + } + + File::File(const std::string name, + const MPI_Comm mpi_communicator, + const Mode mode) + : File(name, true, mpi_communicator, mode) + {} + + File::File(const std::string name, const Mode mode) + : File(name, false, MPI_COMM_NULL, mode) + {} + + + // explicit instantiations of classes + template class DataSet; + template class DataSet>; + + // explicit instantiations of functions + template double + 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, + double 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 void + HDF5Object::write_attr>( + const std::string attr_name, + std::complex 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; + + +} // namespace hdf5 + +DEAL_II_NAMESPACE_CLOSE -- 2.39.5