diff --git a/BaseLib/cpp23.h b/BaseLib/cpp23.h
new file mode 100644
index 0000000000000000000000000000000000000000..3d7b5680da5eaf78e7be161c9cfcdbdbd6953fea
--- /dev/null
+++ b/BaseLib/cpp23.h
@@ -0,0 +1,33 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#pragma once
+
+#include <type_traits>
+
+// ToDo (tm) remove with c++23
+
+namespace BaseLib
+{
+/**
+ *  \brief Converts an enumeration to its underlying type.
+ *  \param e Enumeration value to convert
+ *  \return The integer value of the underlying type of Enum, converted from e.
+ *
+ * https://en.cppreference.com/w/cpp/utility/to_underlying
+ * Implementation from Scott Meyers : Modern effective c++ , Item 10 Scoped enums
+ */
+
+template <typename E>
+constexpr auto to_underlying(E e) noexcept
+{
+    return static_cast<std::underlying_type_t<E>>(e);
+}
+}  // namespace BaseLib
\ No newline at end of file
diff --git a/CMakeLists.txt b/CMakeLists.txt
index f254ad85d301c82e07283cf869da4e6cd7ac7b6e..26b83698bd052b2073a2bccc3b05827855a13e2f 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -102,7 +102,6 @@ if(OGS_USE_PETSC AND MSVC)
     )
 endif()
 option(OGS_USE_NETCDF "Add NetCDF support." OFF)
-option(OGS_USE_XDMF "Add Xdmf file IO support" OFF)
 
 # Eigen
 option(OGS_USE_EIGEN_UNSUPPORTED "Use Eigen unsupported modules" ON)
@@ -195,7 +194,7 @@ include(packaging/Pack)
 # ---- Subdirectories ----
 include_directories(${PROJECT_SOURCE_DIR})
 # xdmfdiff
-if(OGS_USE_XDMF AND OGS_BUILD_TESTING)
+if(OGS_BUILD_TESTING)
     add_subdirectory(Tests/xdmfdiff)
 endif()
 
diff --git a/CMakePresets.json b/CMakePresets.json
index 051de3e0e66d90c656463f9a6446fa73cae0e9d1..6a2b0af8b59d3d95c143db6eea7d6539a29b5426 100644
--- a/CMakePresets.json
+++ b/CMakePresets.json
@@ -134,8 +134,7 @@
       "name": "_all",
       "hidden": true,
       "cacheVariables": {
-        "OGS_USE_MFRONT": "ON",
-        "OGS_USE_XDMF": "ON"
+        "OGS_USE_MFRONT": "ON"
       }
     },
     {
diff --git a/MeshLib/CMakeLists.txt b/MeshLib/CMakeLists.txt
index 69d1cc9aa16a7c8b43aa5ff8f0f0ecb7ec66cf85..fce297b5f462b5d09964b3547cdfa2ed818a8974 100644
--- a/MeshLib/CMakeLists.txt
+++ b/MeshLib/CMakeLists.txt
@@ -13,14 +13,13 @@ append_source_files(SOURCES IO)
 append_source_files(SOURCES IO/Legacy)
 append_source_files(SOURCES IO/VtkIO)
 append_source_files(SOURCES Utils)
-if(OGS_USE_XDMF)
-    append_source_files(SOURCES IO/XDMF)
-    if(OGS_USE_PETSC)
-        append_source_files(SOURCES IO/XDMF/mpi)
-    else()
-        append_source_files(SOURCES IO/XDMF/posix)
-    endif()
+append_source_files(SOURCES IO/XDMF)
+if(OGS_USE_PETSC)
+    append_source_files(SOURCES IO/XDMF/mpi)
+else()
+    append_source_files(SOURCES IO/XDMF/posix)
 endif()
+
 append_source_files(SOURCES MeshQuality)
 append_source_files(SOURCES Vtk)
 
@@ -50,5 +49,4 @@ target_compile_definitions(
     MeshLib
     PUBLIC
         $<$<AND:$<BOOL:$ENV{VTK_USE_64BIT_IDS}>,$<NOT:$<BOOL:VTK_ADDED>>>:VTK_USE_64BIT_IDS>
-    PRIVATE $<$<BOOL:${OGS_USE_XDMF}>:OGS_USE_XDMF>
 )
diff --git a/MeshLib/IO/XDMF/HdfWriter.cpp b/MeshLib/IO/XDMF/HdfWriter.cpp
index 76b1f2812731ce5eb7d8423814d9ab48a361c8b9..fcdd27f465740088d11b9429dec2ccee3f61099b 100644
--- a/MeshLib/IO/XDMF/HdfWriter.cpp
+++ b/MeshLib/IO/XDMF/HdfWriter.cpp
@@ -76,6 +76,9 @@ static hid_t createDataSet(
     std::vector<Hdf5DimType> time_data_global_dims =
         prependDimension(1, max_dims);
 
+    std::vector<Hdf5DimType> time_data_chunk_dims =
+        prependDimension(1, chunk_dims);
+
     hid_t fspace =
         H5Screate_simple(time_dim_local_size, time_data_global_dims.data(),
                          time_max_dims.data());
@@ -85,7 +88,7 @@ static hid_t createDataSet(
     assert(dcpl >= 0);
 
     hid_t status =
-        H5Pset_chunk(dcpl, time_dim_local_size, time_data_global_dims.data());
+        H5Pset_chunk(dcpl, chunk_dims.size() + 1, time_data_chunk_dims.data());
     if (status < 0)
     {
         OGS_FATAL("H5Pset_layout failed for data set: {:s}.", dataset_name);
@@ -163,13 +166,14 @@ namespace MeshLib::IO
 HdfWriter::HdfWriter(std::vector<HdfData> constant_attributes,
                      std::vector<HdfData>
                          variable_attributes,
-                     int const step,
+                     int const initial_step,
                      std::filesystem::path const& filepath,
                      bool const use_compression)
     : _variable_attributes(std::move(variable_attributes)),
       _hdf5_filepath(filepath),
       _use_compression(checkCompression() && use_compression),
-      _file(createFile(filepath))
+      _file(createFile(filepath)),
+      _output_step(initial_step)
 {
     _group = H5Gcreate2(_file, "data", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
 
@@ -183,7 +187,7 @@ HdfWriter::HdfWriter(std::vector<HdfData> constant_attributes,
         writeDataSet(attribute.data_start, attribute.data_type,
                      attribute.data_space, attribute.offsets,
                      attribute.file_space, attribute.chunk_space,
-                     attribute.name, step, dataset);
+                     attribute.name, _output_step, dataset);
         return dataset;
     };
 
@@ -199,6 +203,7 @@ HdfWriter::HdfWriter(std::vector<HdfData> constant_attributes,
         // datasets are kept open
         _datasets.insert({attribute.name, dataset});
     }
+    _output_step++;
 }
 
 HdfWriter::~HdfWriter()
@@ -211,7 +216,7 @@ HdfWriter::~HdfWriter()
     H5Fclose(_file);
 }
 
-bool HdfWriter::writeStep(int const step) const
+void HdfWriter::writeStep()
 {
     for (auto const& attribute : _variable_attributes)
     {
@@ -221,12 +226,11 @@ bool HdfWriter::writeStep(int const step) const
             OGS_FATAL("Writing HDF5 Dataset: {:s} failed.", attribute.name);
         }
 
-        writeDataSet(attribute.data_start, attribute.data_type,
-                     attribute.data_space, attribute.offsets,
-                     attribute.file_space, attribute.chunk_space,
-                     attribute.name, step, _datasets.at(attribute.name));
+        writeDataSet(
+            attribute.data_start, attribute.data_type, attribute.data_space,
+            attribute.offsets, attribute.file_space, attribute.chunk_space,
+            attribute.name, _output_step, _datasets.at(attribute.name));
     }
-
-    return true;
+    _output_step++;
 }
 }  // namespace MeshLib::IO
diff --git a/MeshLib/IO/XDMF/HdfWriter.h b/MeshLib/IO/XDMF/HdfWriter.h
index 6f62beba496c103716e7b4273c9cef5b56146237..23edc9a50ccfa8a7a8d3f9c59add8a5bf0c60fe7 100644
--- a/MeshLib/IO/XDMF/HdfWriter.h
+++ b/MeshLib/IO/XDMF/HdfWriter.h
@@ -33,14 +33,15 @@ public:
      * attributes
      * @param variable_attributes vector of variable attributes (each attribute
      * is a OGS mesh property
-     * @param step number of the step (temporal collection)
+     * @param initial_step number of the step (temporal collection), usually 0,
+     * greater 0 with continuation of simulation
      * @param filepath absolute or relative filepath to the hdf5 file
      * @param use_compression if true gzip compression is enabled
      */
     HdfWriter(std::vector<HdfData> constant_attributes,
               std::vector<HdfData>
                   variable_attributes,
-              int const step,
+              int const initial_step,
               std::filesystem::path const& filepath,
               bool const use_compression);
 
@@ -50,9 +51,8 @@ public:
      * the data holder to not change during writing and HdfData given to
      * constructor to be still valid
      * @param step number of the step (temporal collection)
-     * @return true = success, false = error
      */
-    bool writeStep(int step) const;
+    void writeStep();
     ~HdfWriter();
 
 private:
@@ -62,5 +62,6 @@ private:
     hid_t const _file;
     hid_t _group;
     std::map<std::string, hid_t> _datasets;
+    int _output_step;
 };
 }  // namespace MeshLib::IO
\ No newline at end of file
diff --git a/MeshLib/IO/XDMF/MeshPropertyDataType.h b/MeshLib/IO/XDMF/MeshPropertyDataType.h
index 35d32bf1fa2b1527106a7ef8662d7110df0c0d45..e792809d9624d16e93d7d0172232dfbfb1f8f2f0 100644
--- a/MeshLib/IO/XDMF/MeshPropertyDataType.h
+++ b/MeshLib/IO/XDMF/MeshPropertyDataType.h
@@ -15,7 +15,7 @@
 // TODO (tm) If used on several other places move definition of propertyVector
 enum class MeshPropertyDataType
 {
-    unknown,
+    unknown=0,
     float64,
     float32,
     int32,
@@ -23,5 +23,6 @@ enum class MeshPropertyDataType
     uint32,
     uint64,
     int8,
-    uint8
+    uint8,
+    enum_length
 };
\ No newline at end of file
diff --git a/MeshLib/IO/XDMF/Xdmf3Writer.cpp b/MeshLib/IO/XDMF/Xdmf3Writer.cpp
deleted file mode 100644
index 91556c378c9757f8173be3d8792c5e10f19e8510..0000000000000000000000000000000000000000
--- a/MeshLib/IO/XDMF/Xdmf3Writer.cpp
+++ /dev/null
@@ -1,182 +0,0 @@
-/**
- * \file
- * \copyright
- * Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
- *            Distributed under a Modified BSD License.
- *              See accompanying file LICENSE.txt or
- *              http://www.opengeosys.org/project/license
- */
-
-#include "Xdmf3Writer.h"
-
-#include <Xdmf.hpp>
-#include <XdmfAttribute.hpp>
-#include <XdmfDomain.hpp>
-#include <XdmfGeometryType.hpp>
-#include <XdmfGridCollection.hpp>
-#include <XdmfGridCollectionType.hpp>
-#include <XdmfHDF5Controller.hpp>
-#include <XdmfHeavyDataDescription.hpp>
-#include <XdmfInformation.hpp>
-#include <XdmfTopologyType.hpp>
-#include <XdmfUnstructuredGrid.hpp>
-#include <XdmfWriter.hpp>
-#include <string>
-
-#include "InfoLib/GitInfo.h"
-
-using namespace MeshLib::IO;
-using namespace std::string_literals;
-
-boost::shared_ptr<const XdmfAttributeCenter> elemTypeOGS2XDMF(
-    MeshLib::MeshItemType const elem_type)
-{
-    std::map<MeshLib::MeshItemType,
-             boost::shared_ptr<const XdmfAttributeCenter>>
-        mesh_item_type_ogs2xdmf = {
-            {MeshLib::MeshItemType::Cell, XdmfAttributeCenter::Cell()},
-            {MeshLib::MeshItemType::Edge, XdmfAttributeCenter::Edge()},
-            {MeshLib::MeshItemType::Face, XdmfAttributeCenter::Face()},
-            {MeshLib::MeshItemType::Node, XdmfAttributeCenter::Node()},
-            {MeshLib::MeshItemType::IntegrationPoint,
-             XdmfAttributeCenter::Other()}};
-
-    return mesh_item_type_ogs2xdmf.at(elem_type);
-}
-
-static std::string getDataSection(std::string const& name)
-{
-    return "data/"s + name;
-}
-
-static std::vector<XdmfDimType> prependDimension(
-    XdmfDimType const prepend_value, std::vector<XdmfDimType> const& dimensions)
-{
-    std::vector<XdmfDimType> dims = {prepend_value};
-    dims.insert(dims.end(), dimensions.begin(), dimensions.end());
-    return dims;
-}
-
-static boost::shared_ptr<XdmfGeometry> getLightGeometry(
-    std::string const& hdf5filename, XdmfData const& geometry)
-{
-    auto xdmf_geometry = XdmfGeometry::New();
-    xdmf_geometry->setType(XdmfGeometryType::XYZ());
-    boost::shared_ptr<XdmfHDF5Controller> geometry_controller =
-        XdmfHDF5Controller::New(hdf5filename,
-                                getDataSection("geometry"),
-                                XdmfArrayType::Float64(),
-                                geometry.starts,
-                                geometry.strides,
-                                geometry.global_block_dims,
-                                geometry.global_block_dims);
-    xdmf_geometry->setHeavyDataController(geometry_controller);
-    return xdmf_geometry;
-}
-
-static boost::shared_ptr<XdmfTopology> getLightTopology(
-    std::string const& hdf5filename, XdmfData const& topology)
-{
-    auto xdmf_topology = XdmfTopology::New();
-    xdmf_topology->setType(XdmfTopologyType::Mixed());
-    auto topology_controller =
-        XdmfHDF5Controller::New(hdf5filename,
-                                getDataSection("topology"),
-                                XdmfArrayType::Int32(),
-                                topology.starts,
-                                topology.strides,
-                                topology.global_block_dims,
-                                topology.global_block_dims);
-    xdmf_topology->setHeavyDataController(topology_controller);
-    return xdmf_topology;
-}
-
-static boost::shared_ptr<XdmfAttribute> getLightAttribute(
-    std::string const& hdf5filename, int const step, XdmfData const& attribute)
-{
-    std::vector<XdmfDimType> starts = prependDimension(step, attribute.starts);
-    std::vector<XdmfDimType> strides = prependDimension(1, attribute.strides);
-    std::vector<XdmfDimType> global_block_dims =
-        prependDimension(1, attribute.global_block_dims);
-    std::vector<XdmfDimType> all_global_block_dims =
-        prependDimension(step + 1, attribute.global_block_dims);
-
-    auto const attribute_controller =
-        XdmfHDF5Controller::New(hdf5filename,
-                                getDataSection(attribute.name),
-                                attribute.data_type,
-                                starts,
-                                strides,
-                                global_block_dims,
-                                all_global_block_dims);
-
-    auto const xdmf_attribute = XdmfAttribute::New();
-    auto const center = elemTypeOGS2XDMF(*(attribute.attribute_center));
-    xdmf_attribute->setCenter(center);
-    xdmf_attribute->setName(attribute.name);
-    xdmf_attribute->setHeavyDataController(attribute_controller);
-    return xdmf_attribute;
-}
-
-namespace MeshLib::IO
-{
-Xdmf3Writer::Xdmf3Writer(XdmfData const& geometry, XdmfData const& topology,
-                         std::vector<XdmfData> constant_attributes,
-                         std::vector<XdmfData> variable_attributes,
-                         std::filesystem::path const& filepath,
-                         int const time_step)
-    : _variable_attributes(std::move(variable_attributes)),
-      _hdf5filename(filepath.stem().string() + ".h5")
-{
-    _initial_geometry = getLightGeometry(_hdf5filename, geometry);
-    _initial_topology = getLightTopology(_hdf5filename, topology);
-
-    std::transform(
-        constant_attributes.begin(), constant_attributes.end(),
-        std::back_inserter(_constant_attributes),
-        [&](XdmfData const& attribute) -> boost::shared_ptr<XdmfAttribute> {
-            return getLightAttribute(_hdf5filename, time_step, attribute);
-        });
-
-    _writer = XdmfWriter::New(filepath.string());
-    _writer->setMode(XdmfWriter::DistributedHeavyData);
-
-    auto version = XdmfInformation::New();
-    version->setKey(GitInfoLib::GitInfo::OGS_VERSION);
-    version->setValue(GitInfoLib::GitInfo::ogs_version);
-
-    auto grid_collection = XdmfGridCollection::New();
-    grid_collection->setType(XdmfGridCollectionType::Temporal());
-
-    _root = XdmfDomain::New();
-    _root->insert(version);
-    _root->insert(grid_collection);
-}
-
-Xdmf3Writer::~Xdmf3Writer()
-{
-    _root->accept(_writer);
-}
-
-void Xdmf3Writer::writeStep(int const time_step, double const time)
-{
-    auto grid = XdmfUnstructuredGrid::New();
-    grid->setGeometry(_initial_geometry);
-    grid->setTopology(_initial_topology);
-
-    for (auto const& constant_attribute : _constant_attributes)
-    {
-        grid->insert(constant_attribute);
-    }
-
-    grid->setTime(XdmfTime::New(time));
-
-    for (auto const& attribute : _variable_attributes)
-    {
-        grid->insert(getLightAttribute(_hdf5filename, time_step, attribute));
-    }
-
-    auto grid_collection = _root->getGridCollection(0);
-    grid_collection->insert(grid);
-}
-}  // namespace MeshLib::IO
diff --git a/MeshLib/IO/XDMF/Xdmf3Writer.h b/MeshLib/IO/XDMF/Xdmf3Writer.h
deleted file mode 100644
index 5489609d5fcacd57380a1af5f2c081093a9f16b7..0000000000000000000000000000000000000000
--- a/MeshLib/IO/XDMF/Xdmf3Writer.h
+++ /dev/null
@@ -1,77 +0,0 @@
-/**
- * \file
- * \author Tobias Meisel
- * \date   2020-11-13
- * \brief  XdmfWriter which takes contiguous data and writes 1 xdmf + 1 hdf file
- * \copyright
- * Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
- *            Distributed under a Modified BSD License.
- *              See accompanying file LICENSE.txt or
- *              http://www.opengeosys.org/project/license
- */
-
-#pragma once
-
-#include <filesystem.h>
-
-#include <boost/shared_ptr.hpp>
-
-#include "MeshLib/Mesh.h"
-#include "XdmfData.h"
-
-class XdmfAttribute;
-class XdmfGridCollection;
-class XdmfTopology;
-class XdmfGeometry;
-class XdmfWriter;
-class XdmfDomain;
-
-namespace MeshLib::IO
-{
-struct Geometry;
-struct Topology;
-}  // namespace MeshLib::IO
-
-namespace MeshLib::IO
-{
-class Xdmf3Writer final
-{
-public:
-    /**
-     * \brief Write xdmf and h5 file with geometry and topology data.
-     * The data itself is held by a structure outside of this class.
-     * The writer assumes a constant data holder (data itself can change, memory
-     * address is the same)
-     * @param geometry contains meta data to coordinates
-     * @param topology contains meta data cells
-     * @param constant_attributes vector of constant attributes (each attribute
-     * is a OGS mesh property)
-     * @param variable_attributes vector of variable attributes (each attribute
-     * is a OGS mesh property
-     * @param filepath absolute or relative filepath to the hdf5 file
-     * @param time_step number of the step (temporal collection)
-     */
-    Xdmf3Writer(XdmfData const& geometry, XdmfData const& topology,
-                std::vector<XdmfData> constant_attributes,
-                std::vector<XdmfData> variable_attributes,
-                std::filesystem::path const& filepath, int time_step);
-    /**
-     * \brief Write attribute data that has modified to previous time step or
-     * initial
-     * @param time_step number of the step (temporal collection)
-     * @param time time value of the current time_step
-     */
-    void writeStep(int time_step, double time);
-    ~Xdmf3Writer();
-
-private:
-    boost::shared_ptr<XdmfGridCollection> _gridCollection;
-    boost::shared_ptr<XdmfTopology> _initial_topology;
-    boost::shared_ptr<XdmfGeometry> _initial_geometry;
-    std::vector<boost::shared_ptr<XdmfAttribute>> _constant_attributes;
-    std::vector<XdmfData> const _variable_attributes;
-    boost::shared_ptr<XdmfWriter> _writer;
-    boost::shared_ptr<XdmfDomain> _root;
-    std::string const _hdf5filename;
-};
-}  // namespace MeshLib::IO
diff --git a/MeshLib/IO/XDMF/XdmfData.cpp b/MeshLib/IO/XDMF/XdmfData.cpp
index a2e0462e48e4caeb339fa2363a15ba03eff9c443..2aa9bf64cb2275524d6b3ed2240b8e6aaab9d30e 100644
--- a/MeshLib/IO/XDMF/XdmfData.cpp
+++ b/MeshLib/IO/XDMF/XdmfData.cpp
@@ -9,9 +9,6 @@
 
 #include "XdmfData.h"
 
-#include <XdmfArrayType.hpp>
-#include <XdmfAttributeCenter.hpp>
-#include <XdmfTopologyType.hpp>
 #include <cassert>
 #include <map>
 
@@ -20,38 +17,13 @@
 #include "partition.h"
 
 namespace MeshLib::IO
-
-{
-static boost::shared_ptr<XdmfArrayType const> MeshPropertyDataType2XdmfType(
-    MeshPropertyDataType const ogs_data_type)
 {
-    std::map<MeshPropertyDataType, boost::shared_ptr<XdmfArrayType const>> const
-        ogs_to_xdmf_type = {
-            {MeshPropertyDataType::float64, XdmfArrayType::Float64()},
-            {MeshPropertyDataType::float32, XdmfArrayType::Float32()},
-            {MeshPropertyDataType::int32, XdmfArrayType::Int32()},
-            // TODO (tm) XdmfLib does not support 64 bit data types so far
-            //{MeshPropertyDataType::int64, XdmfArrayType::Int64()},
-            {MeshPropertyDataType::uint32, XdmfArrayType::UInt32()},
-            // TODO (tm) XdmfLib does not support 64 bit data types so far
-            //{MeshPropertyDataType::uint64, XdmfArrayType::UInt64},
-            {MeshPropertyDataType::int8, XdmfArrayType::Int8()},
-            {MeshPropertyDataType::uint8, XdmfArrayType::UInt8()}};
-    try
-    {
-        return ogs_to_xdmf_type.at(ogs_data_type);
-    }
-    catch (const std::exception& e)
-    {
-        OGS_FATAL("No known HDF5 type for XDMF type. {:s}", e.what());
-    }
-}
-
 XdmfData::XdmfData(std::size_t const size_partitioned_dim,
                    std::size_t const size_tuple,
                    MeshPropertyDataType const mesh_property_data_type,
                    std::string const& name,
-                   std::optional<MeshLib::MeshItemType> const attribute_center)
+                   std::optional<MeshLib::MeshItemType> const attribute_center,
+                   unsigned int const index)
     : starts([&size_tuple]() {
           if (size_tuple > 1)
           {
@@ -65,15 +37,17 @@ XdmfData::XdmfData(std::size_t const size_partitioned_dim,
       strides([&size_tuple]() {
           if (size_tuple > 1)
           {
-              return std::vector<XdmfDimType>{1};
+              return std::vector<XdmfDimType>{1, 1};
           }
           else
           {
-              return std::vector<XdmfDimType>{1, 1};
+              return std::vector<XdmfDimType>{1};
           }
       }()),
+      data_type(mesh_property_data_type),
       name(name),
-      attribute_center(attribute_center)
+      attribute_center(attribute_center),
+      index(index)
 {
     auto partition_info = getPartitionInfo(size_partitioned_dim);
     // TODO (tm) XdmfLib does not support 64 bit data types so far
@@ -92,7 +66,6 @@ XdmfData::XdmfData(std::size_t const size_partitioned_dim,
         global_block_dims = {ui_global_components, ui_tuple_size};
     }
 
-    data_type = MeshPropertyDataType2XdmfType(mesh_property_data_type);
     DBUG(
         "XDMF: dataset name: {:s}, offset: {:d} "
         "global_blocks: {:d}, tuples: {:d}",
diff --git a/MeshLib/IO/XDMF/XdmfData.h b/MeshLib/IO/XDMF/XdmfData.h
index 959e7d1536b103bc38d273f22fcdbaca39b513ee..b6a826eeed8c5cd00ab4514a8463fca5440817eb 100644
--- a/MeshLib/IO/XDMF/XdmfData.h
+++ b/MeshLib/IO/XDMF/XdmfData.h
@@ -26,7 +26,7 @@ namespace MeshLib::IO
 {
 using XdmfDimType = unsigned int;
 
-struct XdmfData
+struct XdmfData final
 {
     /**
      * \brief XdmfData contains meta data to be passed to the XdmfWriter - it
@@ -43,19 +43,25 @@ struct XdmfData
      * @param attribute_center XdmfData is used for topology, geometry and
      * attributes. Geometry and topology have never a attribute_center.
      * Attributes have always an  attribute_center
+     * @param index The position of the DataItem parents in a grid
+     * (representing a single step). Convention is: 1=Time, 2=
+     * Geometry, 3=Topology, 4>=Attribute
+     *
      */
     XdmfData(std::size_t size_partitioned_dim, std::size_t size_tuple,
              MeshPropertyDataType mesh_property_data_type,
              std::string const& name,
-             std::optional<MeshLib::MeshItemType> attribute_center);
+             std::optional<MeshLib::MeshItemType> attribute_center,
+             unsigned int const index);
     // a hyperslab is defined by starts and strides see
     // https://www.xdmf.org/index.php/XDMF_Model_and_Format#HyperSlab
     std::vector<XdmfDimType> const starts;
     std::vector<XdmfDimType> const strides;
     std::vector<XdmfDimType> global_block_dims;
-    boost::shared_ptr<const XdmfArrayType> data_type;
+    MeshPropertyDataType const data_type;
     std::string const name;
     std::optional<MeshLib::MeshItemType> const attribute_center;
+    unsigned int index;
 };
 
 }  // namespace MeshLib::IO
diff --git a/MeshLib/IO/XDMF/XdmfHdfWriter.cpp b/MeshLib/IO/XDMF/XdmfHdfWriter.cpp
index a33354b118fca1c5ca1b66b321d9ff352319f295..7ec8740380646eb9d6314fb306556014ebf07ac2 100644
--- a/MeshLib/IO/XDMF/XdmfHdfWriter.cpp
+++ b/MeshLib/IO/XDMF/XdmfHdfWriter.cpp
@@ -12,16 +12,17 @@
 #include <algorithm>
 #include <functional>
 
+#include "InfoLib/GitInfo.h"
 #include "partition.h"
 #include "transformData.h"
+#include "writeXdmf.h"
 
 namespace MeshLib::IO
 {
 XdmfHdfWriter::XdmfHdfWriter(MeshLib::Mesh const& mesh,
                              std::filesystem::path const& filepath,
-                             int const time_step,
-                             std::set<std::string>
-                                 variable_output_names,
+                             int const time_step, double const initial_time,
+                             std::set<std::string> const& variable_output_names,
                              bool const use_compression)
 {
     // transform Data into contiguous data and collect meta data
@@ -87,7 +88,13 @@ XdmfHdfWriter::XdmfHdfWriter(MeshLib::Mesh const& mesh,
     std::vector<XdmfData> xdmf_attributes;
     std::transform(attributes.begin(), attributes.end(),
                    std::back_inserter(xdmf_attributes),
-                   [](AttributeMeta att) -> XdmfData { return att.xdmf; });
+                   [&attributes](AttributeMeta& att) -> XdmfData {
+                       size_t const i = &att - &attributes.front();
+                       // index 1 geo, index 2 topology, attributes start at
+                       // index 3
+                       att.xdmf.index = i + 4;
+                       return att.xdmf;
+                   });
 
     std::function<bool(XdmfData)> is_variable_xdmf_attribute =
         [&variable_output_names](XdmfData data) -> bool {
@@ -100,28 +107,33 @@ XdmfHdfWriter::XdmfHdfWriter(MeshLib::Mesh const& mesh,
     std::copy_if(xdmf_attributes.begin(), xdmf_attributes.end(),
                  back_inserter(xdmf_variable_attributes),
                  is_variable_xdmf_attribute);
-
     std::vector<XdmfData> xdmf_constant_attributes;
     std::copy_if(xdmf_attributes.begin(), xdmf_attributes.end(),
                  back_inserter(xdmf_constant_attributes),
                  std::not_fn(is_variable_xdmf_attribute));
 
-    _xdmf_writer = std::make_unique<Xdmf3Writer>(
-        geometry.xdmf, topology.xdmf, std::move(xdmf_constant_attributes),
-        std::move(xdmf_variable_attributes), xdmf_filepath, time_step);
+    if (isFileManager())
+    {
+        auto xdmf_writer_fn =
+            write_xdmf(geometry.xdmf, topology.xdmf, xdmf_constant_attributes,
+                       xdmf_variable_attributes, hdf_filepath.string(),
+                       GitInfoLib::GitInfo::ogs_version);
+        _xdmf_writer = std::make_unique<XdmfWriter>(xdmf_filepath.string(),
+                                                    xdmf_writer_fn);
+        _xdmf_writer->addTimeStep(initial_time);
+    }
 }
 
-void XdmfHdfWriter::writeStep(int const time_step, double const time) const
+void XdmfHdfWriter::writeStep([[maybe_unused]] int const& time_step,
+                              double const& time)
 {
-    _hdf_writer->writeStep(time_step);
-
-    // XDMF
+    // ToDo (tm) time_step will be used for simulation continuation (restart)
+    _hdf_writer->writeStep();
     // The light data is only written by just one process
-    if (!_xdmf_writer)
+    if (isFileManager())
     {
-        return;
+        _xdmf_writer->addTimeStep(time);
     }
-    _xdmf_writer->writeStep(time_step, time);
 }
 
 }  // namespace MeshLib::IO
diff --git a/MeshLib/IO/XDMF/XdmfHdfWriter.h b/MeshLib/IO/XDMF/XdmfHdfWriter.h
index 311821c098ba4a125aa203f80afa2fd350f20137..d4301f0b7ab1c949930ae2daed7e8ab1896ef968 100644
--- a/MeshLib/IO/XDMF/XdmfHdfWriter.h
+++ b/MeshLib/IO/XDMF/XdmfHdfWriter.h
@@ -19,7 +19,7 @@
 
 #include "HdfWriter.h"
 #include "MeshLib/Mesh.h"
-#include "Xdmf3Writer.h"
+#include "XdmfWriter.h"
 
 namespace MeshLib::IO
 {
@@ -31,6 +31,7 @@ public:
      * @param mesh Mesh or NodePartitionedMesh to be written to file(s)
      * @param filepath absolute or relative filepath to the hdf5 file
      * @param time_step number of the step (temporal collection)
+     * @param initial_time time in seconds of the first time step
      * @param variable_output_names names of all process variables (attributes)
      * that change over time
      * @param use_compression if true, zlib compression in HDFWriter component
@@ -38,19 +39,20 @@ public:
      */
     XdmfHdfWriter(MeshLib::Mesh const& mesh,
                   std::filesystem::path const& filepath, int time_step,
-                  std::set<std::string> variable_output_names,
+                  double initial_time,
+                  std::set<std::string> const& variable_output_names,
                   bool use_compression);
+
     /**
-     * \brief Write attribute data that has modified to previous time step or
-     * initial
+     * \brief Adds data for either lazy (xdmf) or eager (hdf) writing algorithm
      * @param time_step number of the step (temporal collection)
      * @param time time value of the current time_step
      */
-    void writeStep(int time_step, double time) const;
+    void writeStep(int const& time_step, double const& time);
 
 private:
     // hdf_writer must be destructed before xdmf_writer
-    std::unique_ptr<Xdmf3Writer> _xdmf_writer;
     std::unique_ptr<HdfWriter> _hdf_writer;
+    std::unique_ptr<XdmfWriter> _xdmf_writer;
 };
 }  // namespace MeshLib::IO
\ No newline at end of file
diff --git a/MeshLib/IO/XDMF/XdmfWriter.cpp b/MeshLib/IO/XDMF/XdmfWriter.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..fa26367eeb2a4de6d3db4d922ce8aa79a0a68233
--- /dev/null
+++ b/MeshLib/IO/XDMF/XdmfWriter.cpp
@@ -0,0 +1,34 @@
+#include "XdmfWriter.h"
+
+#include <fstream>
+
+#include "writeXdmf.h"
+#include "BaseLib/RunTime.h"
+#include "BaseLib/Logging.h"
+
+namespace MeshLib::IO
+{
+XdmfWriter::XdmfWriter(std::string const& xdmf_filename,
+                       std::function<std::string(std::vector<double>)>
+                           xdmf_writer_fn)
+    : filename(xdmf_filename), xdmf_writer(xdmf_writer_fn)
+{
+}
+
+XdmfWriter::~XdmfWriter()
+{
+    BaseLib::RunTime time_output;
+    time_output.start();
+    std::ofstream fout;
+    fout.open(this->filename);
+    fout << xdmf_writer(times);
+
+    INFO("[time] Output of XDMF took {:g} s.", time_output.elapsed());
+}
+
+void XdmfWriter::addTimeStep(double const& time_step)
+{
+    times.push_back(time_step);
+}
+
+}  // namespace MeshLib::IO
\ No newline at end of file
diff --git a/MeshLib/IO/XDMF/XdmfWriter.h b/MeshLib/IO/XDMF/XdmfWriter.h
new file mode 100644
index 0000000000000000000000000000000000000000..c5b675cc19266c4025fbd95e8e880adea9ac152e
--- /dev/null
+++ b/MeshLib/IO/XDMF/XdmfWriter.h
@@ -0,0 +1,48 @@
+/**
+ * \file
+ * \author Tobias Meisel
+ * \date   2021-06-30
+ * \brief  Collects and holds all metadata for writing XDMF file
+ * \copyright
+ * Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ */
+
+#pragma once
+
+#include <functional>
+#include <string>
+
+namespace MeshLib::IO
+{
+class XdmfWriter final
+{
+public:
+    /**
+     * \brief Writes xdmf string into file on class destruction
+     * @param xdmf_filename absolute or relative filepath to the xdmf file
+     * @param xdmf_writer_fn function that generates xdmf string
+     */
+    XdmfWriter(std::string const& xdmf_filename,
+               std::function<std::string(std::vector<double>)>
+                   xdmf_writer_fn);
+    XdmfWriter(XdmfWriter&&) = default;
+    XdmfWriter& operator=(XdmfWriter&&) = default;
+    XdmfWriter(XdmfWriter const&) = delete;
+    XdmfWriter& operator=(XdmfWriter const&) = delete;
+    ~XdmfWriter();
+
+    /**
+     * \brief Adds data for lazy (xdmf) writing algorithm
+     * @param time_step time value of the current time_step
+     */
+    void addTimeStep(double const& time_step);
+
+private:
+    std::string filename;
+    std::vector<double> times;
+    std::function<std::string(std::vector<double>)> xdmf_writer;
+};
+}  // namespace MeshLib::IO
diff --git a/MeshLib/IO/XDMF/transformData.cpp b/MeshLib/IO/XDMF/transformData.cpp
index 10244d1bb9383e6e2c6a7d95817b1e1e45f8ee1a..33ac0e59859425baf2aa297cc3c978f527e703c8 100644
--- a/MeshLib/IO/XDMF/transformData.cpp
+++ b/MeshLib/IO/XDMF/transformData.cpp
@@ -11,11 +11,13 @@
 
 #include "transformData.h"
 
-#include <XdmfTopologyType.hpp>
 #include <algorithm>
+#include <array>
 #include <optional>
 #include <string>
 
+
+#include "BaseLib/cpp23.h"
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Elements/Element.h"
 #include "MeshLib/Mesh.h"
@@ -24,37 +26,47 @@
 #include "MeshPropertyDataType.h"
 #include "partition.h"
 
+
+using namespace BaseLib;
 namespace MeshLib::IO
 {
-// \TODO (tm) constexpr by other function signature can not be transformed to
-// constexpr shared_ptr is not literal type / has non trivial destructor
-boost::shared_ptr<const XdmfTopologyType> cellTypeOGS2XDMF(
-    MeshLib::CellType const cell_type)
+struct XdmfTopology
+{
+    // https://www.xdmf.org/index.php/XDMF_Model_and_Format#Topology, Section
+    // Arbitrary
+    unsigned int id;
+    unsigned int number_of_nodes;
+};
+
+constexpr auto elemOGSTypeToXDMFType()
+{
+    std::array<XdmfTopology, to_underlying(MeshLib::CellType::enum_length)> elem_type{};
+    elem_type[to_underlying(MeshLib::CellType::POINT1)] = {0x1, 1};
+    elem_type[to_underlying(MeshLib::CellType::LINE2)] = {0x2, 2};
+    elem_type[to_underlying(MeshLib::CellType::LINE3)] = {0x2, 3};
+    elem_type[to_underlying(MeshLib::CellType::TRI3)] = {0x4, 3};
+    elem_type[to_underlying(MeshLib::CellType::TRI6)] = {0x24, 6};
+    elem_type[to_underlying(MeshLib::CellType::QUAD4)] = {0x5, 4};
+    elem_type[to_underlying(MeshLib::CellType::QUAD8)] = {0x25, 8};
+    elem_type[to_underlying(MeshLib::CellType::QUAD9)] = {0x23, 9};
+    elem_type[to_underlying(MeshLib::CellType::TET4)] = {0x6, 4};
+    elem_type[to_underlying(MeshLib::CellType::TET10)] = {0x26, 10};
+    elem_type[to_underlying(MeshLib::CellType::HEX8)] = {0x9, 8};
+    elem_type[to_underlying(MeshLib::CellType::HEX20)] = {0x30, 20};
+    elem_type[to_underlying(MeshLib::CellType::HEX27)] = {0x32, 27};
+    elem_type[to_underlying(MeshLib::CellType::PRISM6)] = {
+        0x8, 6};  // parallel triangle wedge
+    elem_type[to_underlying(MeshLib::CellType::PRISM15)] = {0x28, 15};
+    elem_type[to_underlying(MeshLib::CellType::PRISM18)] = {0x29, 18};
+    elem_type[to_underlying(MeshLib::CellType::PYRAMID5)] = {0x7, 5};
+    elem_type[to_underlying(MeshLib::CellType::PYRAMID13)] = {0x27, 13};
+    return elem_type;
+}
+
+auto cellTypeOGS2XDMF(MeshLib::CellType const& cell_type)
 {
-    static std::map<MeshLib::CellType const,
-                    boost::shared_ptr<const XdmfTopologyType> const>
-        elem_type_ogs2xdmf = {
-            {MeshLib::CellType::POINT1, XdmfTopologyType::Polyvertex()},
-            {MeshLib::CellType::LINE2, XdmfTopologyType::Polyline(2)},
-            {MeshLib::CellType::LINE3, XdmfTopologyType::Polyline(3)},
-            {MeshLib::CellType::QUAD4, XdmfTopologyType::Quadrilateral()},
-            {MeshLib::CellType::QUAD8, XdmfTopologyType::Quadrilateral_8()},
-            {MeshLib::CellType::QUAD9, XdmfTopologyType::Quadrilateral_9()},
-            {MeshLib::CellType::TET4, XdmfTopologyType::Tetrahedron()},
-            {MeshLib::CellType::TET10, XdmfTopologyType::Tetrahedron_10()},
-            {MeshLib::CellType::TRI3, XdmfTopologyType::Triangle()},
-            {MeshLib::CellType::TRI6, XdmfTopologyType::Triangle_6()},
-            {MeshLib::CellType::PRISM6,
-             XdmfTopologyType::Wedge()},  // parallel triangle wedge
-            {MeshLib::CellType::PRISM15, XdmfTopologyType::Wedge_15()},
-            {MeshLib::CellType::PRISM18, XdmfTopologyType::Wedge_18()},
-            {MeshLib::CellType::PYRAMID13, XdmfTopologyType::Pyramid_13()},
-            {MeshLib::CellType::PYRAMID5, XdmfTopologyType::Pyramid()},
-            {MeshLib::CellType::HEX20, XdmfTopologyType::Hexahedron_20()},
-            {MeshLib::CellType::HEX27, XdmfTopologyType::Hexahedron_27()},
-            {MeshLib::CellType::HEX8, XdmfTopologyType::Hexahedron()}};
-
-    return elem_type_ogs2xdmf.at(cell_type);
+    constexpr auto elem_type_ogs2xdmf = elemOGSTypeToXDMFType();
+    return elem_type_ogs2xdmf[to_underlying(cell_type)];
 }
 
 std::optional<AttributeMeta> transformAttribute(
@@ -106,33 +118,31 @@ std::optional<AttributeMeta> transformAttribute(
                           "Signed int has 32-1 bits");
             data_type = MeshPropertyDataType::int32;
         }
-        // \TODO (tm) Reimplement size checks
+        // ToDo (tm) These tests are platform specific and currently fail on Windows
         // else if constexpr (std::is_same_v<long, decltype(basic_type)>)
-        // {
-        //     static_assert((std::numeric_limits<long>::digits == 63),
-        //                   "Signed int has 64-1 bits");
-        //     data_type = XdmfArrayType::Int64();
-        // }
+        //{
+        //    static_assert((std::numeric_limits<long>::digits == 63),
+        //                  "Signed int has 64-1 bits");
+        //    data_type = MeshPropertyDataType::int64;
+        //}
+        // else if constexpr (std::is_same_v<unsigned long,
+        // decltype(basic_type)>)
+        //{
+        //    static_assert((std::numeric_limits<unsigned long>::digits == 64),
+        //                  "Unsigned long has 64 bits");
+        //    data_type = MeshPropertyDataType::uint64;
+        //}
         else if constexpr (std::is_same_v<unsigned int, decltype(basic_type)>)
         {
             static_assert((std::numeric_limits<unsigned int>::digits == 32),
                           "Unsigned int has 32 bits");
             data_type = MeshPropertyDataType::uint32;
         }
-        // else if constexpr (std::is_same_v<unsigned long,
-        // decltype(basic_type)>)
-        // {
-        //     static_assert((std::numeric_limits<unsigned long>::digits == 64),
-        //                   "Unsigned long has 64 bits");
-        //     // \TODO (tm) Extend XdmfLibrary with 64bit data types
-        //     data_type = XdmfArrayType::UInt32();
-        // }
         else if constexpr (std::is_same_v<std::size_t, decltype(basic_type)>)
         {
             static_assert((std::numeric_limits<std::size_t>::digits == 64),
                           "size_t has 64 bits");
-            // \TODO (tm) Extend XdmfLibrary with 64bit data types
-            data_type = MeshPropertyDataType::uint32;
+            data_type = MeshPropertyDataType::uint64;
         }
         else if constexpr (std::is_same_v<char, decltype(basic_type)>)
         {
@@ -182,7 +192,7 @@ std::optional<AttributeMeta> transformAttribute(
         HdfData(data_ptr, num_of_tuples, ui_global_components, name, data_type);
 
     XdmfData xdmf = XdmfData(num_of_tuples, ui_global_components, data_type,
-                             name, mesh_item_type);
+                             name, mesh_item_type, 0);
 
     return AttributeMeta{std::move(hdf), std::move(xdmf)};
 }
@@ -234,11 +244,9 @@ Geometry transformGeometry(MeshLib::Mesh const& mesh)
                           point_size,
                           name,
                           MeshPropertyDataType::float64);
-    XdmfData xdmf = XdmfData(partition_dim,
-                             point_size,
-                             MeshPropertyDataType::float64,
-                             name,
-                             std::nullopt);
+    XdmfData xdmf =
+        XdmfData(partition_dim, point_size, MeshPropertyDataType::float64, name,
+                 std::nullopt, 2);
 
     return Geometry{std::move(values), std::move(hdf), std::move(xdmf)};
 }
@@ -247,18 +255,18 @@ Topology transformTopology(MeshLib::Mesh const& mesh, std::size_t const offset)
 {
     std::string const name = "topology";
     std::vector<MeshLib::Element*> const& elements = mesh.getElements();
-    // \TODO (tm) Precalculate exact size
     std::vector<int> values;
     values.reserve(elements.size());
 
     for (auto const& cell : elements)
     {
-        auto const cell_type = cellTypeOGS2XDMF(cell->getCellType());
-        auto const cell_type_id = cell_type->getID();
-        values.push_back(cell_type_id);
-        if (cell_type_id == 2 || cell_type_id == 3)
+        auto const ogs_cell_type = cell->getCellType();
+        auto const xdmf_cell_id = cellTypeOGS2XDMF(ogs_cell_type).id;
+        values.push_back(xdmf_cell_id);
+        if (ogs_cell_type == MeshLib::CellType::LINE2 ||
+            ogs_cell_type == MeshLib::CellType::LINE3)
         {
-            values.push_back(cell_type->getNodesPerElement());
+            values.push_back(cellTypeOGS2XDMF(ogs_cell_type).number_of_nodes);
         }
 
         for (std::size_t i = 0; i < cell->getNumberOfNodes(); ++i)
@@ -271,7 +279,7 @@ Topology transformTopology(MeshLib::Mesh const& mesh, std::size_t const offset)
     HdfData hdf = HdfData(values.data(), values.size(), 1, name,
                           MeshPropertyDataType::int32);
     XdmfData xdmf = XdmfData(values.size(), 1, MeshPropertyDataType::int32,
-                             name, std::nullopt);
+                             name, std::nullopt, 3);
 
     return Topology{std::move(values), std::move(hdf), std::move(xdmf)};
 }
diff --git a/MeshLib/IO/XDMF/writeXdmf.cpp b/MeshLib/IO/XDMF/writeXdmf.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..c1d1bde0f4d5da1fc26a1312922888edad6cab20
--- /dev/null
+++ b/MeshLib/IO/XDMF/writeXdmf.cpp
@@ -0,0 +1,333 @@
+
+#include "writeXdmf.h"
+
+// ToDo (tm) Remove spdlogfmt includes with c++20
+#include <spdlog/fmt/bundled/core.h>
+#include <spdlog/fmt/bundled/format.h>
+
+#include <algorithm>
+#include <fstream>
+#include <iostream>
+#include <iterator>
+#include <list>
+#include <map>
+#include <optional>
+#include <string_view>
+#include <vector>
+
+#include "BaseLib/cpp23.h"
+#include "MeshLib/Location.h"
+
+// naming convention for local function objects by postfix:
+// _transform: functions that take data (mostly XDMF meta data) and return
+// transformed data (XDMF string)
+// _fn: functions that take a function and return a new function
+
+using namespace fmt::literals;
+using namespace BaseLib;
+namespace MeshLib::IO
+{
+// similar definition in Location.h - XDMF uses capital letters
+// Actually there is no correspondece for "IntegrationPoint" in XDMF Format
+// specification
+constexpr char const* mesh_item_type_strings[] = {"Node", "Edge", "Face",
+                                                  "Cell", "IntegrationPoint"};
+
+// Transforms MeshItemType into string
+static std::string mesh_item_type_string(
+    std::optional<MeshItemType> const& item_type)
+{
+    /// Char array names for all of MeshItemType values.
+    if (item_type)
+    {
+        return mesh_item_type_strings[static_cast<int>(item_type.value())];
+    }
+    else
+    {
+        return "";
+    }
+}
+
+// Transforms MeshPropertyDataType into string
+// constexpr with literal type std::string not yet supported
+static auto meshPropertyDatatypeString()
+{
+    std::array<std::string, to_underlying(MeshPropertyDataType::enum_length)>
+        ogs_to_xdmf_type = {};
+    ogs_to_xdmf_type[to_underlying(MeshPropertyDataType::float64)] = "Float";
+    ogs_to_xdmf_type[to_underlying(MeshPropertyDataType::float32)] = "Float";
+    ogs_to_xdmf_type[to_underlying(MeshPropertyDataType::int32)] = "Int";
+    ogs_to_xdmf_type[to_underlying(MeshPropertyDataType::int64)] = "Int";
+    ogs_to_xdmf_type[to_underlying(MeshPropertyDataType::uint32)] = "UInt";
+    ogs_to_xdmf_type[to_underlying(MeshPropertyDataType::uint64)] = "UInt";
+    ogs_to_xdmf_type[to_underlying(MeshPropertyDataType::int8)] = "Int";
+    ogs_to_xdmf_type[to_underlying(MeshPropertyDataType::uint8)] = "UInt";
+    return ogs_to_xdmf_type;
+}
+
+// Transform MeshPropertyDatatype into string )
+static std::string getPropertyDataTypeString(
+    MeshPropertyDataType const& ogs_data_type)
+{
+    auto ogs_to_xdmf_type = meshPropertyDatatypeString();
+    return ogs_to_xdmf_type[to_underlying(ogs_data_type)];
+}
+
+// Returns MeshPropertyDataType-to-Size (in bytes) lookup-table
+static constexpr auto meshPropertyDatatypeSize()
+{
+    std::array<int, to_underlying(MeshPropertyDataType::enum_length)>
+        property_sizes{};
+    property_sizes[to_underlying(MeshPropertyDataType::float64)] = 8;
+    property_sizes[to_underlying(MeshPropertyDataType::float32)] = 4,
+    property_sizes[to_underlying(MeshPropertyDataType::int32)] = 4,
+    property_sizes[to_underlying(MeshPropertyDataType::int64)] = 8,
+    property_sizes[to_underlying(MeshPropertyDataType::uint32)] = 4,
+    property_sizes[to_underlying(MeshPropertyDataType::uint64)] = 8,
+    property_sizes[to_underlying(MeshPropertyDataType::int8)] = 1,
+    property_sizes[to_underlying(MeshPropertyDataType::uint8)] = 1;
+    return property_sizes;
+}
+
+// Transform MeshPropertyDataType into type_sizes (in bytes)
+static std::string getPropertyDataTypeSize(
+    MeshPropertyDataType const& ogs_data_type)
+{
+    constexpr auto ogs_to_xdmf_type = meshPropertyDatatypeSize();
+    return fmt::format("{}", ogs_to_xdmf_type[to_underlying(ogs_data_type)]);
+}
+
+// Collects all known data and return a function that takes all unknown data
+// known: XdmfData of mesh to be written, unknown: timesteps
+std::function<std::string(std::vector<double>)> write_xdmf(
+    XdmfData const& geometry, XdmfData const& topology,
+    std::vector<XdmfData> const& constant_attributes,
+    std::vector<XdmfData> const& variable_attributes,
+    std::string const& h5filename, std::string const& ogs_version)
+{
+    // Generates function that writes <DataItem>. Late binding needed because
+    // maximum number of steps unknown. Time step and h5filename are captured to
+    // return a unary function
+    // _a suffix is a user defined literal for fmt named arguments
+    auto const time_dataitem_genfn = [](int const time_step, int const max_step,
+                                        std::string const h5filename) {
+        return [time_step, max_step, h5filename](auto const& xdmfdata) {
+            return fmt::format(
+                "\n\t\t<DataItem DataType=\"{datatype}\" "
+                "Dimensions=\"{local_dimensions}\" "
+                "Format=\"HDF\" "
+                "Precision=\"{precision}\">{filename}:data/{datasetname}|"
+                "{time_step} {starts}:1 {strides}:1 "
+                "{local_dimensions}:{max_step} "
+                "{global_dimensions}</"
+                "DataItem>",
+                "datatype"_a = getPropertyDataTypeString(xdmfdata.data_type),
+                "local_dimensions"_a =
+                    fmt::join(xdmfdata.global_block_dims, " "),
+                "precision"_a = getPropertyDataTypeSize(xdmfdata.data_type),
+                "filename"_a = h5filename,
+                "datasetname"_a = xdmfdata.name,
+                "starts"_a = fmt::join(xdmfdata.starts, " "),
+                "strides"_a = fmt::join(xdmfdata.strides, " "),
+                "global_dimensions"_a =
+                    fmt::join(xdmfdata.global_block_dims, " "),
+                "time_step"_a = fmt::format("{}", time_step),
+                "max_step"_a = fmt::format("{}", max_step));
+        };
+    };
+
+    // string_join could be moved to ogs lib if used more
+    auto const string_join_fn = [](auto const& collection) {
+        return fmt::to_string(fmt::join(collection, ""));
+    };
+
+    // m_bind could be moved to ogs lib if used more
+    auto const m_bind_fn = [](auto const& transform, auto const& join) {
+        return [join, transform](auto const& collection) {
+            std::vector<std::string> temp;
+            temp.resize(collection.size());
+            std::transform(collection.begin(), collection.end(),
+                           std::back_inserter(temp), transform);
+            return join(temp);
+        };
+    };
+
+    // XDMF if part of the data that is already written in any previous step
+    // subsequent steps can refer to this data collection of elements navigates
+    // the xml tree (take first child, then in child take 2nd child ...)
+    auto const pointer_transfrom = [](auto const& elements) {
+        return fmt::format(
+            "\n\t<xi:include xpointer=\"element(/{elements})\"/>",
+            "elements"_a = fmt::join(elements, "/"));
+    };
+
+    // Defines content of <Attribute> in XDMF, child of Attribute is a single
+    // <DataItem>
+    auto const attribute_transform = [](XdmfData const& attribute,
+                                        auto const& dataitem_transform) {
+        return fmt::format(
+            "\n\t<Attribute Center=\"{center}\" ElementCell=\"\" "
+            "ElementDegree=\"0\" "
+            "ElementFamily=\"\" ItemType=\"\" Name=\"{name}\" "
+            "Type=\"None\">{dataitem}\n\t</Attribute>",
+            "center"_a = mesh_item_type_string(attribute.attribute_center),
+            "name"_a = attribute.name,
+            "dataitem"_a = dataitem_transform(attribute));
+    };
+
+    // Define content of <Geometry> in XDMF, same as attribute_transform
+    auto const geometry_transform = [](XdmfData const& geometry,
+                                       auto const& dataitem_transform) {
+        return fmt::format(
+            "\n\t<Geometry Origin=\"\" Type=\"XYZ\">{dataitem}\n\t</Geometry>",
+            "dataitem"_a = dataitem_transform(geometry));
+    };
+
+    // Define content of <Topology> in XDMF, same as attribute_transform
+    auto const topology_transform = [](XdmfData const& topology,
+                                       auto const& dataitem_transform) {
+        return fmt::format(
+            "\n\t<Topology Dimensions=\"{dimensions}\" "
+            "Type=\"Mixed\">{dataitem}\n\t</Topology>",
+            "dataitem"_a = dataitem_transform(topology),
+            "dimensions"_a = fmt::join(topology.global_block_dims, " "));
+    };
+
+    // Defines content of <Grid> of a single time step, takes specific transform
+    // functions for all of its child elements
+    auto const grid_transform = [](double const& time_value,
+                                   auto const& geometry, auto const& topology,
+                                   auto const& constant_attributes,
+                                   auto const& variable_attributes) {
+        return fmt::format(
+            "\n<Grid Name=\"Grid\">\n\t<Time "
+            "Value=\"{time_value}\"/"
+            ">{geometry}{topology}{fix_attributes}{variable_attributes}\n</"
+            "Grid>",
+            "time_value"_a = time_value, "geometry"_a = geometry,
+            "topology"_a = topology, "fix_attributes"_a = constant_attributes,
+            "variable_attributes"_a = variable_attributes);
+    };
+
+    // An attribute may change over time (variable) or stay constant
+    enum class time_attribute
+    {
+        constant,
+        variable
+    };
+
+    // Generates a function that either writes the data or points to existing
+    // data
+    auto const time_step_fn = [time_dataitem_genfn, pointer_transfrom,
+                               h5filename](auto const& component_transform,
+                                           int const time_step,
+                                           int const max_step,
+                                           time_attribute const variable) {
+        return [component_transform, time_dataitem_genfn, pointer_transfrom,
+                variable, time_step, max_step,
+                h5filename](XdmfData const& attr) {
+            // new data arrived
+            bool changed =
+                ((time_step > 0 && (variable == time_attribute::variable)) ||
+                 time_step == 0);
+            if (changed)
+            {
+                auto dataitem =
+                    time_dataitem_genfn(time_step, max_step, h5filename);
+                return component_transform(attr, dataitem);
+            }
+            else
+            {
+                std::vector<unsigned int> d = {1, 1, 2, 1, attr.index};
+                return pointer_transfrom(d);
+            };
+        };
+    };
+
+    // Top-Level transform function (take spatial and temporal transform
+    // functions) and return the time depended grid transform function
+    // ToDo (tm) Remove capturing m_bind and string_join as helper function
+
+    auto const time_grid_transform =
+        [time_step_fn, m_bind_fn, string_join_fn, grid_transform,
+         geometry_transform, topology_transform, attribute_transform](
+            double const time, int const step, int const max_step,
+            auto const& geometry, auto const& topology,
+            auto const& constant_attributes, auto const& variable_attributes) {
+            auto const time_step_geometry_transform = time_step_fn(
+                geometry_transform, step, max_step, time_attribute::constant);
+            auto const time_step_topology_transform = time_step_fn(
+                topology_transform, step, max_step, time_attribute::constant);
+            auto const time_step_const_attribute_fn = time_step_fn(
+                attribute_transform, step, max_step, time_attribute::constant);
+            auto const time_step_variable_attribute_fn = time_step_fn(
+                attribute_transform, step, max_step, time_attribute::variable);
+
+            // lift both functions from handling single attributes to work with
+            // collection of attributes
+            auto const variable_attributes_transform =
+                m_bind_fn(time_step_variable_attribute_fn, string_join_fn);
+            auto const constant_attributes_transform =
+                m_bind_fn(time_step_const_attribute_fn, string_join_fn);
+
+            return grid_transform(
+                time, time_step_geometry_transform(geometry),
+                time_step_topology_transform(topology),
+                constant_attributes_transform(constant_attributes),
+                variable_attributes_transform(variable_attributes));
+        };
+
+    // Function that combines all the data from spatial (geometry, topology,
+    // attribute) and temporal domain (time) with the time domain
+    // (time_step_fn). And writes <Grid CollectionType="Temporal">
+    auto const temporal_grid_collection_transform =
+        [time_grid_transform](
+            auto const& times, auto const& geometry, auto const& topology,
+            auto const& constant_attributes, auto const& variable_attributes) {
+            // c++20 ranges zip missing
+            auto const max_step = times.size();
+            std::vector<std::string> grids;
+            grids.resize(max_step);
+            for (size_t time_step = 0; time_step < max_step; ++time_step)
+            {
+                grids.push_back(time_grid_transform(
+                    times[time_step], time_step, max_step, geometry, topology,
+                    constant_attributes, variable_attributes));
+            }
+            return fmt::format(
+                "\n<Grid CollectionType=\"Temporal\" GridType=\"Collection\" "
+                "Name=\"Collection\">{grids}\n</Grid>\n",
+                "grids"_a = fmt::join(grids, ""));
+        };
+
+    // Generator function that return a function that represents complete xdmf
+    // structure
+    auto const xdmf_writer_fn = [temporal_grid_collection_transform](
+                                    auto const& ogs_version,
+                                    auto const& geometry, auto const& topology,
+                                    auto const& constant_attributes,
+                                    auto const& variable_attributes) {
+        // This function will be the return of the surrounding named function
+        // capture by copy - it is the function that will survive this scope!!
+        // Use generalized lambda capture to move for optimization
+        return [temporal_grid_collection_transform, ogs_version,
+                geometry = std::move(geometry), topology = std::move(topology),
+                constant_attributes = std::move(constant_attributes),
+                variable_attributes = std::move(variable_attributes)](
+                   std::vector<double> const& times) {
+            return fmt::format(
+                "<?xml version=\"1.0\" encoding=\"utf-8\"?>\n<Xdmf "
+                "xmlns:xi=\"http://www.w3.org/2001/XInclude\" "
+                "Version=\"3.0\">\n<Domain>\n<Information Name=\"OGS_VERSION\" "
+                "Value=\"{ogs_version}\"/>{grid_collection}</Domain>\n</Xdmf>",
+                "ogs_version"_a = ogs_version,
+                "grid_collection"_a = temporal_grid_collection_transform(
+                    times, geometry, topology, constant_attributes,
+                    variable_attributes));
+        };
+    };
+
+    return xdmf_writer_fn(ogs_version, geometry, topology, constant_attributes,
+                          variable_attributes);
+}
+}  // namespace MeshLib::IO
\ No newline at end of file
diff --git a/MeshLib/IO/XDMF/writeXdmf.h b/MeshLib/IO/XDMF/writeXdmf.h
new file mode 100644
index 0000000000000000000000000000000000000000..d713a4fb6605379b1386ae3598f63daab0d76935
--- /dev/null
+++ b/MeshLib/IO/XDMF/writeXdmf.h
@@ -0,0 +1,43 @@
+/**
+ * \file
+ * \author Tobias Meisel
+ * \date   2021-07-13
+ * \brief  write_xdmf generates a function based on spatial mesh data. The
+ * generated function finally generates an XDMF string when temporal data is
+ * known
+ * \copyright
+ * Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ */
+
+#pragma once
+
+#include <vector>
+
+#include "MeshLib/IO/XDMF/XdmfData.h"
+#include "MeshLib/Location.h"
+#include "MeshPropertyDataType.h"
+
+using XdmfDimType = unsigned int;
+
+namespace MeshLib::IO
+{
+/**
+ * \brief Generator function that creates a function capturing the spatial data
+ * of a mesh Temporal data can later be passed as argument
+ * @param geometry Metadata for the geometry (points) of the mesh
+ * @param topology Metadata for the topology of the mesh
+ * @param variable_attributes Meta data for attributes changing over time
+ * @param constant_attributes Meta data for attributes NOT changing over time
+ * @param h5filename Name of the file where the actual data was written
+ * @param ogs_version OGS Version to be added to XdmfInformation tag
+ * @return unary function with vector of time step values, returning XDMF string
+ */
+std::function<std::string(std::vector<double>)> write_xdmf(
+    XdmfData const& geometry, XdmfData const& topology,
+    std::vector<XdmfData> const& constant_attributes,
+    std::vector<XdmfData> const& variable_attributes,
+    std::string const& h5filename, std::string const& ogs_version);
+}  // namespace MeshLib::IO
\ No newline at end of file
diff --git a/MeshLib/IO/writeMeshToFile.cpp b/MeshLib/IO/writeMeshToFile.cpp
index 81d9eb6d83965c6c39f22ecf0566a5270b56ae79..34af711a71166304c0ca847a24225de75d16d2a9 100644
--- a/MeshLib/IO/writeMeshToFile.cpp
+++ b/MeshLib/IO/writeMeshToFile.cpp
@@ -14,9 +14,7 @@
 #include "BaseLib/StringTools.h"
 #include "MeshLib/IO/Legacy/MeshIO.h"
 #include "MeshLib/IO/VtkIO/VtuInterface.h"
-#ifdef OGS_USE_XDMF
 #include "MeshLib/IO/XDMF/XdmfHdfWriter.h"
-#endif
 #include "MeshLib/Mesh.h"
 
 namespace MeshLib::IO
@@ -45,18 +43,15 @@ int writeMeshToFile(const MeshLib::Mesh& mesh,
         }
         return 0;
     }
-#ifdef OGS_USE_XDMF
     if (file_path.extension().string() == ".xdmf")
     {
-        MeshLib::IO::XdmfHdfWriter(
-            mesh, file_path, 0, variable_output_names, true);
-
+        MeshLib::IO::XdmfHdfWriter(mesh, file_path, 0, 0.0,
+                                   variable_output_names, true);
         return 0;
     }
-#endif
-    ERR("writeMeshToFile(): Unknown mesh file format in file {:s}.",
-        file_path.string());
-    return -1;
+    ERR("writeMeshToFile(): Unknown file extension '{:s}'. Can not write file "
+        "'{'s}'.",
+        file_path.extension().string(), file_path.string());
+    return 0;
 }
-
 }  // namespace MeshLib::IO
diff --git a/MeshLib/MeshEnums.h b/MeshLib/MeshEnums.h
index 9a137657584169429ad62ae7f7a7cdb3752a492d..658b9d93faa3c58cacbf4aa1477ab792aaeedf64 100644
--- a/MeshLib/MeshEnums.h
+++ b/MeshLib/MeshEnums.h
@@ -59,7 +59,8 @@ enum class CellType
     PRISM15 = 15,
     PRISM18 = 16,
     PYRAMID5 = 17,
-    PYRAMID13 = 18
+    PYRAMID13 = 18,
+    enum_length
 };
 
 /**
diff --git a/ProcessLib/CMakeLists.txt b/ProcessLib/CMakeLists.txt
index a84b271af04213852f56b59ffa61cd3374422b3f..2c417a1d6d8e1aaa246deb9b984e299074daec8a 100644
--- a/ProcessLib/CMakeLists.txt
+++ b/ProcessLib/CMakeLists.txt
@@ -38,14 +38,12 @@ target_link_libraries(
         $<$<TARGET_EXISTS:ProcessLibSourceTermPython>:ProcessLibSourceTermPython>
         $<$<TARGET_EXISTS:petsc>:petsc>
         nlohmann_json
-    PRIVATE ParameterLib GitInfoLib
-            $<$<TARGET_EXISTS:InSituLib>:InSituLib>
+    PRIVATE ParameterLib GitInfoLib $<$<TARGET_EXISTS:InSituLib>:InSituLib>
 )
 
 target_compile_definitions(
     ProcessLib
-    PUBLIC $<$<BOOL:${OGS_USE_XDMF}>:OGS_USE_XDMF>
-           # Enabled elements
+    PUBLIC # Enabled elements
            OGS_MAX_ELEMENT_DIM=${OGS_MAX_ELEMENT_DIM}
            OGS_MAX_ELEMENT_ORDER=${OGS_MAX_ELEMENT_ORDER}
            $<$<BOOL:${OGS_ENABLE_ELEMENT_SIMPLEX}>:OGS_ENABLE_ELEMENT_SIMPLEX>
diff --git a/ProcessLib/LiquidFlow/Tests.cmake b/ProcessLib/LiquidFlow/Tests.cmake
index ca84d1cc97e2bb3cd22e55c5633954e7a9955f3b..a302aca454319a6fdada998be8a693b33b576a9d 100644
--- a/ProcessLib/LiquidFlow/Tests.cmake
+++ b/ProcessLib/LiquidFlow/Tests.cmake
@@ -482,7 +482,7 @@ AddTest(
     WRAPPER time
     TESTER xdmfdiff
     # See https://gitlab.opengeosys.org/ogs/ogs/-/merge_requests/3184#note_85104
-    REQUIREMENTS NOT OGS_USE_MPI AND OGS_USE_XDMF AND NOT COMPILER_IS_APPLE_CLANG
+    REQUIREMENTS NOT OGS_USE_MPI AND NOT COMPILER_IS_APPLE_CLANG
     DIFF_DATA
     square_5x5_tris_32.xdmf square_5x5_tris_32.xdmf pressure pressure 1e-7 1e-13
     square_5x5_tris_32.xdmf square_5x5_tris_32.xdmf HydraulicFlow HydraulicFlow 1e-7 1e-13
diff --git a/ProcessLib/Output/Output.cpp b/ProcessLib/Output/Output.cpp
index 790bae18b580defcf71ec9c16528fc601e0e3252..73ee24acc75dd77991320a0142b93125f3d25122 100644
--- a/ProcessLib/Output/Output.cpp
+++ b/ProcessLib/Output/Output.cpp
@@ -93,8 +93,8 @@ bool Output::shallDoOutput(int timestep, double const t)
     return false;
 }
 
-Output::Output(std::string output_directory, OutputType output_file_type,
-               std::string output_file_prefix, std::string output_file_suffix,
+Output::Output(std::string directory, OutputType file_type,
+               std::string file_prefix, std::string file_suffix,
                bool const compress_output, std::string const& data_mode,
                bool const output_nonlinear_iteration_results,
                std::vector<PairRepeatEachSteps> repeats_each_steps,
@@ -102,10 +102,10 @@ Output::Output(std::string output_directory, OutputType output_file_type,
                OutputDataSpecification&& output_data_specification,
                std::vector<std::string>&& mesh_names_for_output,
                std::vector<std::unique_ptr<MeshLib::Mesh>> const& meshes)
-    : _output_directory(std::move(output_directory)),
-      _output_file_type(output_file_type),
-      _output_file_prefix(std::move(output_file_prefix)),
-      _output_file_suffix(std::move(output_file_suffix)),
+    : _output_directory(std::move(directory)),
+      _output_file_type(file_type),
+      _output_file_prefix(std::move(file_prefix)),
+      _output_file_suffix(std::move(file_suffix)),
       _output_file_compression(compress_output),
       _output_file_data_mode(convertVtkDataMode(data_mode)),
       _output_nonlinear_iteration_results(output_nonlinear_iteration_results),
@@ -232,7 +232,6 @@ struct Output::OutputFile
     }
 };
 
-#ifdef OGS_USE_XDMF
 void Output::outputMeshXdmf(OutputFile const& output_file,
                             MeshLib::Mesh const& mesh,
                             int const timestep,
@@ -245,12 +244,15 @@ void Output::outputMeshXdmf(OutputFile const& output_file,
         std::filesystem::path path(output_file.path);
         _mesh_xdmf_hdf_writer = std::make_unique<MeshLib::IO::XdmfHdfWriter>(
         MeshLib::IO::XdmfHdfWriter(
-                mesh, path, timestep,
+                mesh, path, timestep, t,
                 _output_data_specification.output_variables, output_file.compression));
     }
+    else
+    {
     _mesh_xdmf_hdf_writer->writeStep(timestep, t);
+    };
 }
-#endif
+
 
 void Output::outputMesh(OutputFile const& output_file,
                         MeshLib::IO::PVDFile* const pvd_file,
@@ -320,7 +322,7 @@ void Output::doOutputAlways(Process const& process,
         }
         else if (_output_file_type == ProcessLib::OutputType::xdmf)
         {
-#ifdef OGS_USE_XDMF
+
             OutputFile const file(
                 _output_directory, _output_file_type, _output_file_prefix, "",
                 mesh.getName(), timestep, t, iteration, _output_file_data_mode,
@@ -328,11 +330,9 @@ void Output::doOutputAlways(Process const& process,
                 _output_data_specification.output_variables);
 
             outputMeshXdmf(file, mesh, timestep, t);
-#else
-            OGS_FATAL(
-                "Trying to write Xdmf file but OGS was not built with "
-                "Xdmf-support.");
-#endif
+
+
+
         }
     };
 
diff --git a/ProcessLib/Output/Output.h b/ProcessLib/Output/Output.h
index 057ce57f0a44e38fab3830f7b7903b9e27f5db93..1cf0709012c3c815e27d423d2771ce6d08444103 100644
--- a/ProcessLib/Output/Output.h
+++ b/ProcessLib/Output/Output.h
@@ -14,9 +14,7 @@
 #include <utility>
 
 #include "MeshLib/IO/VtkIO/PVDFile.h"
-#ifdef OGS_USE_XDMF
 #include "MeshLib/IO/XDMF/XdmfHdfWriter.h"
-#endif
 #include "ProcessOutput.h"
 
 namespace ProcessLib
@@ -40,7 +38,7 @@ public:
     };
 
 public:
-    Output(std::string output_directory, OutputType const type,
+    Output(std::string directory, OutputType const type,
            std::string prefix, std::string suffix, bool const compress_output,
            std::string const& data_mode,
            bool const output_nonlinear_iteration_results,
@@ -93,17 +91,17 @@ private:
                     MeshLib::IO::PVDFile* const pvd_file,
                     MeshLib::Mesh const& mesh,
                     double const t);
-#ifdef OGS_USE_XDMF
+
     void outputMeshXdmf(OutputFile const& output_file,
                         MeshLib::Mesh const& mesh,
                         int const timestep,
                         double const t);
-#endif
+
 
 private:
-#ifdef OGS_USE_XDMF
+
     std::unique_ptr<MeshLib::IO::XdmfHdfWriter> _mesh_xdmf_hdf_writer;
-#endif
+
     std::string const _output_directory;
     OutputType const _output_file_type;
     std::string const _output_file_prefix;
diff --git a/ProcessLib/Output/ProcessOutput.cpp b/ProcessLib/Output/ProcessOutput.cpp
index f223937c5c4a62210d645846375faeaf21c82639..d34ed006228cec6ef273e71afd615a2bd64e9bfa 100644
--- a/ProcessLib/Output/ProcessOutput.cpp
+++ b/ProcessLib/Output/ProcessOutput.cpp
@@ -139,7 +139,7 @@ void addProcessDataToMesh(
     const double t, std::vector<GlobalVector*> const& x, int const process_id,
     MeshLib::Mesh& mesh,
     [[maybe_unused]] std::vector<NumLib::LocalToGlobalIndexMap const*> const&
-        bulk_dof_table,
+        bulk_dof_tables,
     std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
     std::vector<std::reference_wrapper<ProcessVariable>> const&
         process_variables,
@@ -212,7 +212,7 @@ void addProcessDataToMesh(
             for (auto const* node : mesh_subset.getNodes())
             {
 #ifdef USE_PETSC
-                if (bulk_dof_table[process_id] != dof_table[process_id])
+                if (bulk_dof_tables[process_id] != dof_table[process_id])
                 {
                     if (!mesh.getProperties().existsPropertyVector<std::size_t>(
                             "bulk_node_ids"))
@@ -234,7 +234,7 @@ void addProcessDataToMesh(
                     {
                         auto const bulk_node_id =
                             bulk_node_id_map[node->getID()];
-                        // use bulk_dof_table to find information
+                        // use bulk_dof_tables to find information
                         std::size_t const bulk_mesh_id = 0;
                         MeshLib::Location const l(bulk_mesh_id,
                                                   MeshLib::MeshItemType::Node,
@@ -242,7 +242,7 @@ void addProcessDataToMesh(
                         auto const global_component_id =
                             global_component_offset + component_id;
                         auto const index =
-                            bulk_dof_table[process_id]->getLocalIndex(
+                            bulk_dof_tables[process_id]->getLocalIndex(
                                 l, global_component_id,
                                 x[process_id]->getRangeBegin(),
                                 x[process_id]->getRangeEnd());
diff --git a/ProcessLib/Output/ProcessOutput.h b/ProcessLib/Output/ProcessOutput.h
index 07cb8ef6dbf70b9606b3705aebe84987ab2dce20..71541c5877d9e0207017a25138f75abd921bb6b5 100644
--- a/ProcessLib/Output/ProcessOutput.h
+++ b/ProcessLib/Output/ProcessOutput.h
@@ -41,7 +41,7 @@ void addProcessDataToMesh(
     bool const output_secondary_variable,
     std::vector<std::unique_ptr<IntegrationPointWriter>> const*
         integration_point_writer,
-    OutputDataSpecification const& output_data_specification);
+    OutputDataSpecification const& process_output);
 
 //! Writes output to the given \c file_name using the specified file format.
 ///
diff --git a/ProcessLib/SteadyStateDiffusion/Tests.cmake b/ProcessLib/SteadyStateDiffusion/Tests.cmake
index 2d30d35fcc09c9e38f34fe94a4d0eddbb7b8f853..9dd07c7fcee15b08bd0024f83744ee67eeba8672 100644
--- a/ProcessLib/SteadyStateDiffusion/Tests.cmake
+++ b/ProcessLib/SteadyStateDiffusion/Tests.cmake
@@ -423,7 +423,7 @@ AddTest(
     WRAPPER mpirun
     WRAPPER_ARGS -np 3
     TESTER xdmfdiff
-    REQUIREMENTS (OGS_USE_MPI AND OGS_USE_XDMF)
+    REQUIREMENTS OGS_USE_MPI
     DIFF_DATA
     cube_1e3_np3.xdmf cube_1e3_np3.xdmf pressure pressure 1e-3 1e-3
     cube_1e3_np3.xdmf cube_1e3_np3.xdmf v v 1e-3 1e-3
@@ -437,7 +437,7 @@ AddTest(
     WRAPPER mpirun
     WRAPPER_ARGS -np 2
     TESTER xdmfdiff
-    REQUIREMENTS (OGS_USE_MPI AND OGS_USE_XDMF)
+    REQUIREMENTS OGS_USE_MPI
     DIFF_DATA
     cube_1e3_np2.xdmf cube_1e3_np2.xdmf pressure pressure 1e-3 1e-3
     cube_1e3_np2.xdmf cube_1e3_np2.xdmf v v 1e-3 1e-3
diff --git a/scripts/ci/jobs/build-linux-petsc.yml b/scripts/ci/jobs/build-linux-petsc.yml
index 8af9b11b906b28e65d11323042cdfb17634074ca..225c04c09c89106b218b70e9647903b7e839d951 100644
--- a/scripts/ci/jobs/build-linux-petsc.yml
+++ b/scripts/ci/jobs/build-linux-petsc.yml
@@ -9,5 +9,3 @@ build linux petsc:
   variables:
     BUILD_CTEST_LARGE_ON_MASTER: "true"
     CMAKE_PRESET: release-petsc
-    CMAKE_ARGS: >-
-      -DOGS_USE_XDMF=ON
diff --git a/scripts/ci/jobs/build-linux.yml b/scripts/ci/jobs/build-linux.yml
index c3a0a1daab9f02b1eccd292ebc87f23adde2fc36..8a1dcf52956729ed41a84de58fe8be35c6d4a72d 100644
--- a/scripts/ci/jobs/build-linux.yml
+++ b/scripts/ci/jobs/build-linux.yml
@@ -14,7 +14,6 @@ build linux:
       -DOGS_USE_MFRONT=ON
       -DOGS_$USE_PYTHON
       -DOGS_INSTALL_DEPENDENCIES=ON
-      -DOGS_USE_XDMF=ON
   parallel:
     matrix:
       - USE_PYTHON: ["USE_PYTHON=ON", "USE_PYTHON=OFF"]
@@ -36,7 +35,6 @@ linux ctest large:
       -DOGS_USE_CONAN=OFF
       -DOGS_USE_MFRONT=ON
       -DOGS_USE_PYTHON=ON
-      -DOGS_USE_XDMF=ON
 
 build linux (no unity):
   image: $CONTAINER_GCC_IMAGE
@@ -57,7 +55,6 @@ build linux (no unity):
       -DOGS_USE_CONAN=OFF
       -DOGS_USE_MFRONT=ON
       -DOGS_USE_UNITY_BUILDS=OFF
-      -DOGS_USE_XDMF=ON
 
 build linux (no deps, no procs):
   image: $CONTAINER_GCC_IMAGE
diff --git a/scripts/ci/jobs/build-mac.yml b/scripts/ci/jobs/build-mac.yml
index 1c3e79ee2feff28f6385f1c6ee6b041ffc268562..5c46c156f61eb87cf62b412c3fc0a3e877027a7b 100644
--- a/scripts/ci/jobs/build-mac.yml
+++ b/scripts/ci/jobs/build-mac.yml
@@ -11,5 +11,4 @@ build mac:
     CMAKE_PRESET: release
     CMAKE_ARGS: >-
       -DOGS_INSTALL_DEPENDENCIES=ON
-      -DOGS_USE_XDMF=ON
       -DOGS_USE_MFRONT=ON
diff --git a/scripts/ci/jobs/build-win.yml b/scripts/ci/jobs/build-win.yml
index 9c2ed063ede0a1ba66d3403d8767a8421aea293e..ce3ca3d0946315016c64556da598dd5dbedf2354 100644
--- a/scripts/ci/jobs/build-win.yml
+++ b/scripts/ci/jobs/build-win.yml
@@ -1,7 +1,7 @@
 build win:
   extends: .template-build-win
   rules:
-    - if: '$USE_PYTHON =~ /ON$/'
+    - if: "$USE_PYTHON =~ /ON$/"
     - if: '$CI_COMMIT_BRANCH == "master"'
     - if: $CI_COMMIT_BRANCH =~ /^v[0-9]\.[0-9]\.[0-9]/
     - when: manual
@@ -14,7 +14,6 @@ build win:
       -DOGS_CI_TESTRUNNER_REPEAT=1
       -DOGS_$USE_PYTHON
       -DOGS_USE_CONAN=AUTO
-      -DOGS_USE_XDMF=ON
   parallel:
     matrix:
       - USE_PYTHON: ["USE_PYTHON=ON", "USE_PYTHON=OFF"]
diff --git a/scripts/ci/jobs/code-quality.yml b/scripts/ci/jobs/code-quality.yml
index 6967592e56f598d4991548748eb8514cbf33f9bb..58b21d27a99d388859bf63aee74eba5824dfeee8 100644
--- a/scripts/ci/jobs/code-quality.yml
+++ b/scripts/ci/jobs/code-quality.yml
@@ -1,7 +1,7 @@
 cppcheck:
   stage: check
   image: $CONTAINER_GCC_IMAGE
-  needs: [ ci_images, meta ]
+  needs: [ci_images, meta]
   before_script:
     - mkdir -p build
     - cd build
@@ -9,7 +9,7 @@ cppcheck:
     - >
       cmake .. -G Ninja -DCMAKE_BUILD_TYPE=Release -DOGS_USE_CONAN=OFF
       -DOGS_USE_PYTHON=OFF -DOGS_BUILD_UTILS=OFF
-      -DOGS_USE_UNITY_BUILDS=OFF -DOGS_USE_XDMF=ON -DCMAKE_EXPORT_COMPILE_COMMANDS=ON
+      -DOGS_USE_UNITY_BUILDS=OFF -DCMAKE_EXPORT_COMPILE_COMMANDS=ON
     - bash cppcheck.sh
   artifacts:
     reports:
diff --git a/scripts/cmake/Dependencies.cmake b/scripts/cmake/Dependencies.cmake
index 1109abb8529279fa233fccbbffa82fe6e9b0bca1..a1a61fae4d44279b5972c39cf85687c95a848302 100644
--- a/scripts/cmake/Dependencies.cmake
+++ b/scripts/cmake/Dependencies.cmake
@@ -169,85 +169,83 @@ CPMAddPackage(
     OPTIONS "BUILD_SHARED_LIBS OFF"
 )
 
-if(OGS_USE_XDMF)
-    # ZLIB is a HDF5 dependency
-    CPMFindPackage(
-        NAME ZLIB
-        GITHUB_REPOSITORY madler/zlib
-        VERSION 1.2.11
-        EXCLUDE_FROM_ALL YES
-    )
-    if(ZLIB_ADDED)
-        add_library(ZLIB::ZLIB ALIAS zlibstatic)
-    endif()
+# ZLIB is a HDF5 dependency
+CPMFindPackage(
+    NAME ZLIB
+    GITHUB_REPOSITORY madler/zlib
+    VERSION 1.2.11
+    EXCLUDE_FROM_ALL YES
+)
+if(ZLIB_ADDED)
+    add_library(ZLIB::ZLIB ALIAS zlibstatic)
+endif()
 
-    if(OGS_USE_MPI)
-        set(_hdf5_options "HDF5_ENABLE_PARALLEL ON")
-    endif()
+if(OGS_USE_MPI)
+    set(_hdf5_options "HDF5_ENABLE_PARALLEL ON")
+endif()
 
-    string(REPLACE "." "_" HDF5_TAG ${ogs.minimum_version.hdf5})
-    CPMFindPackage(
-        NAME HDF5
-        GITHUB_REPOSITORY HDFGroup/hdf5
-        GIT_TAG hdf5-${HDF5_TAG}
-        VERSION ${ogs.minimum_version.hdf5}
-        OPTIONS "HDF5_EXTERNALLY_CONFIGURED 1"
-                "HDF5_GENERATE_HEADERS OFF"
-                "HDF5_BUILD_TOOLS OFF"
-                "HDF5_BUILD_EXAMPLES OFF"
-                "HDF5_BUILD_HL_LIB OFF"
-                "HDF5_BUILD_FORTRAN OFF"
-                "HDF5_BUILD_CPP_LIB OFF"
-                "HDF5_BUILD_JAVA OFF"
-                ${_hdf5_options}
-        EXCLUDE_FROM_ALL YES
-    )
-    if(HDF5_ADDED)
-        target_include_directories(hdf5-static INTERFACE ${HDF5_BINARY_DIR})
-        list(APPEND DISABLE_WARNINGS_TARGETS hdf5-static)
-        set(HDF5_LIBRARIES hdf5-static)
-        set(HDF5_C_INCLUDE_DIR ${HDF5_SOURCE_DIR})
-        set(HDF5_INCLUDE_DIR ${HDF5_SOURCE_DIR})
-    endif()
+string(REPLACE "." "_" HDF5_TAG ${ogs.minimum_version.hdf5})
+CPMFindPackage(
+    NAME HDF5
+    GITHUB_REPOSITORY HDFGroup/hdf5
+    GIT_TAG hdf5-${HDF5_TAG}
+    VERSION ${ogs.minimum_version.hdf5}
+    OPTIONS "HDF5_EXTERNALLY_CONFIGURED 1"
+            "HDF5_GENERATE_HEADERS OFF"
+            "HDF5_BUILD_TOOLS OFF"
+            "HDF5_BUILD_EXAMPLES OFF"
+            "HDF5_BUILD_HL_LIB OFF"
+            "HDF5_BUILD_FORTRAN OFF"
+            "HDF5_BUILD_CPP_LIB OFF"
+            "HDF5_BUILD_JAVA OFF"
+            ${_hdf5_options}
+    EXCLUDE_FROM_ALL YES
+)
+if(HDF5_ADDED)
+    target_include_directories(hdf5-static INTERFACE ${HDF5_BINARY_DIR})
+    list(APPEND DISABLE_WARNINGS_TARGETS hdf5-static)
+    set(HDF5_LIBRARIES hdf5-static)
+    set(HDF5_C_INCLUDE_DIR ${HDF5_SOURCE_DIR})
+    set(HDF5_INCLUDE_DIR ${HDF5_SOURCE_DIR})
+endif()
 
-    set(XDMF_LIBNAME OgsXdmf CACHE STRING "")
-    CPMAddPackage(
-        NAME xdmf
-        VERSION 3.0.0
-        GIT_REPOSITORY https://gitlab.opengeosys.org/ogs/xdmflib.git
-        GIT_TAG 8d5ae1e1cbf506b8ca2160745fc914e25690c8a4
-        OPTIONS "XDMF_LIBNAME OgsXdmf"
+set(XDMF_LIBNAME OgsXdmf CACHE STRING "")
+CPMAddPackage(
+    NAME xdmf
+    VERSION 3.0.0
+    GIT_REPOSITORY https://gitlab.opengeosys.org/ogs/xdmflib.git
+    GIT_TAG 8d5ae1e1cbf506b8ca2160745fc914e25690c8a4
+    OPTIONS "XDMF_LIBNAME OgsXdmf"
+)
+if(xdmf_ADDED)
+    target_include_directories(
+        OgsXdmf PUBLIC ${xdmf_SOURCE_DIR} ${xdmf_BINARY_DIR}
     )
-    if(xdmf_ADDED)
-        target_include_directories(
-            OgsXdmf PUBLIC ${xdmf_SOURCE_DIR} ${xdmf_BINARY_DIR}
-        )
 
-        target_link_libraries(OgsXdmf Boost::boost ZLIB::ZLIB)
-        target_include_directories(
-            OgsXdmfCore PUBLIC ${xdmf_SOURCE_DIR}/core ${xdmf_BINARY_DIR}/core
-            PRIVATE ${xdmf_SOURCE_DIR}/CMake/VersionSuite
-        )
-        target_link_libraries(
-            OgsXdmfCore PUBLIC Boost::boost LibXml2::LibXml2 ${HDF5_LIBRARIES}
-        )
+    target_link_libraries(OgsXdmf Boost::boost ZLIB::ZLIB)
+    target_include_directories(
+        OgsXdmfCore PUBLIC ${xdmf_SOURCE_DIR}/core ${xdmf_BINARY_DIR}/core
+        PRIVATE ${xdmf_SOURCE_DIR}/CMake/VersionSuite
+    )
+    target_link_libraries(
+        OgsXdmfCore PUBLIC Boost::boost LibXml2::LibXml2 ${HDF5_LIBRARIES}
+    )
 
-        set_target_properties(
-            OgsXdmf OgsXdmfCore
-            PROPERTIES RUNTIME_OUTPUT_DIRECTORY
-                       ${PROJECT_BINARY_DIR}/${CMAKE_INSTALL_BINDIR}
-                       LIBRARY_OUTPUT_DIRECTORY
-                       ${PROJECT_BINARY_DIR}/${CMAKE_INSTALL_LIBDIR}
-                       ARCHIVE_OUTPUT_DIRECTORY
-                       ${PROJECT_BINARY_DIR}/${CMAKE_INSTALL_LIBDIR}
+    set_target_properties(
+        OgsXdmf OgsXdmfCore
+        PROPERTIES RUNTIME_OUTPUT_DIRECTORY
+                    ${PROJECT_BINARY_DIR}/${CMAKE_INSTALL_BINDIR}
+                    LIBRARY_OUTPUT_DIRECTORY
+                    ${PROJECT_BINARY_DIR}/${CMAKE_INSTALL_LIBDIR}
+                    ARCHIVE_OUTPUT_DIRECTORY
+                    ${PROJECT_BINARY_DIR}/${CMAKE_INSTALL_LIBDIR}
+    )
+    if(BUILD_SHARED_LIBS)
+        install(TARGETS OgsXdmf OgsXdmfCore
+                LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR}
         )
-        if(BUILD_SHARED_LIBS)
-            install(TARGETS OgsXdmf OgsXdmfCore
-                    LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR}
-            )
-        endif()
-        list(APPEND DISABLE_WARNINGS_TARGETS OgsXdmf OgsXdmfCore)
     endif()
+    list(APPEND DISABLE_WARNINGS_TARGETS OgsXdmf OgsXdmfCore)
 endif()
 
 if(OGS_BUILD_SWMM)
diff --git a/scripts/cmake/DocumentationSetup.cmake b/scripts/cmake/DocumentationSetup.cmake
index 5397b05341fe5788b5c865ca5bf5b6be2077e90d..a9ce78bab59f666d7c4fc7d00851e20697872a5f 100644
--- a/scripts/cmake/DocumentationSetup.cmake
+++ b/scripts/cmake/DocumentationSetup.cmake
@@ -62,7 +62,6 @@ set(DOXYGEN_PREDEFINED
     OGS_USE_MFRONT
     OGS_USE_NETCDF
     OGS_USE_PYTHON
-    OGS_USE_XDMF
     USE_INSITU
     USE_LIS
     USE_MKL