From 2d637ba15848aa8147be354c8c785e230f3fe6c5 Mon Sep 17 00:00:00 2001
From: Daniel Garcia-Sanchez <daniel.garcia-sanchez@insp.upmc.fr>
Date: Sun, 19 Aug 2018 13:53:26 +0200
Subject: [PATCH] Improve write_hyperslab()

---
 include/deal.II/base/hdf5.h |  64 ++++++++-
 source/base/hdf5.cc         | 267 +++++++++++++++++++++++++++++++++---
 2 files changed, 306 insertions(+), 25 deletions(-)

diff --git a/include/deal.II/base/hdf5.h b/include/deal.II/base/hdf5.h
index 768135eeb9..3705a01983 100644
--- a/include/deal.II/base/hdf5.h
+++ b/include/deal.II/base/hdf5.h
@@ -29,7 +29,7 @@
 DEAL_II_NAMESPACE_OPEN
 
 /**
- * Namespace containing the HDF5 wrappers.
+ * Namespace containing the HDF5 interface.
  *
  * This set of classes can be used to store the data in the hdf5 file format.
  *
@@ -247,10 +247,56 @@ namespace HDF5
                     const std::vector<hsize_t> coordinates);
 
     /**
-     * Writes data to a subset of the dataset. T can be double, int, unsigned
-     * int, bool or std::complex<double>.
+     * Writes a data hyperslab to the dataset. The parameters are summarized
+     * below:
+     *  - Offset: The starting location for the hyperslab.
+     *  - Count: The number of elements to select along each dimension.
+     */
+    // clang-format off
+    /**
+     *
+     * stride and block are set to NULL. For more complex hyperslabs see
+     * write_hyperslab(const Container<T> &data, const std::vector<hsize_t> &data_dimensions, const std::vector<hsize_t> &offset, const std::vector<hsize_t> &stride, const std::vector<hsize_t> &count, const std::vector<hsize_t> &block).
+     */
+    // clang-format on
+    /**
+     *
+     * See the <a
+     * href="https://support.hdfgroup.org/HDF5/doc/UG/HDF5_Users_Guide-Responsive%20HTML5/HDF5_Users_Guide/Dataspaces/HDF5_Dataspaces_and_Partial_I_O.htm?rhtocid=7.2#TOC_7_4_Dataspaces_and_Databc-6">"Dataspaces
+     * and Data Transfer"</a>  section in the HDF5 User's Guide. See as well the
+     * <a
+     * href="https://support.hdfgroup.org/HDF5/doc1.8/RM/RM_H5S.html#Dataspace-SelectHyperslab">"H5Sselect_hyperslab
+     * definition"</a>.
+     *
+     * Datatype conversion takes place at the time of a read or write and is
+     * automatic. See the <a
+     * href="https://support.hdfgroup.org/HDF5/doc/UG/HDF5_Users_Guide-Responsive%20HTML5/index.html#t=HDF5_Users_Guide%2FDatatypes%2FHDF5_Datatypes.htm%23TOC_6_10_Data_Transferbc-26&rhtocid=6.5_2">"Data
+     * Transfer: Datatype Conversion and Selection"</a>  section in the HDF5
+     * User's Guide.
+     */
+    template <template <class...> class Container, typename T>
+    void
+    write_hyperslab(const Container<T> &        data,
+                    const std::vector<hsize_t> &offset,
+                    const std::vector<hsize_t> &count);
+
+    /**
+     * Writes a data hyperslab to the dataset. The parameters are summarized
+     * below:
+     *  - Dataset_dimensions: the dimensions of the data memory block.
+     *  - Offset: The starting location for the hyperslab.
+     *  - Stride: The number of elements to separate each element or block to be
+     * selected.
+     *  - Count: The number of elements or blocks to select along each
+     * dimension.
+     *  - Block: The size of the block selected from the dataspace.
      *
-     * The selected elements form a hyperslab in the dataset.
+     * See the <a
+     * href="https://support.hdfgroup.org/HDF5/doc/UG/HDF5_Users_Guide-Responsive%20HTML5/HDF5_Users_Guide/Dataspaces/HDF5_Dataspaces_and_Partial_I_O.htm?rhtocid=7.2#TOC_7_4_Dataspaces_and_Databc-6">"Dataspaces
+     * and Data Transfer"</a>  section in the HDF5 User's Guide. See as well the
+     * <a
+     * href="https://support.hdfgroup.org/HDF5/doc1.8/RM/RM_H5S.html#Dataspace-SelectHyperslab">"H5Sselect_hyperslab
+     * definition"</a>.
      *
      * Datatype conversion takes place at the time of a read or write and is
      * automatic. See the <a
@@ -260,9 +306,13 @@ namespace HDF5
      */
     template <template <class...> class Container, typename T>
     void
-    write_hyperslab(const Container<T> &       data,
-                    const std::vector<hsize_t> offset,
-                    const std::vector<hsize_t> count);
+    write_hyperslab(const Container<T> &        data,
+                    const std::vector<hsize_t> &data_dimensions,
+                    const std::vector<hsize_t> &offset,
+                    const std::vector<hsize_t> &stride,
+                    const std::vector<hsize_t> &count,
+                    const std::vector<hsize_t> &block);
+
 
     /**
      * This function does not write any data, but can contribute to a collective
diff --git a/source/base/hdf5.cc b/source/base/hdf5.cc
index 69b53fa2fa..519166e14e 100644
--- a/source/base/hdf5.cc
+++ b/source/base/hdf5.cc
@@ -598,57 +598,131 @@ namespace HDF5
 
   template <template <class...> class Container, typename T>
   void
-  DataSet::write_hyperslab(const Container<T> &       data,
-                           const std::vector<hsize_t> offset,
-                           const std::vector<hsize_t> count)
+  DataSet::write_hyperslab(const Container<T> &        data,
+                           const std::vector<hsize_t> &offset,
+                           const std::vector<hsize_t> &count)
   {
     AssertDimension(std::accumulate(count.begin(),
                                     count.end(),
                                     1,
                                     std::multiplies<unsigned int>()),
-                    data.size());
-    std::shared_ptr<hid_t> t_type          = internal::get_hdf5_datatype<T>();
-    std::vector<hsize_t>   data_dimensions = {data.size()};
+                    internal::get_container_size(data));
+    std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<T>();
+    // In this particular overload of write_hyperslab the data_dimensions are
+    // the same as count
+    std::vector<hsize_t> data_dimensions = count;
 
     hid_t  memory_dataspace;
     hid_t  plist;
     herr_t ret;
 
-
-    memory_dataspace = H5Screate_simple(1, data_dimensions.data(), NULL);
-    H5Sselect_hyperslab(
+    memory_dataspace =
+      H5Screate_simple(data_dimensions.size(), data_dimensions.data(), NULL);
+    Assert(memory_dataspace >= 0, ExcMessage("Error at H5Screate_simple"));
+    ret = H5Sselect_hyperslab(
       *dataspace, H5S_SELECT_SET, offset.data(), NULL, count.data(), NULL);
+    Assert(ret >= 0, ExcMessage("Error at H5Sselect_hyperslab"));
 
     if (mpi)
       {
         plist = H5Pcreate(H5P_DATASET_XFER);
-        H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE);
+        Assert(plist >= 0, ExcMessage("Error at H5Pcreate"));
+        ret = H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE);
+        Assert(ret >= 0, ExcMessage("Error at H5Pset_dxpl_mpio"));
       }
     else
       {
         plist = H5P_DEFAULT;
       }
-    H5Dwrite(*hdf5_reference,
-             *t_type,
-             memory_dataspace,
-             *dataspace,
-             plist,
-             internal::get_container_const_pointer<Container, T>(data));
+    ret = H5Dwrite(*hdf5_reference,
+                   *t_type,
+                   memory_dataspace,
+                   *dataspace,
+                   plist,
+                   internal::get_container_const_pointer<Container, T>(data));
+    Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
 
     if (mpi)
       {
         if (_check_io_mode)
           {
             ret = H5Pget_mpio_actual_io_mode(plist, &_io_mode);
-            Assert(ret >= 0, ExcInternalError());
+            Assert(ret >= 0, ExcMessage("Error at H5Pget_mpio_actual_io_mode"));
             ret = H5Pget_mpio_no_collective_cause(plist,
                                                   &_local_no_collective_cause,
                                                   &_global_no_collective_cause);
-            Assert(ret >= 0, ExcInternalError());
+            Assert(ret >= 0,
+                   ExcMessage("Error at H5Pget_mpio_no_collective_cause"));
           }
         H5Pclose(plist);
+        Assert(ret >= 0, ExcMessage("Error at H5Pclose"));
       }
     H5Sclose(memory_dataspace);
+    Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
+  }
+
+  template <template <class...> class Container, typename T>
+  void
+  DataSet::write_hyperslab(const Container<T> &        data,
+                           const std::vector<hsize_t> &data_dimensions,
+                           const std::vector<hsize_t> &offset,
+                           const std::vector<hsize_t> &stride,
+                           const std::vector<hsize_t> &count,
+                           const std::vector<hsize_t> &block)
+  {
+    std::shared_ptr<hid_t> t_type = internal::get_hdf5_datatype<T>();
+
+    hid_t  memory_dataspace;
+    hid_t  plist;
+    herr_t ret;
+
+    memory_dataspace =
+      H5Screate_simple(data_dimensions.size(), data_dimensions.data(), NULL);
+    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"));
+
+    if (mpi)
+      {
+        plist = H5Pcreate(H5P_DATASET_XFER);
+        Assert(plist >= 0, ExcMessage("Error at H5Pcreate"));
+        ret = H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE);
+        Assert(ret >= 0, ExcMessage("Error at H5Pset_dxpl_mpio"));
+      }
+    else
+      {
+        plist = H5P_DEFAULT;
+      }
+    ret = H5Dwrite(*hdf5_reference,
+                   *t_type,
+                   memory_dataspace,
+                   *dataspace,
+                   plist,
+                   internal::get_container_const_pointer<Container, T>(data));
+    Assert(ret >= 0, ExcMessage("Error at H5Dwrite"));
+
+    if (mpi)
+      {
+        if (_check_io_mode)
+          {
+            ret = H5Pget_mpio_actual_io_mode(plist, &_io_mode);
+            Assert(ret >= 0, ExcMessage("Error at H5Pget_mpio_actual_io_mode"));
+            ret = H5Pget_mpio_no_collective_cause(plist,
+                                                  &_local_no_collective_cause,
+                                                  &_global_no_collective_cause);
+            Assert(ret >= 0,
+                   ExcMessage("Error at H5Pget_mpio_no_collective_cause"));
+          }
+        H5Pclose(plist);
+        Assert(ret >= 0, ExcMessage("Error at H5Pclose"));
+      }
+    H5Sclose(memory_dataspace);
+    Assert(ret >= 0, ExcMessage("Error at H5Sclose"));
   }
 
   template <typename T>
@@ -1209,6 +1283,163 @@ namespace HDF5
     const std::vector<unsigned int> &data,
     const std::vector<hsize_t>       coordinates);
 
+  template void
+  DataSet::write_hyperslab<std::vector, float>(
+    const std::vector<float> &  data,
+    const std::vector<hsize_t> &offset,
+    const std::vector<hsize_t> &count);
+  template void
+  DataSet::write_hyperslab<std::vector, double>(
+    const std::vector<double> & data,
+    const std::vector<hsize_t> &offset,
+    const std::vector<hsize_t> &count);
+  template void
+  DataSet::write_hyperslab<std::vector, long double>(
+    const std::vector<long double> &data,
+    const std::vector<hsize_t> &    offset,
+    const std::vector<hsize_t> &    count);
+  template void
+  DataSet::write_hyperslab<std::vector, std::complex<float>>(
+    const std::vector<std::complex<float>> &data,
+    const std::vector<hsize_t> &            offset,
+    const std::vector<hsize_t> &            count);
+  template void
+  DataSet::write_hyperslab<std::vector, std::complex<double>>(
+    const std::vector<std::complex<double>> &data,
+    const std::vector<hsize_t> &             offset,
+    const std::vector<hsize_t> &             count);
+  template void
+  DataSet::write_hyperslab<std::vector, std::complex<long double>>(
+    const std::vector<std::complex<long double>> &data,
+    const std::vector<hsize_t> &                  offset,
+    const std::vector<hsize_t> &                  count);
+  template void
+  DataSet::write_hyperslab<std::vector, int>(const std::vector<int> &    data,
+                                             const std::vector<hsize_t> &offset,
+                                             const std::vector<hsize_t> &count);
+  template void
+  DataSet::write_hyperslab<std::vector, unsigned int>(
+    const std::vector<unsigned int> &data,
+    const std::vector<hsize_t> &     offset,
+    const std::vector<hsize_t> &     count);
+  template void
+  DataSet::write_hyperslab<FullMatrix, float>(
+    const FullMatrix<float> &   data,
+    const std::vector<hsize_t> &offset,
+    const std::vector<hsize_t> &count);
+  template void
+  DataSet::write_hyperslab<FullMatrix, double>(
+    const FullMatrix<double> &  data,
+    const std::vector<hsize_t> &offset,
+    const std::vector<hsize_t> &count);
+  template void
+  DataSet::write_hyperslab<FullMatrix, std::complex<float>>(
+    const FullMatrix<std::complex<float>> &data,
+    const std::vector<hsize_t> &           offset,
+    const std::vector<hsize_t> &           count);
+  template void
+  DataSet::write_hyperslab<FullMatrix, std::complex<double>>(
+    const FullMatrix<std::complex<double>> &data,
+    const std::vector<hsize_t> &            offset,
+    const std::vector<hsize_t> &            count);
+
+  template void
+  DataSet::write_hyperslab<std::vector, float>(
+    const std::vector<float> &  data,
+    const std::vector<hsize_t> &data_dimensions,
+    const std::vector<hsize_t> &offset,
+    const std::vector<hsize_t> &stride,
+    const std::vector<hsize_t> &count,
+    const std::vector<hsize_t> &block);
+  template void
+  DataSet::write_hyperslab<std::vector, double>(
+    const std::vector<double> & data,
+    const std::vector<hsize_t> &data_dimensions,
+    const std::vector<hsize_t> &offset,
+    const std::vector<hsize_t> &stride,
+    const std::vector<hsize_t> &count,
+    const std::vector<hsize_t> &block);
+  template void
+  DataSet::write_hyperslab<std::vector, long double>(
+    const std::vector<long double> &data,
+    const std::vector<hsize_t> &    data_dimensions,
+    const std::vector<hsize_t> &    offset,
+    const std::vector<hsize_t> &    stride,
+    const std::vector<hsize_t> &    count,
+    const std::vector<hsize_t> &    block);
+  template void
+  DataSet::write_hyperslab<std::vector, std::complex<float>>(
+    const std::vector<std::complex<float>> &data,
+    const std::vector<hsize_t> &            data_dimensions,
+    const std::vector<hsize_t> &            offset,
+    const std::vector<hsize_t> &            stride,
+    const std::vector<hsize_t> &            count,
+    const std::vector<hsize_t> &            block);
+  template void
+  DataSet::write_hyperslab<std::vector, std::complex<double>>(
+    const std::vector<std::complex<double>> &data,
+    const std::vector<hsize_t> &             data_dimensions,
+    const std::vector<hsize_t> &             offset,
+    const std::vector<hsize_t> &             stride,
+    const std::vector<hsize_t> &             count,
+    const std::vector<hsize_t> &             block);
+  template void
+  DataSet::write_hyperslab<std::vector, std::complex<long double>>(
+    const std::vector<std::complex<long double>> &data,
+    const std::vector<hsize_t> &                  data_dimensions,
+    const std::vector<hsize_t> &                  offset,
+    const std::vector<hsize_t> &                  stride,
+    const std::vector<hsize_t> &                  count,
+    const std::vector<hsize_t> &                  block);
+  template void
+  DataSet::write_hyperslab<std::vector, int>(
+    const std::vector<int> &    data,
+    const std::vector<hsize_t> &data_dimensions,
+    const std::vector<hsize_t> &offset,
+    const std::vector<hsize_t> &stride,
+    const std::vector<hsize_t> &count,
+    const std::vector<hsize_t> &block);
+  template void
+  DataSet::write_hyperslab<std::vector, unsigned int>(
+    const std::vector<unsigned int> &data,
+    const std::vector<hsize_t> &     data_dimensions,
+    const std::vector<hsize_t> &     offset,
+    const std::vector<hsize_t> &     stride,
+    const std::vector<hsize_t> &     count,
+    const std::vector<hsize_t> &     block);
+  template void
+  DataSet::write_hyperslab<FullMatrix, float>(
+    const FullMatrix<float> &   data,
+    const std::vector<hsize_t> &data_dimensions,
+    const std::vector<hsize_t> &offset,
+    const std::vector<hsize_t> &stride,
+    const std::vector<hsize_t> &count,
+    const std::vector<hsize_t> &block);
+  template void
+  DataSet::write_hyperslab<FullMatrix, double>(
+    const FullMatrix<double> &  data,
+    const std::vector<hsize_t> &data_dimensions,
+    const std::vector<hsize_t> &offset,
+    const std::vector<hsize_t> &stride,
+    const std::vector<hsize_t> &count,
+    const std::vector<hsize_t> &block);
+  template void
+  DataSet::write_hyperslab<FullMatrix, std::complex<float>>(
+    const FullMatrix<std::complex<float>> &data,
+    const std::vector<hsize_t> &           data_dimensions,
+    const std::vector<hsize_t> &           offset,
+    const std::vector<hsize_t> &           stride,
+    const std::vector<hsize_t> &           count,
+    const std::vector<hsize_t> &           block);
+  template void
+  DataSet::write_hyperslab<FullMatrix, std::complex<double>>(
+    const FullMatrix<std::complex<double>> &data,
+    const std::vector<hsize_t> &            data_dimensions,
+    const std::vector<hsize_t> &            offset,
+    const std::vector<hsize_t> &            stride,
+    const std::vector<hsize_t> &            count,
+    const std::vector<hsize_t> &            block);
+
   template void
   DataSet::write_none<float>();
   template void
-- 
2.39.5