From: Daniel Garcia-Sanchez Date: Wed, 12 Aug 2020 08:42:30 +0000 (+0200) Subject: Move template definitions from hdf5.cc to hdf5.h X-Git-Tag: v9.3.0-rc1~1169^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ec591cbf2c335b5e43ad5f7255b6d21ee6702abf;p=dealii.git Move template definitions from hdf5.cc to hdf5.h --- diff --git a/include/deal.II/base/hdf5.h b/include/deal.II/base/hdf5.h index 5e8be8e6c5..52becf1764 100644 --- a/include/deal.II/base/hdf5.h +++ b/include/deal.II/base/hdf5.h @@ -20,11 +20,13 @@ #ifdef DEAL_II_WITH_HDF5 +# include + # include # include -# include +# include DEAL_II_NAMESPACE_OPEN @@ -1103,6 +1105,1049 @@ namespace HDF5 const bool mpi, const MPI_Comm mpi_communicator); }; + + + + // definitions + + 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. + // std::shared_ptr is used instead of std::unique_ptr because the destructor + // of std::shared_ptr doesn't have to be defined in the template argument. + // In the other hand, the destructor of std::unique has to be defined in the + // template argument. Native types such as H5T_NATIVE_DOUBLE do not + // require a destructor, but compound types such as std::complex + // require a destructor to free the HDF5 resources. + template + std::shared_ptr + get_hdf5_datatype() + { + std::shared_ptr t_type; + if (std::is_same::value) + { + return std::make_shared(H5T_NATIVE_FLOAT); + } + else if (std::is_same::value) + { + return std::make_shared(H5T_NATIVE_DOUBLE); + } + else if (std::is_same::value) + { + return std::make_shared(H5T_NATIVE_INT); + } + else if (std::is_same::value) + { + return std::make_shared(H5T_NATIVE_UINT); + } + else if (std::is_same>::value) + { + t_type = std::shared_ptr(new hid_t, [](hid_t *pointer) { + // Release the HDF5 resource + const herr_t ret = H5Tclose(*pointer); + AssertNothrow(ret >= 0, ExcInternalError()); + (void)ret; + 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. + herr_t ret = H5Tinsert(*t_type, "r", 0, H5T_NATIVE_FLOAT); + Assert(ret >= 0, ExcInternalError()); + ret = H5Tinsert(*t_type, "i", sizeof(float), H5T_NATIVE_FLOAT); + Assert(ret >= 0, ExcInternalError()); + (void)ret; + return t_type; + } + else if (std::is_same>::value) + { + t_type = std::shared_ptr(new hid_t, [](hid_t *pointer) { + // Release the HDF5 resource + const herr_t ret = H5Tclose(*pointer); + AssertNothrow(ret >= 0, ExcInternalError()); + (void)ret; + 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. + herr_t ret = H5Tinsert(*t_type, "r", 0, H5T_NATIVE_DOUBLE); + Assert(ret >= 0, ExcInternalError()); + ret = H5Tinsert(*t_type, "i", sizeof(double), H5T_NATIVE_DOUBLE); + Assert(ret >= 0, ExcInternalError()); + (void)ret; + return t_type; + } + + // The function should not reach this point + Assert(false, ExcInternalError()); + return t_type; + } + + + + // Several HDF5 functions such as H5Screate_simple() require a + // one-dimensional array that specifies the size of each dimension of the + // container, see: + // https://support.hdfgroup.org/HDF5/doc1.8/RM/RM_H5S.html#Dataspace-CreateSimple + // The function get_container_dimensions returns a vector with the + // dimensions of the container. + // For a std::vector the function returns std::vector{vector_size} + // For a Vector the function returns std::vector{vector_size} + // For a FullMatrix the function returns std::vector{rows, columns} + template + std::vector + get_container_dimensions(const std::vector &data) + { + std::vector dimensions = {data.size()}; + return dimensions; + } + + + + template + std::vector + get_container_dimensions(const Vector &data) + { + std::vector dimensions = {data.size()}; + return dimensions; + } + + + + template + std::vector + get_container_dimensions(const FullMatrix &data) + { + std::vector dimensions = {data.m(), data.n()}; + return dimensions; + } + + + + // This function returns the total size of the container. + // For a std::vector the function returns int(vector_size) + // For a Vector the function returns int(vector_size) + // For a FullMatrix the function returns int(rows*columns) + template + unsigned int + get_container_size(const std::vector &data) + { + return static_cast(data.size()); + } + + + + template + unsigned int + get_container_size(const Vector &data) + { + return static_cast(data.size()); + } + + + + template + unsigned int + get_container_size(const FullMatrix &data) + { + return static_cast(data.m() * data.n()); + } + + + + // This function initializes and returns a container of type std::vector, + // Vector or FullMatrix. The function does not set the values of the + // elements of the container. The container can store data of a HDF5 dataset + // or a HDF5 selection. The dimensions parameter holds the dimensions of the + // HDF5 dataset or selection. + // + // In the case of a std::vector, the size of the vector will be the total + // size given by dimensions. For example in the case of a dataset of rank 3, + // the dimensions are std::vector{dim_0,dim_1,dim_2}. The size of + // the returned std::vector will be dim_0*dim_1*dim_2 + // + // In the case of a dealii::Vector, the size of the returned dealii::Vector + // will be as well dim_0*dim_1*dim_2 + // + // A FullMatrix can store only data of HDF5 datasets with rank 2. The size + // of the FullMatrix will be FullMatrix(dim_0,dim_2) + template + typename std::enable_if< + std::is_same>::value, + Container>::type + initialize_container(const std::vector &dimensions) + { + return Container(std::accumulate( + dimensions.begin(), dimensions.end(), 1, std::multiplies())); + } + + + + template + typename std::enable_if< + std::is_same>::value, + Container>::type + initialize_container(const std::vector &dimensions) + { + return Container(std::accumulate( + dimensions.begin(), dimensions.end(), 1, std::multiplies())); + } + + + + template + typename std::enable_if< + std::is_same>::value, + Container>::type + initialize_container(const std::vector &dimensions) + { + // If the rank is higher than 2, then remove single-dimensional entries + // from the shape defined by dimensions. This is equivalent to the squeeze + // function of python/numpy. For example the following code would convert + // the vector {1,3,1,2} to {3,2} + std::vector squeezed_dimensions; + + if (dimensions.size() > 2) + { + for (const auto &dimension : dimensions) + { + if (dimension > 1) + squeezed_dimensions.push_back(dimension); + } + } + else + { + squeezed_dimensions = dimensions; + } + + AssertDimension(squeezed_dimensions.size(), 2); + return Container(squeezed_dimensions[0], squeezed_dimensions[1]); + } + + + // This helper function sets the property list of the read and write + // operations of DataSet. A property list has to be created for the MPI + // driver. For the serial driver the default H5P_DEFAULT can be used. In + // addition H5Pset_dxpl_mpio is used to set the MPI mode to collective. + void + set_plist(hid_t &plist, const bool mpi) + { + if (mpi) + { +# ifdef DEAL_II_WITH_MPI + plist = H5Pcreate(H5P_DATASET_XFER); + Assert(plist >= 0, ExcInternalError()); + const herr_t ret = H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE); + (void)ret; + Assert(ret >= 0, ExcInternalError()); +# else + AssertThrow(false, ExcNotImplemented()); +# endif + } + else + { + plist = H5P_DEFAULT; + } + + (void)plist; + (void)mpi; + } + + + // This helper function release the property list handler of the read and + // write operations of DataSet. For the serial there is no need to release + // the property list handler because H5P_DEFAULT has been used. If + // query_io_mode is True then H5Pget_mpio_actual_io_mode and + // H5Pget_mpio_no_collective_cause are used to check if the operation has + // been collective. + void + release_plist(hid_t & plist, + H5D_mpio_actual_io_mode_t &io_mode, + uint32_t & local_no_collective_cause, + uint32_t & global_no_collective_cause, + const bool mpi, + const bool query_io_mode) + { + if (mpi) + { +# ifdef DEAL_II_WITH_MPI + herr_t ret; + (void)ret; + if (query_io_mode) + { + ret = H5Pget_mpio_actual_io_mode(plist, &io_mode); + Assert(ret >= 0, ExcInternalError()); + ret = + H5Pget_mpio_no_collective_cause(plist, + &local_no_collective_cause, + &global_no_collective_cause); + Assert(ret >= 0, ExcInternalError()); + } + ret = H5Pclose(plist); + Assert(ret >= 0, ExcInternalError()); +# else + AssertThrow(false, ExcNotImplemented()); +# endif + } + + (void)plist; + (void)io_mode; + (void)local_no_collective_cause; + (void)global_no_collective_cause; + (void)mpi; + (void)query_io_mode; + } + + + // Convert a HDF5 no_collective_cause code to a human readable string + std::string + no_collective_cause_to_string(const uint32_t no_collective_cause) + { + std::string message; + + auto append_to_message = [&message](const char *p) { + if (message.length() > 0) + message += ", "; + message += p; + }; + + // The first is not a bitmask comparison, the rest are bitmask + // comparisons. + // https://support.hdfgroup.org/HDF5/doc/RM/RM_H5P.html#Property-GetMpioNoCollectiveCause + // See H5Ppublic.h + // Hex codes are used because the HDF5 Group can deprecate some of the + // enum codes. For example the enum code H5D_MPIO_FILTERS is not defined + // in 1.10.2 because it is possible to use compressed datasets with the + // MPI/IO driver. + + // H5D_MPIO_COLLECTIVE + if (no_collective_cause == 0x00) + { + append_to_message("H5D_MPIO_COLLECTIVE"); + } + // H5D_MPIO_SET_INDEPENDENT + if ((no_collective_cause & 0x01) == 0x01) + { + append_to_message("H5D_MPIO_SET_INDEPENDENT"); + } + // H5D_MPIO_DATATYPE_CONVERSION + if ((no_collective_cause & 0x02) == 0x02) + { + append_to_message("H5D_MPIO_DATATYPE_CONVERSION"); + } + // H5D_MPIO_DATA_TRANSFORMS + if ((no_collective_cause & 0x04) == 0x04) + { + append_to_message("H5D_MPIO_DATA_TRANSFORMS"); + } + // H5D_MPIO_NOT_SIMPLE_OR_SCALAR_DATASPACES + if ((no_collective_cause & 0x10) == 0x10) + { + append_to_message("H5D_MPIO_NOT_SIMPLE_OR_SCALAR_DATASPACES"); + } + // H5D_MPIO_NOT_CONTIGUOUS_OR_CHUNKED_DATASET + if ((no_collective_cause & 0x20) == 0x20) + { + append_to_message("H5D_MPIO_NOT_CONTIGUOUS_OR_CHUNKED_DATASET"); + } + // H5D_MPIO_FILTERS + if ((no_collective_cause & 0x40) == 0x40) + { + append_to_message("H5D_MPIO_FILTERS"); + } + return message; + } + } // namespace internal + + + template + T + HDF5Object::get_attribute(const std::string &attr_name) const + { + const std::shared_ptr t_type = internal::get_hdf5_datatype(); + T value; + hid_t attr; + herr_t ret; + + attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT); + Assert(attr >= 0, ExcMessage("Error at H5Aopen")); + (void)ret; + ret = H5Aread(attr, *t_type, &value); + Assert(ret >= 0, ExcMessage("Error at H5Aread")); + (void)ret; + ret = H5Aclose(attr); + Assert(ret >= 0, ExcMessage("Error at H5Aclose")); + + return value; + } + + + + template <> + bool + HDF5Object::get_attribute(const std::string &attr_name) const + { + // The enum field generated by h5py can be casted to int + int int_value; + hid_t attr; + herr_t ret; + + attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT); + Assert(attr >= 0, ExcMessage("Error at H5Aopen")); + ret = H5Aread(attr, H5T_NATIVE_INT, &int_value); + (void)ret; + Assert(ret >= 0, ExcMessage("Error at H5Aread")); + ret = H5Aclose(attr); + (void)ret; + Assert(ret >= 0, ExcMessage("Error at H5Aclose")); + + return (int_value != 0); + } + + + + template <> + std::string + HDF5Object::get_attribute(const std::string &attr_name) const + { + // Reads a UTF8 variable string + // + // code inspired 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. H5Aread will reserve the memory and the + // user has to free the memory. + // + // Todo: + // - Use H5Dvlen_reclaim instead of free + + 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); + ret = H5Tclose(type); + Assert(ret >= 0, ExcInternalError()); + + ret = H5Aclose(attr); + Assert(ret >= 0, ExcInternalError()); + + (void)ret; + return string_value; + } + + + + template + void + HDF5Object::set_attribute(const std::string &attr_name, const T value) + { + hid_t attr; + hid_t aid; + herr_t ret; + + const std::shared_ptr t_type = internal::get_hdf5_datatype(); + + + /* + * Create scalar attribute. + */ + aid = H5Screate(H5S_SCALAR); + Assert(aid >= 0, ExcMessage("Error at H5Screate")); + attr = H5Acreate2(*hdf5_reference, + attr_name.data(), + *t_type, + aid, + H5P_DEFAULT, + H5P_DEFAULT); + Assert(attr >= 0, ExcMessage("Error at H5Acreate2")); + + /* + * Write scalar attribute. + */ + ret = H5Awrite(attr, *t_type, &value); + Assert(ret >= 0, ExcMessage("Error at H5Awrite")); + + ret = H5Sclose(aid); + Assert(ret >= 0, ExcMessage("Error at H5Sclose")); + ret = H5Aclose(attr); + Assert(ret >= 0, ExcMessage("Error at H5Aclose")); + + (void)ret; + } + + + + template <> + void + HDF5Object::set_attribute(const std::string &attr_name, + const std::string value) // NOLINT + { + // Writes a UTF8 variable string + // + // code inspired from + // https://support.hdfgroup.org/ftp/HDF5/examples/misc-examples/vlstratt.c + + hid_t attr; + hid_t aid; + hid_t t_type; + herr_t ret; + + /* Create a datatype to refer to. */ + t_type = H5Tcopy(H5T_C_S1); + Assert(t_type >= 0, ExcInternalError()); + + // Python strings are encoded in UTF8 + ret = H5Tset_cset(t_type, H5T_CSET_UTF8); + Assert(t_type >= 0, ExcInternalError()); + + ret = H5Tset_size(t_type, H5T_VARIABLE); + Assert(ret >= 0, ExcInternalError()); + + /* + * Create scalar attribute. + */ + aid = H5Screate(H5S_SCALAR); + Assert(aid >= 0, ExcMessage("Error at H5Screate")); + attr = H5Acreate2( + *hdf5_reference, attr_name.data(), t_type, aid, H5P_DEFAULT, H5P_DEFAULT); + Assert(attr >= 0, ExcMessage("Error at H5Acreate2")); + + /* + * Write scalar attribute. + * In most of the cases H5Awrite and H5Dwrite take a pointer to the data. + * But in the particular case of a variable length string, H5Awrite takes + * the address of the pointer of the string. + */ + const char *c_string_value = value.c_str(); + ret = H5Awrite(attr, t_type, &c_string_value); + Assert(ret >= 0, ExcInternalError()); + + ret = H5Sclose(aid); + Assert(ret >= 0, ExcMessage("Error at H5Sclose")); + ret = H5Aclose(attr); + Assert(ret >= 0, ExcMessage("Error at H5Aclose")); + + (void)ret; + } + + + + template + Container + DataSet::read() + { + const std::shared_ptr t_type = + internal::get_hdf5_datatype(); + hid_t plist; + herr_t ret; + + Container data = internal::initialize_container(dimensions); + + internal::set_plist(plist, mpi); + + ret = H5Dread(*hdf5_reference, + *t_type, + H5S_ALL, + H5S_ALL, + plist, + make_array_view(data).data()); + Assert(ret >= 0, ExcInternalError()); + + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + (void)ret; + return data; + } + + + + template + Container + DataSet::read_selection(const std::vector &coordinates) + { + Assert(coordinates.size() % rank == 0, + ExcMessage( + "The dimension of coordinates has to be divisible by the rank")); + const std::shared_ptr t_type = + internal::get_hdf5_datatype(); + hid_t plist; + hid_t memory_dataspace; + herr_t ret; + + std::vector data_dimensions{ + static_cast(coordinates.size() / rank)}; + + Container data = internal::initialize_container(data_dimensions); + + memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr); + Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); + ret = H5Sselect_elements(*dataspace, + H5S_SELECT_SET, + data.size(), + coordinates.data()); + Assert(ret >= 0, ExcMessage("Error at H5Sselect_elements")); + + internal::set_plist(plist, mpi); + + ret = H5Dread(*hdf5_reference, + *t_type, + memory_dataspace, + *dataspace, + plist, + make_array_view(data).data()); + Assert(ret >= 0, ExcMessage("Error at H5Dread")); + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + ret = H5Sclose(memory_dataspace); + Assert(ret >= 0, ExcMessage("Error at H5SClose")); + + (void)ret; + return data; + } + + + + template + Container + DataSet::read_hyperslab(const std::vector &offset, + const std::vector &count) + { + const std::shared_ptr t_type = + internal::get_hdf5_datatype(); + hid_t plist; + hid_t memory_dataspace; + herr_t ret; + + // In this particular overload of read_hyperslab the data_dimensions are + // the same as count + std::vector data_dimensions = count; + + Container data = internal::initialize_container(data_dimensions); + + memory_dataspace = + H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr); + Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); + ret = H5Sselect_hyperslab(*dataspace, + H5S_SELECT_SET, + offset.data(), + nullptr, + count.data(), + nullptr); + Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab")); + + internal::set_plist(plist, mpi); + + ret = H5Dread(*hdf5_reference, + *t_type, + memory_dataspace, + *dataspace, + plist, + make_array_view(data).data()); + Assert(ret >= 0, ExcMessage("Error at H5Dread")); + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + ret = H5Sclose(memory_dataspace); + Assert(ret >= 0, ExcMessage("Error at H5SClose")); + + (void)ret; + return data; + } + + + + template + Container + DataSet::read_hyperslab(const std::vector &data_dimensions, + const std::vector &offset, + const std::vector &stride, + const std::vector &count, + const std::vector &block) + { + const std::shared_ptr t_type = + internal::get_hdf5_datatype(); + hid_t plist; + hid_t memory_dataspace; + herr_t ret; + + Container data = internal::initialize_container(data_dimensions); + + memory_dataspace = + H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr); + Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); + ret = H5Sselect_hyperslab(*dataspace, + H5S_SELECT_SET, + offset.data(), + stride.data(), + count.data(), + block.data()); + Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab")); + + internal::set_plist(plist, mpi); + + ret = H5Dread(*hdf5_reference, + *t_type, + memory_dataspace, + *dataspace, + plist, + make_array_view(data).data()); + Assert(ret >= 0, ExcMessage("Error at H5Dread")); + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + ret = H5Sclose(memory_dataspace); + Assert(ret >= 0, ExcMessage("Error at H5SClose")); + + (void)ret; + return data; + } + + + + template + void + DataSet::read_none() + { + const std::shared_ptr t_type = internal::get_hdf5_datatype(); + const std::vector data_dimensions = {0}; + + hid_t memory_dataspace; + hid_t plist; + herr_t ret; + + memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr); + Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); + ret = H5Sselect_none(*dataspace); + Assert(ret >= 0, ExcMessage("H5Sselect_none")); + + internal::set_plist(plist, mpi); + + // The pointer of data can safely be nullptr, see the discussion at the HDF5 + // forum: + // https://forum.hdfgroup.org/t/parallel-i-o-does-not-support-filters-yet/884/17 + ret = H5Dread( + *hdf5_reference, *t_type, memory_dataspace, *dataspace, plist, nullptr); + Assert(ret >= 0, ExcMessage("Error at H5Dread")); + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + ret = H5Sclose(memory_dataspace); + Assert(ret >= 0, ExcMessage("Error at H5SClose")); + + (void)ret; + } + + + + template + void + DataSet::write(const Container &data) + { + AssertDimension(size, internal::get_container_size(data)); + const std::shared_ptr t_type = + internal::get_hdf5_datatype(); + hid_t plist; + herr_t ret; + + internal::set_plist(plist, mpi); + + ret = H5Dwrite(*hdf5_reference, + *t_type, + H5S_ALL, + H5S_ALL, + plist, + make_array_view(data).data()); + Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + (void)ret; + } + + + + template + void + DataSet::write_selection(const Container & data, + const std::vector &coordinates) + { + AssertDimension(coordinates.size(), data.size() * rank); + const std::shared_ptr t_type = + internal::get_hdf5_datatype(); + const std::vector data_dimensions = + internal::get_container_dimensions(data); + + hid_t memory_dataspace; + hid_t plist; + herr_t ret; + + + memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr); + Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); + ret = H5Sselect_elements(*dataspace, + H5S_SELECT_SET, + data.size(), + coordinates.data()); + Assert(ret >= 0, ExcMessage("Error at H5Sselect_elements")); + + internal::set_plist(plist, mpi); + + ret = H5Dwrite(*hdf5_reference, + *t_type, + memory_dataspace, + *dataspace, + plist, + make_array_view(data).data()); + Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + ret = H5Sclose(memory_dataspace); + Assert(ret >= 0, ExcMessage("Error at H5SClose")); + + (void)ret; + } + + + + template + void + DataSet::write_hyperslab(const Container & data, + const std::vector &offset, + const std::vector &count) + { + AssertDimension(std::accumulate(count.begin(), + count.end(), + 1, + std::multiplies()), + internal::get_container_size(data)); + const std::shared_ptr t_type = + internal::get_hdf5_datatype(); + // In this particular overload of write_hyperslab the data_dimensions are + // the same as count + const std::vector &data_dimensions = count; + + hid_t memory_dataspace; + hid_t plist; + herr_t ret; + + memory_dataspace = + H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr); + Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); + ret = H5Sselect_hyperslab(*dataspace, + H5S_SELECT_SET, + offset.data(), + nullptr, + count.data(), + nullptr); + Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab")); + + internal::set_plist(plist, mpi); + + ret = H5Dwrite(*hdf5_reference, + *t_type, + memory_dataspace, + *dataspace, + plist, + make_array_view(data).data()); + Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + ret = H5Sclose(memory_dataspace); + Assert(ret >= 0, ExcMessage("Error at H5Sclose")); + + (void)ret; + } + + + + template + void + DataSet::write_hyperslab(const Container & data, + const std::vector &data_dimensions, + const std::vector &offset, + const std::vector &stride, + const std::vector &count, + const std::vector &block) + { + const std::shared_ptr t_type = + internal::get_hdf5_datatype(); + + hid_t memory_dataspace; + hid_t plist; + herr_t ret; + + memory_dataspace = + H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr); + Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); + ret = H5Sselect_hyperslab(*dataspace, + H5S_SELECT_SET, + offset.data(), + stride.data(), + count.data(), + block.data()); + Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab")); + + internal::set_plist(plist, mpi); + + ret = H5Dwrite(*hdf5_reference, + *t_type, + memory_dataspace, + *dataspace, + plist, + make_array_view(data).data()); + Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + ret = H5Sclose(memory_dataspace); + Assert(ret >= 0, ExcMessage("Error at H5Sclose")); + + (void)ret; + } + + + + template + void + DataSet::write_none() + { + std::shared_ptr t_type = internal::get_hdf5_datatype(); + std::vector data_dimensions = {0}; + + hid_t memory_dataspace; + hid_t plist; + herr_t ret; + + memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr); + Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); + ret = H5Sselect_none(*dataspace); + Assert(ret >= 0, ExcMessage("Error at H5PSselect_none")); + + internal::set_plist(plist, mpi); + + // The pointer of data can safely be nullptr, see the discussion at the HDF5 + // forum: + // https://forum.hdfgroup.org/t/parallel-i-o-does-not-support-filters-yet/884/17 + ret = H5Dwrite( + *hdf5_reference, *t_type, memory_dataspace, *dataspace, plist, nullptr); + Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); + + internal::release_plist(plist, + io_mode, + local_no_collective_cause, + global_no_collective_cause, + mpi, + query_io_mode); + + ret = H5Sclose(memory_dataspace); + Assert(ret >= 0, ExcMessage("Error at H5Sclose")); + + (void)ret; + } + + + + template + DataSet + Group::create_dataset(const std::string & name, + const std::vector &dimensions) const + { + std::shared_ptr t_type = internal::get_hdf5_datatype(); + return {name, *hdf5_reference, dimensions, t_type, mpi}; + } + + + + template + void + Group::write_dataset(const std::string &name, const Container &data) const + { + std::vector dimensions = internal::get_container_dimensions(data); + auto dataset = + create_dataset(name, dimensions); + dataset.write(data); + } } // namespace HDF5 DEAL_II_NAMESPACE_CLOSE diff --git a/source/base/hdf5.cc b/source/base/hdf5.cc index 5e6fab8e0c..ae14b47241 100644 --- a/source/base/hdf5.cc +++ b/source/base/hdf5.cc @@ -17,382 +17,14 @@ #ifdef DEAL_II_WITH_HDF5 -# include # include -# include -# include - # include -# include -# include -# include - DEAL_II_NAMESPACE_OPEN 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. - // std::shared_ptr is used instead of std::unique_ptr because the destructor - // of std::shared_ptr doesn't have to be defined in the template argument. - // In the other hand, the destructor of std::unique has to be defined in the - // template argument. Native types such as H5T_NATIVE_DOUBLE do not - // require a destructor, but compound types such as std::complex - // require a destructor to free the HDF5 resources. - template - std::shared_ptr - get_hdf5_datatype() - { - std::shared_ptr t_type; - if (std::is_same::value) - { - return std::make_shared(H5T_NATIVE_FLOAT); - } - else if (std::is_same::value) - { - return std::make_shared(H5T_NATIVE_DOUBLE); - } - else if (std::is_same::value) - { - return std::make_shared(H5T_NATIVE_INT); - } - else if (std::is_same::value) - { - return std::make_shared(H5T_NATIVE_UINT); - } - else if (std::is_same>::value) - { - t_type = std::shared_ptr(new hid_t, [](hid_t *pointer) { - // Release the HDF5 resource - const herr_t ret = H5Tclose(*pointer); - AssertNothrow(ret >= 0, ExcInternalError()); - (void)ret; - 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. - herr_t ret = H5Tinsert(*t_type, "r", 0, H5T_NATIVE_FLOAT); - Assert(ret >= 0, ExcInternalError()); - ret = H5Tinsert(*t_type, "i", sizeof(float), H5T_NATIVE_FLOAT); - Assert(ret >= 0, ExcInternalError()); - (void)ret; - return t_type; - } - else if (std::is_same>::value) - { - t_type = std::shared_ptr(new hid_t, [](hid_t *pointer) { - // Release the HDF5 resource - const herr_t ret = H5Tclose(*pointer); - AssertNothrow(ret >= 0, ExcInternalError()); - (void)ret; - 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. - herr_t ret = H5Tinsert(*t_type, "r", 0, H5T_NATIVE_DOUBLE); - Assert(ret >= 0, ExcInternalError()); - ret = H5Tinsert(*t_type, "i", sizeof(double), H5T_NATIVE_DOUBLE); - Assert(ret >= 0, ExcInternalError()); - (void)ret; - return t_type; - } - - // The function should not reach this point - Assert(false, ExcInternalError()); - return t_type; - } - - - - // Several HDF5 functions such as H5Screate_simple() require a - // one-dimensional array that specifies the size of each dimension of the - // container, see: - // https://support.hdfgroup.org/HDF5/doc1.8/RM/RM_H5S.html#Dataspace-CreateSimple - // The function get_container_dimensions returns a vector with the - // dimensions of the container. - // For a std::vector the function returns std::vector{vector_size} - // For a Vector the function returns std::vector{vector_size} - // For a FullMatrix the function returns std::vector{rows, columns} - template - std::vector - get_container_dimensions(const std::vector &data) - { - std::vector dimensions = {data.size()}; - return dimensions; - } - - - - template - std::vector - get_container_dimensions(const Vector &data) - { - std::vector dimensions = {data.size()}; - return dimensions; - } - - - - template - std::vector - get_container_dimensions(const FullMatrix &data) - { - std::vector dimensions = {data.m(), data.n()}; - return dimensions; - } - - - - // This function returns the total size of the container. - // For a std::vector the function returns int(vector_size) - // For a Vector the function returns int(vector_size) - // For a FullMatrix the function returns int(rows*columns) - template - unsigned int - get_container_size(const std::vector &data) - { - return static_cast(data.size()); - } - - - - template - unsigned int - get_container_size(const Vector &data) - { - return static_cast(data.size()); - } - - - - template - unsigned int - get_container_size(const FullMatrix &data) - { - return static_cast(data.m() * data.n()); - } - - - - // This function initializes and returns a container of type std::vector, - // Vector or FullMatrix. The function does not set the values of the - // elements of the container. The container can store data of a HDF5 dataset - // or a HDF5 selection. The dimensions parameter holds the dimensions of the - // HDF5 dataset or selection. - // - // In the case of a std::vector, the size of the vector will be the total - // size given by dimensions. For example in the case of a dataset of rank 3, - // the dimensions are std::vector{dim_0,dim_1,dim_2}. The size of - // the returned std::vector will be dim_0*dim_1*dim_2 - // - // In the case of a dealii::Vector, the size of the returned dealii::Vector - // will be as well dim_0*dim_1*dim_2 - // - // A FullMatrix can store only data of HDF5 datasets with rank 2. The size - // of the FullMatrix will be FullMatrix(dim_0,dim_2) - template - typename std::enable_if< - std::is_same>::value, - Container>::type - initialize_container(const std::vector &dimensions) - { - return Container(std::accumulate( - dimensions.begin(), dimensions.end(), 1, std::multiplies())); - } - - - - template - typename std::enable_if< - std::is_same>::value, - Container>::type - initialize_container(const std::vector &dimensions) - { - return Container(std::accumulate( - dimensions.begin(), dimensions.end(), 1, std::multiplies())); - } - - - - template - typename std::enable_if< - std::is_same>::value, - Container>::type - initialize_container(const std::vector &dimensions) - { - // If the rank is higher than 2, then remove single-dimensional entries - // from the shape defined by dimensions. This is equivalent to the squeeze - // function of python/numpy. For example the following code would convert - // the vector {1,3,1,2} to {3,2} - std::vector squeezed_dimensions; - - if (dimensions.size() > 2) - { - for (const auto &dimension : dimensions) - { - if (dimension > 1) - squeezed_dimensions.push_back(dimension); - } - } - else - { - squeezed_dimensions = dimensions; - } - - AssertDimension(squeezed_dimensions.size(), 2); - return Container(squeezed_dimensions[0], squeezed_dimensions[1]); - } - - - // This helper function sets the property list of the read and write - // operations of DataSet. A property list has to be created for the MPI - // driver. For the serial driver the default H5P_DEFAULT can be used. In - // addition H5Pset_dxpl_mpio is used to set the MPI mode to collective. - void - set_plist(hid_t &plist, const bool mpi) - { - if (mpi) - { -# ifdef DEAL_II_WITH_MPI - plist = H5Pcreate(H5P_DATASET_XFER); - Assert(plist >= 0, ExcInternalError()); - const herr_t ret = H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE); - (void)ret; - Assert(ret >= 0, ExcInternalError()); -# else - AssertThrow(false, ExcNotImplemented()); -# endif - } - else - { - plist = H5P_DEFAULT; - } - - (void)plist; - (void)mpi; - } - - - // This helper function release the property list handler of the read and - // write operations of DataSet. For the serial there is no need to release - // the property list handler because H5P_DEFAULT has been used. If - // query_io_mode is True then H5Pget_mpio_actual_io_mode and - // H5Pget_mpio_no_collective_cause are used to check if the operation has - // been collective. - void - release_plist(hid_t & plist, - H5D_mpio_actual_io_mode_t &io_mode, - uint32_t & local_no_collective_cause, - uint32_t & global_no_collective_cause, - const bool mpi, - const bool query_io_mode) - { - if (mpi) - { -# ifdef DEAL_II_WITH_MPI - herr_t ret; - (void)ret; - if (query_io_mode) - { - ret = H5Pget_mpio_actual_io_mode(plist, &io_mode); - Assert(ret >= 0, ExcInternalError()); - ret = - H5Pget_mpio_no_collective_cause(plist, - &local_no_collective_cause, - &global_no_collective_cause); - Assert(ret >= 0, ExcInternalError()); - } - ret = H5Pclose(plist); - Assert(ret >= 0, ExcInternalError()); -# else - AssertThrow(false, ExcNotImplemented()); -# endif - } - - (void)plist; - (void)io_mode; - (void)local_no_collective_cause; - (void)global_no_collective_cause; - (void)mpi; - (void)query_io_mode; - } - - - // Convert a HDF5 no_collective_cause code to a human readable string - std::string - no_collective_cause_to_string(const uint32_t no_collective_cause) - { - std::string message; - - auto append_to_message = [&message](const char *p) { - if (message.length() > 0) - message += ", "; - message += p; - }; - - // The first is not a bitmask comparison, the rest are bitmask - // comparisons. - // https://support.hdfgroup.org/HDF5/doc/RM/RM_H5P.html#Property-GetMpioNoCollectiveCause - // See H5Ppublic.h - // Hex codes are used because the HDF5 Group can deprecate some of the - // enum codes. For example the enum code H5D_MPIO_FILTERS is not defined - // in 1.10.2 because it is possible to use compressed datasets with the - // MPI/IO driver. - - // H5D_MPIO_COLLECTIVE - if (no_collective_cause == 0x00) - { - append_to_message("H5D_MPIO_COLLECTIVE"); - } - // H5D_MPIO_SET_INDEPENDENT - if ((no_collective_cause & 0x01) == 0x01) - { - append_to_message("H5D_MPIO_SET_INDEPENDENT"); - } - // H5D_MPIO_DATATYPE_CONVERSION - if ((no_collective_cause & 0x02) == 0x02) - { - append_to_message("H5D_MPIO_DATATYPE_CONVERSION"); - } - // H5D_MPIO_DATA_TRANSFORMS - if ((no_collective_cause & 0x04) == 0x04) - { - append_to_message("H5D_MPIO_DATA_TRANSFORMS"); - } - // H5D_MPIO_NOT_SIMPLE_OR_SCALAR_DATASPACES - if ((no_collective_cause & 0x10) == 0x10) - { - append_to_message("H5D_MPIO_NOT_SIMPLE_OR_SCALAR_DATASPACES"); - } - // H5D_MPIO_NOT_CONTIGUOUS_OR_CHUNKED_DATASET - if ((no_collective_cause & 0x20) == 0x20) - { - append_to_message("H5D_MPIO_NOT_CONTIGUOUS_OR_CHUNKED_DATASET"); - } - // H5D_MPIO_FILTERS - if ((no_collective_cause & 0x40) == 0x40) - { - append_to_message("H5D_MPIO_FILTERS"); - } - return message; - } - } // namespace internal - - - HDF5Object::HDF5Object(const std::string &name, const bool mpi) : name(name) , mpi(mpi) @@ -400,201 +32,6 @@ namespace HDF5 - template - T - HDF5Object::get_attribute(const std::string &attr_name) const - { - const std::shared_ptr t_type = internal::get_hdf5_datatype(); - T value; - hid_t attr; - herr_t ret; - - attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT); - Assert(attr >= 0, ExcMessage("Error at H5Aopen")); - (void)ret; - ret = H5Aread(attr, *t_type, &value); - Assert(ret >= 0, ExcMessage("Error at H5Aread")); - (void)ret; - ret = H5Aclose(attr); - Assert(ret >= 0, ExcMessage("Error at H5Aclose")); - - return value; - } - - - - template <> - bool - HDF5Object::get_attribute(const std::string &attr_name) const - { - // The enum field generated by h5py can be casted to int - int int_value; - hid_t attr; - herr_t ret; - - attr = H5Aopen(*hdf5_reference, attr_name.data(), H5P_DEFAULT); - Assert(attr >= 0, ExcMessage("Error at H5Aopen")); - ret = H5Aread(attr, H5T_NATIVE_INT, &int_value); - (void)ret; - Assert(ret >= 0, ExcMessage("Error at H5Aread")); - ret = H5Aclose(attr); - (void)ret; - Assert(ret >= 0, ExcMessage("Error at H5Aclose")); - - return (int_value != 0); - } - - - - template <> - std::string - HDF5Object::get_attribute(const std::string &attr_name) const - { - // Reads a UTF8 variable string - // - // code inspired 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. H5Aread will reserve the memory and the - // user has to free the memory. - // - // Todo: - // - Use H5Dvlen_reclaim instead of free - - 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); - ret = H5Tclose(type); - Assert(ret >= 0, ExcInternalError()); - - ret = H5Aclose(attr); - Assert(ret >= 0, ExcInternalError()); - - (void)ret; - return string_value; - } - - - - template - void - HDF5Object::set_attribute(const std::string &attr_name, const T value) - { - hid_t attr; - hid_t aid; - herr_t ret; - - const std::shared_ptr t_type = internal::get_hdf5_datatype(); - - - /* - * Create scalar attribute. - */ - aid = H5Screate(H5S_SCALAR); - Assert(aid >= 0, ExcMessage("Error at H5Screate")); - attr = H5Acreate2(*hdf5_reference, - attr_name.data(), - *t_type, - aid, - H5P_DEFAULT, - H5P_DEFAULT); - Assert(attr >= 0, ExcMessage("Error at H5Acreate2")); - - /* - * Write scalar attribute. - */ - ret = H5Awrite(attr, *t_type, &value); - Assert(ret >= 0, ExcMessage("Error at H5Awrite")); - - ret = H5Sclose(aid); - Assert(ret >= 0, ExcMessage("Error at H5Sclose")); - ret = H5Aclose(attr); - Assert(ret >= 0, ExcMessage("Error at H5Aclose")); - - (void)ret; - } - - - - template <> - void - HDF5Object::set_attribute(const std::string &attr_name, - const std::string value) // NOLINT - { - // Writes a UTF8 variable string - // - // code inspired from - // https://support.hdfgroup.org/ftp/HDF5/examples/misc-examples/vlstratt.c - - hid_t attr; - hid_t aid; - hid_t t_type; - herr_t ret; - - /* Create a datatype to refer to. */ - t_type = H5Tcopy(H5T_C_S1); - Assert(t_type >= 0, ExcInternalError()); - - // Python strings are encoded in UTF8 - ret = H5Tset_cset(t_type, H5T_CSET_UTF8); - Assert(t_type >= 0, ExcInternalError()); - - ret = H5Tset_size(t_type, H5T_VARIABLE); - Assert(ret >= 0, ExcInternalError()); - - /* - * Create scalar attribute. - */ - aid = H5Screate(H5S_SCALAR); - Assert(aid >= 0, ExcMessage("Error at H5Screate")); - attr = H5Acreate2( - *hdf5_reference, attr_name.data(), t_type, aid, H5P_DEFAULT, H5P_DEFAULT); - Assert(attr >= 0, ExcMessage("Error at H5Acreate2")); - - /* - * Write scalar attribute. - * In most of the cases H5Awrite and H5Dwrite take a pointer to the data. - * But in the particular case of a variable length string, H5Awrite takes - * the address of the pointer of the string. - */ - const char *c_string_value = value.c_str(); - ret = H5Awrite(attr, t_type, &c_string_value); - Assert(ret >= 0, ExcInternalError()); - - ret = H5Sclose(aid); - Assert(ret >= 0, ExcMessage("Error at H5Sclose")); - ret = H5Aclose(attr); - Assert(ret >= 0, ExcMessage("Error at H5Aclose")); - - (void)ret; - } - - - std::string HDF5Object::get_name() const { @@ -703,470 +140,6 @@ namespace HDF5 - template - Container - DataSet::read() - { - const std::shared_ptr t_type = - internal::get_hdf5_datatype(); - hid_t plist; - herr_t ret; - - Container data = internal::initialize_container(dimensions); - - internal::set_plist(plist, mpi); - - ret = H5Dread(*hdf5_reference, - *t_type, - H5S_ALL, - H5S_ALL, - plist, - make_array_view(data).data()); - Assert(ret >= 0, ExcInternalError()); - - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - (void)ret; - return data; - } - - - - template - Container - DataSet::read_selection(const std::vector &coordinates) - { - Assert(coordinates.size() % rank == 0, - ExcMessage( - "The dimension of coordinates has to be divisible by the rank")); - const std::shared_ptr t_type = - internal::get_hdf5_datatype(); - hid_t plist; - hid_t memory_dataspace; - herr_t ret; - - std::vector data_dimensions{ - static_cast(coordinates.size() / rank)}; - - Container data = internal::initialize_container(data_dimensions); - - memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr); - Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); - ret = H5Sselect_elements(*dataspace, - H5S_SELECT_SET, - data.size(), - coordinates.data()); - Assert(ret >= 0, ExcMessage("Error at H5Sselect_elements")); - - internal::set_plist(plist, mpi); - - ret = H5Dread(*hdf5_reference, - *t_type, - memory_dataspace, - *dataspace, - plist, - make_array_view(data).data()); - Assert(ret >= 0, ExcMessage("Error at H5Dread")); - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - ret = H5Sclose(memory_dataspace); - Assert(ret >= 0, ExcMessage("Error at H5SClose")); - - (void)ret; - return data; - } - - - - template - Container - DataSet::read_hyperslab(const std::vector &offset, - const std::vector &count) - { - const std::shared_ptr t_type = - internal::get_hdf5_datatype(); - hid_t plist; - hid_t memory_dataspace; - herr_t ret; - - // In this particular overload of read_hyperslab the data_dimensions are - // the same as count - std::vector data_dimensions = count; - - Container data = internal::initialize_container(data_dimensions); - - memory_dataspace = - H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr); - Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); - ret = H5Sselect_hyperslab(*dataspace, - H5S_SELECT_SET, - offset.data(), - nullptr, - count.data(), - nullptr); - Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab")); - - internal::set_plist(plist, mpi); - - ret = H5Dread(*hdf5_reference, - *t_type, - memory_dataspace, - *dataspace, - plist, - make_array_view(data).data()); - Assert(ret >= 0, ExcMessage("Error at H5Dread")); - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - ret = H5Sclose(memory_dataspace); - Assert(ret >= 0, ExcMessage("Error at H5SClose")); - - (void)ret; - return data; - } - - - - template - Container - DataSet::read_hyperslab(const std::vector &data_dimensions, - const std::vector &offset, - const std::vector &stride, - const std::vector &count, - const std::vector &block) - { - const std::shared_ptr t_type = - internal::get_hdf5_datatype(); - hid_t plist; - hid_t memory_dataspace; - herr_t ret; - - Container data = internal::initialize_container(data_dimensions); - - memory_dataspace = - H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr); - Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); - ret = H5Sselect_hyperslab(*dataspace, - H5S_SELECT_SET, - offset.data(), - stride.data(), - count.data(), - block.data()); - Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab")); - - internal::set_plist(plist, mpi); - - ret = H5Dread(*hdf5_reference, - *t_type, - memory_dataspace, - *dataspace, - plist, - make_array_view(data).data()); - Assert(ret >= 0, ExcMessage("Error at H5Dread")); - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - ret = H5Sclose(memory_dataspace); - Assert(ret >= 0, ExcMessage("Error at H5SClose")); - - (void)ret; - return data; - } - - - - template - void - DataSet::read_none() - { - const std::shared_ptr t_type = internal::get_hdf5_datatype(); - const std::vector data_dimensions = {0}; - - hid_t memory_dataspace; - hid_t plist; - herr_t ret; - - memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr); - Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); - ret = H5Sselect_none(*dataspace); - Assert(ret >= 0, ExcMessage("H5Sselect_none")); - - internal::set_plist(plist, mpi); - - // The pointer of data can safely be nullptr, see the discussion at the HDF5 - // forum: - // https://forum.hdfgroup.org/t/parallel-i-o-does-not-support-filters-yet/884/17 - ret = H5Dread( - *hdf5_reference, *t_type, memory_dataspace, *dataspace, plist, nullptr); - Assert(ret >= 0, ExcMessage("Error at H5Dread")); - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - ret = H5Sclose(memory_dataspace); - Assert(ret >= 0, ExcMessage("Error at H5SClose")); - - (void)ret; - } - - - - template - void - DataSet::write(const Container &data) - { - AssertDimension(size, internal::get_container_size(data)); - const std::shared_ptr t_type = - internal::get_hdf5_datatype(); - hid_t plist; - herr_t ret; - - internal::set_plist(plist, mpi); - - ret = H5Dwrite(*hdf5_reference, - *t_type, - H5S_ALL, - H5S_ALL, - plist, - make_array_view(data).data()); - Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - (void)ret; - } - - - - template - void - DataSet::write_selection(const Container & data, - const std::vector &coordinates) - { - AssertDimension(coordinates.size(), data.size() * rank); - const std::shared_ptr t_type = - internal::get_hdf5_datatype(); - const std::vector data_dimensions = - internal::get_container_dimensions(data); - - hid_t memory_dataspace; - hid_t plist; - herr_t ret; - - - memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr); - Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); - ret = H5Sselect_elements(*dataspace, - H5S_SELECT_SET, - data.size(), - coordinates.data()); - Assert(ret >= 0, ExcMessage("Error at H5Sselect_elements")); - - internal::set_plist(plist, mpi); - - ret = H5Dwrite(*hdf5_reference, - *t_type, - memory_dataspace, - *dataspace, - plist, - make_array_view(data).data()); - Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - ret = H5Sclose(memory_dataspace); - Assert(ret >= 0, ExcMessage("Error at H5SClose")); - - (void)ret; - } - - - - template - void - DataSet::write_hyperslab(const Container & data, - const std::vector &offset, - const std::vector &count) - { - AssertDimension(std::accumulate(count.begin(), - count.end(), - 1, - std::multiplies()), - internal::get_container_size(data)); - const std::shared_ptr t_type = - internal::get_hdf5_datatype(); - // In this particular overload of write_hyperslab the data_dimensions are - // the same as count - const std::vector &data_dimensions = count; - - hid_t memory_dataspace; - hid_t plist; - herr_t ret; - - memory_dataspace = - H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr); - Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); - ret = H5Sselect_hyperslab(*dataspace, - H5S_SELECT_SET, - offset.data(), - nullptr, - count.data(), - nullptr); - Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab")); - - internal::set_plist(plist, mpi); - - ret = H5Dwrite(*hdf5_reference, - *t_type, - memory_dataspace, - *dataspace, - plist, - make_array_view(data).data()); - Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - ret = H5Sclose(memory_dataspace); - Assert(ret >= 0, ExcMessage("Error at H5Sclose")); - - (void)ret; - } - - - - template - void - DataSet::write_hyperslab(const Container & data, - const std::vector &data_dimensions, - const std::vector &offset, - const std::vector &stride, - const std::vector &count, - const std::vector &block) - { - const std::shared_ptr t_type = - internal::get_hdf5_datatype(); - - hid_t memory_dataspace; - hid_t plist; - herr_t ret; - - memory_dataspace = - H5Screate_simple(data_dimensions.size(), data_dimensions.data(), nullptr); - Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); - ret = H5Sselect_hyperslab(*dataspace, - H5S_SELECT_SET, - offset.data(), - stride.data(), - count.data(), - block.data()); - Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab")); - - internal::set_plist(plist, mpi); - - ret = H5Dwrite(*hdf5_reference, - *t_type, - memory_dataspace, - *dataspace, - plist, - make_array_view(data).data()); - Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - ret = H5Sclose(memory_dataspace); - Assert(ret >= 0, ExcMessage("Error at H5Sclose")); - - (void)ret; - } - - - - template - void - DataSet::write_none() - { - std::shared_ptr t_type = internal::get_hdf5_datatype(); - std::vector data_dimensions = {0}; - - hid_t memory_dataspace; - hid_t plist; - herr_t ret; - - memory_dataspace = H5Screate_simple(1, data_dimensions.data(), nullptr); - Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple")); - ret = H5Sselect_none(*dataspace); - Assert(ret >= 0, ExcMessage("Error at H5PSselect_none")); - - internal::set_plist(plist, mpi); - - // The pointer of data can safely be nullptr, see the discussion at the HDF5 - // forum: - // https://forum.hdfgroup.org/t/parallel-i-o-does-not-support-filters-yet/884/17 - ret = H5Dwrite( - *hdf5_reference, *t_type, memory_dataspace, *dataspace, plist, nullptr); - Assert(ret >= 0, ExcMessage("Error at H5Dwrite")); - - internal::release_plist(plist, - io_mode, - local_no_collective_cause, - global_no_collective_cause, - mpi, - query_io_mode); - - ret = H5Sclose(memory_dataspace); - Assert(ret >= 0, ExcMessage("Error at H5Sclose")); - - (void)ret; - } - - - void DataSet::set_query_io_mode(const bool new_query_io_mode) { @@ -1365,29 +338,6 @@ namespace HDF5 - template - DataSet - Group::create_dataset(const std::string & name, - const std::vector &dimensions) const - { - std::shared_ptr t_type = internal::get_hdf5_datatype(); - return {name, *hdf5_reference, dimensions, t_type, mpi}; - } - - - - template - void - Group::write_dataset(const std::string &name, const Container &data) const - { - std::vector dimensions = internal::get_container_dimensions(data); - auto dataset = - create_dataset(name, dimensions); - dataset.write(data); - } - - - File::File(const std::string &name, const FileAccessMode mode) : File(name, mode, @@ -1475,7 +425,7 @@ namespace HDF5 (void)ret; (void)mpi_communicator; } -} +} // namespace HDF5 DEAL_II_NAMESPACE_CLOSE