]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Enable output of compressed hdf5 files 14958/head
authorChristoph Schmidt <christoph.schmidt@tum.de>
Thu, 16 Feb 2023 17:58:15 +0000 (18:58 +0100)
committerChristoph Schmidt <christoph.schmidt@tum.de>
Tue, 28 Mar 2023 15:13:30 +0000 (17:13 +0200)
- Hdf5Flags are introduced to enable setting the desired DataOutBase::CompressionLevel
- HDF5Flags are added to the corresponding write_hdf5_parallel and do_write_hdf5 methods
- if deal.II is build with zlib support, all necessary calls to perform compression using
  zlib (i.e. setting the compression level, setting the chunk size) are shielded behind
  'ifdef DEAL_II_WITH_ZLIB'
- if deal.II is build without zlib support, but compression_level is not set to no_compression
  an exception is raised.

cmake/config/template-arguments.in
doc/news/changes/incompatibilities/20230324Schmidt [new file with mode: 0644]
include/deal.II/base/data_out_base.h
source/base/data_out_base.cc
source/base/data_out_base.inst.in
tests/mpi/data_out_hdf5_02.cc

index 96654f2d3244615ee587e138ce3467bee39701be..f7d58dec3d1d8e527ed670c1ddfb336e7719c6c8 100644 (file)
@@ -266,7 +266,7 @@ SYM_RANKS := { 2; 4 }
 
 // Flags that are allowed in DataOutInterface::set_flags
 OUTPUT_FLAG_TYPES := { DXFlags; UcdFlags; GnuplotFlags; PovrayFlags; EpsFlags;
-                       GmvFlags; TecplotFlags; VtkFlags; SvgFlags;
+                       GmvFlags; Hdf5Flags; TecplotFlags; VtkFlags; SvgFlags;
                        Deal_II_IntermediateFlags }
 
 // CGAL Kernels
diff --git a/doc/news/changes/incompatibilities/20230324Schmidt b/doc/news/changes/incompatibilities/20230324Schmidt
new file mode 100644 (file)
index 0000000..50eabf4
--- /dev/null
@@ -0,0 +1,3 @@
+Changed: The default behavior of the hdf5 output is changed from "not compressed" to "compressed" (default compression_level: best_speed).
+<br>
+(Christoph Schmidt, 2023/03/24)
index 7624876bece749f477b4207a214f143134ad7a37..5479903734369211009241107b15e8406ce1bfcc 100644 (file)
@@ -1078,6 +1078,23 @@ namespace DataOutBase
   struct GmvFlags : public OutputFlagsBase<GmvFlags>
   {};
 
+  /**
+   * Flags controlling the details of output in HDF5 format.
+   *
+   * @ingroup output
+   */
+  struct Hdf5Flags : public OutputFlagsBase<Hdf5Flags>
+  {
+    /**
+     * Flag determining the compression level at which zlib, if available, is
+     * run. The default is <tt>best_speed</tt>.
+     */
+    DataOutBase::CompressionLevel compression_level;
+
+    explicit Hdf5Flags(
+      const CompressionLevel compression_level = CompressionLevel::best_speed);
+  };
+
   /**
    * Flags controlling the details of output in Tecplot format.
    *
@@ -2396,6 +2413,7 @@ namespace DataOutBase
   void
   write_hdf5_parallel(const std::vector<Patch<dim, spacedim>> &patches,
                       const DataOutFilter &                    data_filter,
+                      const DataOutBase::Hdf5Flags &           flags,
                       const std::string &                      filename,
                       const MPI_Comm &                         comm);
 
@@ -2410,6 +2428,7 @@ namespace DataOutBase
   void
   write_hdf5_parallel(const std::vector<Patch<dim, spacedim>> &patches,
                       const DataOutFilter &                    data_filter,
+                      const DataOutBase::Hdf5Flags &           flags,
                       const bool                               write_mesh_file,
                       const std::string &                      mesh_filename,
                       const std::string &solution_filename,
@@ -3158,6 +3177,12 @@ private:
    */
   DataOutBase::GmvFlags gmv_flags;
 
+  /**
+   * Flags to be used upon output of hdf5 data in one space dimension. Can be
+   * changed by using the <tt>set_flags</tt> function.
+   */
+  DataOutBase::Hdf5Flags hdf5_flags;
+
   /**
    * Flags to be used upon output of Tecplot data in one space dimension. Can
    * be changed by using the <tt>set_flags</tt> function.
index 8365a7b10524b626b3b1d90d3b2191518a234eea..f1ab7a59153cb4169e6f3d6806eec7cd91963419 100644 (file)
@@ -2405,6 +2405,10 @@ namespace DataOutBase
   }
 
 
+  Hdf5Flags::Hdf5Flags(const CompressionLevel compression_level)
+    : compression_level(compression_level)
+  {}
+
 
   TecplotFlags::TecplotFlags(const char *zone_name, const double solution_time)
     : zone_name(zone_name)
@@ -8584,6 +8588,7 @@ namespace
   void
   do_write_hdf5(const std::vector<DataOutBase::Patch<dim, spacedim>> &patches,
                 const DataOutBase::DataOutFilter &data_filter,
+                const DataOutBase::Hdf5Flags &    flags,
                 const bool                        write_mesh_file,
                 const std::string &               mesh_filename,
                 const std::string &               solution_filename,
@@ -8591,7 +8596,7 @@ namespace
   {
     hid_t h5_mesh_file_id = -1, h5_solution_file_id, file_plist_id, plist_id;
     hid_t node_dataspace, node_dataset, node_file_dataspace,
-      node_memory_dataspace;
+      node_memory_dataspace, node_dataset_id;
     hid_t cell_dataspace, cell_dataset, cell_file_dataspace,
       cell_memory_dataspace;
     hid_t pt_data_dataspace, pt_data_dataset, pt_data_file_dataspace,
@@ -8617,6 +8622,10 @@ namespace
     status = H5Pset_fapl_mpio(file_plist_id, comm, MPI_INFO_NULL);
     AssertThrow(status >= 0, ExcIO());
 #    endif
+#  endif
+    // if zlib support is disabled flags are unused
+#  ifndef DEAL_II_WITH_ZLIB
+    (void)flags;
 #  endif
 
     // Compute the global total number of nodes/cells and determine the offset
@@ -8685,13 +8694,20 @@ namespace
                                  node_dataspace,
                                  H5P_DEFAULT);
 #  else
-        node_dataset    = H5Dcreate(h5_mesh_file_id,
+        node_dataset_id = H5Pcreate(H5P_DATASET_CREATE);
+#    ifdef DEAL_II_WITH_ZLIB
+        H5Pset_deflate(node_dataset_id,
+                       get_zlib_compression_level(flags.compression_level));
+        H5Pset_chunk(node_dataset_id, 2, node_ds_dim);
+#    endif
+        node_dataset = H5Dcreate(h5_mesh_file_id,
                                  "nodes",
                                  H5T_NATIVE_DOUBLE,
                                  node_dataspace,
                                  H5P_DEFAULT,
-                                 H5P_DEFAULT,
+                                 node_dataset_id,
                                  H5P_DEFAULT);
+        H5Pclose(node_dataset_id);
 #  endif
         AssertThrow(node_dataset >= 0, ExcIO());
 #  if H5Gcreate_vers == 1
@@ -8701,13 +8717,20 @@ namespace
                                  cell_dataspace,
                                  H5P_DEFAULT);
 #  else
-        cell_dataset    = H5Dcreate(h5_mesh_file_id,
+        node_dataset_id = H5Pcreate(H5P_DATASET_CREATE);
+#    ifdef DEAL_II_WITH_ZLIB
+        H5Pset_deflate(node_dataset_id,
+                       get_zlib_compression_level(flags.compression_level));
+        H5Pset_chunk(node_dataset_id, 2, cell_ds_dim);
+#    endif
+        cell_dataset = H5Dcreate(h5_mesh_file_id,
                                  "cells",
                                  H5T_NATIVE_UINT,
                                  cell_dataspace,
                                  H5P_DEFAULT,
-                                 H5P_DEFAULT,
+                                 node_dataset_id,
                                  H5P_DEFAULT);
+        H5Pclose(node_dataset_id);
 #  endif
         AssertThrow(cell_dataset >= 0, ExcIO());
 
@@ -8836,13 +8859,20 @@ namespace
                                     pt_data_dataspace,
                                     H5P_DEFAULT);
 #  else
+        node_dataset_id = H5Pcreate(H5P_DATASET_CREATE);
+#    ifdef DEAL_II_WITH_ZLIB
+        H5Pset_deflate(node_dataset_id,
+                       get_zlib_compression_level(flags.compression_level));
+        H5Pset_chunk(node_dataset_id, 2, node_ds_dim);
+#    endif
         pt_data_dataset = H5Dcreate(h5_solution_file_id,
                                     vector_name.c_str(),
                                     H5T_NATIVE_DOUBLE,
                                     pt_data_dataspace,
                                     H5P_DEFAULT,
-                                    H5P_DEFAULT,
+                                    node_dataset_id,
                                     H5P_DEFAULT);
+        H5Pclose(node_dataset_id);
 #  endif
         AssertThrow(pt_data_dataset >= 0, ExcIO());
 
@@ -9033,7 +9063,8 @@ DataOutInterface<dim, spacedim>::write_hdf5_parallel(
   const std::string &               filename,
   const MPI_Comm &                  comm) const
 {
-  DataOutBase::write_hdf5_parallel(get_patches(), data_filter, filename, comm);
+  DataOutBase::write_hdf5_parallel(
+    get_patches(), data_filter, hdf5_flags, filename, comm);
 }
 
 
@@ -9049,6 +9080,7 @@ DataOutInterface<dim, spacedim>::write_hdf5_parallel(
 {
   DataOutBase::write_hdf5_parallel(get_patches(),
                                    data_filter,
+                                   hdf5_flags,
                                    write_mesh_file,
                                    mesh_filename,
                                    solution_filename,
@@ -9062,10 +9094,12 @@ void
 DataOutBase::write_hdf5_parallel(
   const std::vector<Patch<dim, spacedim>> &patches,
   const DataOutBase::DataOutFilter &       data_filter,
+  const DataOutBase::Hdf5Flags &           flags,
   const std::string &                      filename,
   const MPI_Comm &                         comm)
 {
-  write_hdf5_parallel(patches, data_filter, true, filename, filename, comm);
+  write_hdf5_parallel(
+    patches, data_filter, flags, true, filename, filename, comm);
 }
 
 
@@ -9075,6 +9109,7 @@ void
 DataOutBase::write_hdf5_parallel(
   const std::vector<Patch<dim, spacedim>> &patches,
   const DataOutBase::DataOutFilter &       data_filter,
+  const DataOutBase::Hdf5Flags &           flags,
   const bool                               write_mesh_file,
   const std::string &                      mesh_filename,
   const std::string &                      solution_filename,
@@ -9091,6 +9126,7 @@ DataOutBase::write_hdf5_parallel(
   // the now unused function arguments
   (void)patches;
   (void)data_filter;
+  (void)flags;
   (void)write_mesh_file;
   (void)mesh_filename;
   (void)solution_filename;
@@ -9136,6 +9172,7 @@ DataOutBase::write_hdf5_parallel(
     {
       do_write_hdf5<dim, spacedim>(patches,
                                    data_filter,
+                                   flags,
                                    write_mesh_file,
                                    mesh_filename,
                                    solution_filename,
@@ -9242,6 +9279,8 @@ DataOutInterface<dim, spacedim>::set_flags(const FlagType &flags)
     eps_flags = *reinterpret_cast<const DataOutBase::EpsFlags *>(&flags);
   else if (typeid(flags) == typeid(gmv_flags))
     gmv_flags = *reinterpret_cast<const DataOutBase::GmvFlags *>(&flags);
+  else if (typeid(flags) == typeid(hdf5_flags))
+    hdf5_flags = *reinterpret_cast<const DataOutBase::Hdf5Flags *>(&flags);
   else if (typeid(flags) == typeid(tecplot_flags))
     tecplot_flags =
       *reinterpret_cast<const DataOutBase::TecplotFlags *>(&flags);
@@ -9311,6 +9350,10 @@ DataOutInterface<dim, spacedim>::declare_parameters(ParameterHandler &prm)
   DataOutBase::GmvFlags::declare_parameters(prm);
   prm.leave_subsection();
 
+  prm.enter_subsection("HDF5 output parameters");
+  DataOutBase::Hdf5Flags::declare_parameters(prm);
+  prm.leave_subsection();
+
   prm.enter_subsection("Tecplot output parameters");
   DataOutBase::TecplotFlags::declare_parameters(prm);
   prm.leave_subsection();
@@ -9359,6 +9402,10 @@ DataOutInterface<dim, spacedim>::parse_parameters(ParameterHandler &prm)
   gmv_flags.parse_parameters(prm);
   prm.leave_subsection();
 
+  prm.enter_subsection("HDF5 output parameters");
+  hdf5_flags.parse_parameters(prm);
+  prm.leave_subsection();
+
   prm.enter_subsection("Tecplot output parameters");
   tecplot_flags.parse_parameters(prm);
   prm.leave_subsection();
@@ -9385,6 +9432,7 @@ DataOutInterface<dim, spacedim>::memory_consumption() const
           MemoryConsumption::memory_consumption(povray_flags) +
           MemoryConsumption::memory_consumption(eps_flags) +
           MemoryConsumption::memory_consumption(gmv_flags) +
+          MemoryConsumption::memory_consumption(hdf5_flags) +
           MemoryConsumption::memory_consumption(tecplot_flags) +
           MemoryConsumption::memory_consumption(vtk_flags) +
           MemoryConsumption::memory_consumption(svg_flags) +
index 35077900a8fb0818d429a11ad5e31a6d71fc63d3..7df4210f7858c1d1151e1c9dfe5ab0bdd9cef2f4 100644 (file)
@@ -222,10 +222,11 @@ for (deal_II_dimension : OUTPUT_DIMENSIONS;
       template void
       write_hdf5_parallel(
         const std::vector<Patch<deal_II_dimension, deal_II_space_dimension>>
-          &                  patches,
-        const DataOutFilter &data_filter,
-        const std::string &  filename,
-        const MPI_Comm &     comm);
+          &                           patches,
+        const DataOutFilter &         data_filter,
+        const DataOutBase::Hdf5Flags &flags,
+        const std::string &           filename,
+        const MPI_Comm &              comm);
 
       template void
       write_filtered_data(
index 05003706ddf87d02d923febb0caf5f653d12edf3..39cd1eefeaa4eeb5d1a4cebbfff32940853ec81f 100644 (file)
@@ -54,10 +54,11 @@ check()
 
   std::string output_basename = std::to_string(dim) + std::to_string(spacedim);
 
-  DataOutBase::write_hdf5_parallel(patches,
-                                   data_filter,
-                                   output_basename + ".h5",
-                                   MPI_COMM_WORLD);
+  DataOutBase::Hdf5Flags hdf5Flags;
+  hdf5Flags.compression_level = DataOutBase::CompressionLevel::no_compression;
+
+  DataOutBase::write_hdf5_parallel(
+    patches, data_filter, hdf5Flags, output_basename + ".h5", MPI_COMM_WORLD);
 
   const double current_time = 0.0;
   XDMFEntry    entry(output_basename + ".h5",

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.