Xdmf io in Idefix (#13)

* [ENH] Added `HighFive` hdf5 io submodule

* [ENH] Added `Xdmf` IO feature

* [RFC] Removed redundant prints

* [DOC] Updated documentation describing `xdmf` dumps

* [RFC] Native `HDF5` dump feature

* [FIX] Fixed datatype of `vect3D` to match dump datatype

* [ENH] Added more info to attributes in hdf5 file

* [FIX] Fixed missing hyperslab select in serial io

* [FIX] Fixed output for all dimensions

* [ADD] Added a 3D test problem with unequal extent

* [FIX] Made `hdf5` includes optional

* [FIX] Optional `hdf5` include

* [FIX] More optional `hdf5` fixes

* [FIX] Linter fixes

* [FIX] Fixed hdf5 optional include pre-processor flags

* [FIX] Fixed hdf5 optional include pre-processor flags

* Used `IDEFIX_ERROR` macro for missing hdf5 linking

* [RFC] Incorporated changes suggested by @glesur

Author: dutta-alankar

.gitignore

@@ -25,6 +25,8 @@ Makefile.local
 # test artifacts
 test/**/*.o
 test/**/*.vtk
+test/**/*.h5
+test/**/*.xmf
 test/**/*.dmp
 test/**/*.idfx
 test/**/*.tmp

CMakeLists.txt

@@ -82,6 +82,19 @@ if(Idefix_MPI)
   target_link_libraries(idefix MPI::MPI_CXX)
 endif()
 
+if(Idefix_HDF5)
+  add_compile_definitions("WITH_HDF5")
+  if(Idefix_MPI)
+    set(HDF5_PREFER_PARALLEL TRUE)
+  endif()
+  find_package(HDF5 REQUIRED)
+  target_link_libraries(idefix "${HDF5_LIBRARIES}")
+  target_include_directories(idefix PUBLIC "${HDF5_INCLUDE_DIRS}")
+  message(STATUS "XDMF (hdf5+xmf) dumps enabled")
+else()
+  set(Idefix_HDF5 OFF)
+endif()
+
 if(Idefix_DEBUG)
   add_compile_definitions("DEBUG")
   set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -O0 -g")
@@ -166,11 +179,12 @@ target_link_libraries(idefix Kokkos::kokkos)
 
 message(STATUS "Idefix final configuration")
 if(Idefix_EVOLVE_VECTOR_POTENTIAL)
-  message(STATUS "    MHD: ${Idefix_MHD} (Vector potential)")
+  message(STATUS "    MHD:  ${Idefix_MHD} (Vector potential)")
 else()
-  message(STATUS "    MHD: ${Idefix_MHD}")
+  message(STATUS "    MHD:  ${Idefix_MHD}")
 endif()
-message(STATUS "    MPI: ${Idefix_MPI}")
+message(STATUS "    MPI:  ${Idefix_MPI}")
+message(STATUS "    HDF5: ${Idefix_HDF5}")
 message(STATUS "    Reconstruction: ${Idefix_RECONSTRUCTION}")
 message(STATUS "    Precision: ${Idefix_PRECISION}")
 message(STATUS "    Version: ${GIT_SHA1}")

README.md

@@ -57,6 +57,9 @@ cmake $IDEFIX_DIR
 
 Several options can be enabled from the command line (a complete list is available with `cmake $IDEFIX_DIR -LH`). For instance: `-DIdefix_RECONSTRUCTION=Parabolic` (enable PPM reconstruction), `-DIdefix_MPI=ON` (enable mpi), `-DKokkos_ENABLE_OPENMP=ON` (enable openmp parallelisation), etc... For more complex target architectures, it is recommended to use cmake GUI launching `ccmake $IDEFIX_DIR` in place of `cmake` and then switching on the required options.
 
+Optional xdmf(hdf5+xmf) file dumping feature has been added to `Idefix`. This uses either serial or parallel implementation of `hdf5` library which needs to be made available. These xdmf file pairs can be easily visualized in `ParaView` or `VisIt` by loading the `xmf` files. The `hdf5` files can also be loaded easily in `python` (using `h5py`) for post-processing and post-run analysis. One can turn on `xdmf` data dumps by using `-DIdefix_HDF5=ON`. The `[Output]` block of `.ini` file is checked during runtime for a `xdmf` entry whih controls the frequency of xdmf file dumps during code execution.
+<!-- TODO: HDF5 Chunking and Compression filters -->
+
 One can then compile the code:
 
 ```shell

src/output/CMakeLists.txt

@@ -6,3 +6,9 @@ target_sources(idefix
   PUBLIC ${CMAKE_CURRENT_LIST_DIR}/vtk.cpp
   PUBLIC ${CMAKE_CURRENT_LIST_DIR}/vtk.hpp
   )
+if(Idefix_HDF5)
+target_sources(idefix
+  PUBLIC ${CMAKE_CURRENT_LIST_DIR}/xdmf.cpp
+  PUBLIC ${CMAKE_CURRENT_LIST_DIR}/xdmf.hpp
+  )
+endif()

src/output/output.cpp

@@ -7,7 +7,6 @@
 
 #include "output.hpp"
 
-
 Output::Output(Input &input, DataBlock &data) {
   idfx::pushRegion("Output::Output");
   // initialise the output objects for each format
@@ -25,6 +24,23 @@ Output::Output(Input &input, DataBlock &data) {
   }
   vtk.Init(input,data); // Always initialised in case of emergency vtk output
 
+  // Initialise xdmf outputs
+  if(input.CheckEntry("Output","xdmf")>0) {
+    xdmfPeriod = input.Get<real>("Output","xdmf",0);
+    if(xdmfPeriod>=0.0) {  // backward compatibility (negative value means no file)
+      xdmfLast = data.t - xdmfPeriod; // write something in the next CheckForWrite()
+      #ifdef WITH_HDF5
+      xdmfEnabled = true;
+      #else
+      xdmfEnabled = false;
+      IDEFIX_ERROR("Attention: HDF5 library not linked when building Idefix!");
+      #endif
+    }
+  }
+  #ifdef WITH_HDF5
+  xdmf.Init(input,data); // Always initialised in case of emergency xdmf output
+  #endif
+
   // intialise dump outputs
   if(input.CheckEntry("Output","dmp")>0) {
     dumpPeriod = input.Get<real>("Output","dmp",0);
@@ -116,6 +132,37 @@ int Output::CheckForWrites(DataBlock &data) {
       }
     }
   }
+  #ifdef WITH_HDF5
+  // Do we need a XDMF output?
+  if(xdmfEnabled) {
+    if(data.t >= xdmfLast + xdmfPeriod) {
+      elapsedTime -= timer.seconds();
+      if(userDefVariablesEnabled) {
+        if(haveUserDefVariablesFunc) {
+          // Call user-def function to fill the userdefined variable arrays
+          idfx::pushRegion("UserDef::User-defined variables function");
+          userDefVariablesFunc(data, userDefVariables);
+          idfx::popRegion();
+        } else {
+          IDEFIX_ERROR("Cannot output user-defined variables without "
+                        "enrollment of your user-defined variables function");
+        }
+      }
+      xdmfLast += xdmfPeriod;
+      xdmf.Write(data, *this);
+      nfiles++;
+      elapsedTime += timer.seconds();
+
+      // Check if our next predicted output should already have happened
+      if((xdmfLast+xdmfPeriod <= data.t) && xdmfPeriod>0.0) {
+        // Move forward xdmfLast
+        while(xdmfLast <= data.t - xdmfPeriod) {
+          xdmfLast += xdmfPeriod;
+        }
+      }
+    }
+  }
+  #endif
 
   // Do we need an analysis ?
   if(analysisEnabled) {
@@ -184,6 +231,29 @@ void Output::ForceWriteVtk(DataBlock &data) {
   idfx::popRegion();
 }
 
+#ifdef WITH_HDF5
+void Output::ForceWriteXdmf(DataBlock &data) {
+  idfx::pushRegion("Output::ForceWriteXdmf");
+
+  if(!forceNoWrite) {
+    if(userDefVariablesEnabled) {
+        if(haveUserDefVariablesFunc) {
+          // Call user-def function to fill the userdefined variable arrays
+          idfx::pushRegion("UserDef::User-defined variables function");
+          userDefVariablesFunc(data, userDefVariables);
+          idfx::popRegion();
+        } else {
+          IDEFIX_ERROR("Cannot output user-defined variables without "
+                        "enrollment of your user-defined variables function");
+        }
+      }
+      xdmfLast += xdmfPeriod;
+      xdmf.Write(data, *this);
+  }
+  idfx::popRegion();
+}
+#endif
+
 void Output::EnrollAnalysis(AnalysisFunc myFunc) {
   idfx::pushRegion("Output::EnrollAnalysis");
   if(!analysisEnabled) {

src/output/output.hpp

@@ -13,6 +13,9 @@
 #include "input.hpp"
 #include "dataBlock.hpp"
 #include "vtk.hpp"
+#ifdef WITH_HDF5
+#include "xdmf.hpp"
+#endif
 #include "dump.hpp"
 
 
@@ -25,13 +28,19 @@ using UserDefVariablesFunc = void (*) (DataBlock &, UserDefVariablesContainer &)
 class Output {
   friend class Dump;    // Allow dump to have R/W access to private variables
   friend class Vtk;     // Allow VTK to have access to user-defined variables
+  #ifdef WITH_HDF5
+  friend class Xdmf;    // Allow XDMF to have access to user-defined variables
+  #endif
   friend class DumpImage; // Allow dumpimag to have access to dump API
  public:
   Output(Input &, DataBlock &);           // Create Output Object
   int CheckForWrites(DataBlock &);        // Check if outputs are needed at this stage
   void RestartFromDump(DataBlock &, int);  // Restart from a dump file.
   void ForceWriteDump(DataBlock &);            // Force write dumps (needed during an abort)
   void ForceWriteVtk(DataBlock &);            // Force write vtks
+  #ifdef WITH_HDF5
+  void ForceWriteXdmf(DataBlock &);          // Force write xdmfs
+  #endif
   void ResetTimer();                      // Reset internal timer
   double GetTimer();
   void EnrollAnalysis(AnalysisFunc);
@@ -40,7 +49,9 @@ class Output {
  private:
   Vtk vtk;          // local instance of Vtk class
   Dump dump;        // local instance of Dump class
-
+  #ifdef WITH_HDF5
+  Xdmf xdmf;        // local instance of Xdmf class
+  #endif
   bool forceNoWrite = false;    //< explicitely disable all writes
   bool vtkEnabled = false;
   real vtkPeriod = 0.0;   // periodicity of vtk outputs
@@ -50,6 +61,10 @@ class Output {
   real dumpPeriod = 0.0;
   real dumpLast = 0.0;
 
+  bool xdmfEnabled = false;
+  real xdmfPeriod = 0.0;   // periodicity of xdmf outputs
+  real xdmfLast = 0.0;
+
   bool analysisEnabled = false;
   real analysisPeriod = 0.0;
   real analysisLast = 0.0;
@@ -65,6 +80,4 @@ class Output {
   Kokkos::Timer timer;
   double elapsedTime{0.0};
 };
-
-
 #endif // OUTPUT_OUTPUT_HPP_

src/output/vtk.cpp

@@ -26,7 +26,7 @@
   #endif
 #endif
 
-// Whether of not we write the time in the VTK file
+// Whether or not we write the time in the VTK file
 #define WRITE_TIME
 
 void Vtk::WriteHeaderString(const char* header, IdfxFileHandler fvtk) {
@@ -216,7 +216,7 @@ void Vtk::Init(Input &input, DataBlock &datain) {
 
 #endif
 
-  // Creat MPI view when using MPI I/O
+  // Create MPI view when using MPI I/O
 #ifdef WITH_MPI
   int start[3];
   int size[3];
@@ -346,7 +346,7 @@ void Vtk::WriteHeader(IdfxFileHandler fvtk, real time) {
 * Write VTK header in parallel or serial mode.
 *
 * \param [in]  fvtk  pointer to file
-* \param [in]  grid  pointer to an array of Grid structures
+* \param [in]  time  time of the simulation in code units
 *
 * \todo  Write the grid using several processors.
 *********************************************************************** */