diff --git a/NumLib/Function/Interpolation.h b/NumLib/Function/Interpolation.h
index 8999d16d0b49b5d9c93e7e3567872fd5acc0fcc2..cc9f5a4fa3ec566b31a3fa38a99e851b9079b738 100644
--- a/NumLib/Function/Interpolation.h
+++ b/NumLib/Function/Interpolation.h
@@ -12,6 +12,9 @@
 #include<array>
 #include<cassert>
 
+#include "NumLib/Fem/CoordinatesMapping/NaturalNodeCoordinates.h"
+#include "NumLib/Fem/FiniteElement/TemplateIsoparametric.h"
+
 namespace NumLib
 {
 
@@ -85,4 +88,55 @@ void shapeFunctionInterpolate(
                                         interpolated_values...);
 }
 
+/// Interpolates scalar \c node_values given in lower order element nodes (e.g.
+/// the base nodes) to higher order element's nodes (e.g. quadratic nodes) and
+/// writes the result into the global property vector.
+///
+/// The base nodes' values are copied.  For each higher order node the shape
+/// matrices are evaluated for the lower order element (the base nodes), and
+/// used for the the scalar quantity interpolation.
+template <typename LowerOrderShapeFunction, typename HigherOrderMeshElementType,
+          int GlobalDim, typename EigenMatrixType>
+void interpolateToHigherOrderNodes(
+    MeshLib::Element const& element, bool const is_axially_symmetric,
+    Eigen::MatrixBase<EigenMatrixType> const& node_values,
+    MeshLib::PropertyVector<double>& interpolated_values_global_vector)
+{
+    assert(dynamic_cast<HigherOrderMeshElementType const*>(&element));
+    assert(node_values.cols() == 1);  // Scalar quantity only.
+
+    using SF = LowerOrderShapeFunction;
+    using ShapeMatricesType = ShapeMatrixPolicyType<SF, GlobalDim>;
+    using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices;
+    using FemType = TemplateIsoparametric<SF, ShapeMatricesType>;
+
+    FemType fe(*static_cast<const typename SF::MeshElement*>(&element));
+    int const number_base_nodes = element.getNumberOfBaseNodes();
+    int const number_all_nodes = element.getNumberOfNodes();
+
+    // Copy the values for linear nodes.
+    for (int n = 0; n < number_base_nodes; ++n)
+    {
+        std::size_t const global_index = element.getNodeIndex(n);
+        interpolated_values_global_vector[global_index] = node_values[n];
+    }
+
+    // Shape matrices storage reused in the interpolation loop.
+    ShapeMatrices shape_matrices{SF::DIM, GlobalDim, SF::NPOINTS};
+
+    // Interpolate values for higher order nodes.
+    for (int n = number_base_nodes; n < number_all_nodes; ++n)
+    {
+        // Evaluated at higher order nodes' coordinates.
+        fe.template computeShapeFunctions<ShapeMatrixType::N>(
+            NaturalCoordinates<HigherOrderMeshElementType>::coordinates[n]
+                .data(),
+            shape_matrices, GlobalDim, is_axially_symmetric);
+
+        std::size_t const global_index = element.getNodeIndex(n);
+        interpolated_values_global_vector[global_index] =
+            shape_matrices.N * node_values;
+    }
+}
+
 } // namespace NumLib
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
index ae05e4465c717fbee315070d85e7cf1e59a3bd2e..723d320358259e38ecd5d7e8ad6a5ef3d09629c3 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
@@ -586,5 +586,22 @@ void HydroMechanicsLocalAssembler<ShapeFunctionDisplacement,
     }
 }
 
+template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
+          typename IntegrationMethod, int DisplacementDim>
+void HydroMechanicsLocalAssembler<ShapeFunctionDisplacement,
+                                  ShapeFunctionPressure, IntegrationMethod,
+                                  DisplacementDim>::
+    computeSecondaryVariableConcrete(double const /*t*/,
+                                     std::vector<double> const& local_x)
+{
+    auto p = Eigen::Map<typename ShapeMatricesTypePressure::template VectorType<
+        pressure_size> const>(local_x.data() + pressure_index, pressure_size);
+
+    NumLib::interpolateToHigherOrderNodes<
+        ShapeFunctionPressure, typename ShapeFunctionDisplacement::MeshElement,
+        DisplacementDim>(_element, _is_axially_symmetric, p,
+                         *_process_data.pressure_interpolated);
+}
+
 }  // namespace HydroMechanics
 }  // namespace ProcessLib
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM.h
index f48950a68cf664d725a0123f859a316870afdfbd..a00bdc1a0b23c3f96e1cc5c55ba0f0c6b0d79275 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsFEM.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM.h
@@ -166,6 +166,8 @@ public:
         }
     }
 
+    void computeSecondaryVariableConcrete(
+        double const t, std::vector<double> const& local_x) override;
     void postNonLinearSolverConcrete(std::vector<double> const& local_x,
                                      double const t,
                                      bool const use_monolithic_scheme) override;
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp b/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp
index 68427b9bcca6015decef96af44cf883caca1ae6c..f6af65800383f77e0a9b0df697ceab0731c7082a 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp
+++ b/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp
@@ -226,6 +226,11 @@ void HydroMechanicsProcess<DisplacementDim>::initializeConcreteProcess(
         makeExtrapolator(mesh.getDimension(), getExtrapolator(),
                          _local_assemblers,
                          &LocalAssemblerInterface::getIntPtDarcyVelocity));
+
+    _process_data.pressure_interpolated =
+        MeshLib::getOrCreateMeshProperty<double>(
+            const_cast<MeshLib::Mesh&>(mesh), "pressure_interpolated",
+            MeshLib::MeshItemType::Node, 1);
 }
 
 template <int DisplacementDim>
@@ -380,6 +385,16 @@ void HydroMechanicsProcess<DisplacementDim>::postNonLinearSolverConcreteProcess(
         getDOFTable(process_id), x, t, _use_monolithic_scheme);
 }
 
+template <int DisplacementDim>
+void HydroMechanicsProcess<DisplacementDim>::computeSecondaryVariableConcrete(
+    const double t, GlobalVector const& x)
+{
+    DBUG("Compute the secondary variables for HydroMechanicsProcess.");
+    GlobalExecutor::executeMemberOnDereferenced(
+        &LocalAssemblerInterface::computeSecondaryVariable, _local_assemblers,
+        *_local_to_global_index_map, t, x, _coupled_solutions);
+}
+
 template <int DisplacementDim>
 std::tuple<NumLib::LocalToGlobalIndexMap*, bool>
 HydroMechanicsProcess<DisplacementDim>::getDOFTableForExtrapolatorData() const
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsProcess.h b/ProcessLib/HydroMechanics/HydroMechanicsProcess.h
index 2d6c630e76fc1cbad39c04203d930bac50d0401e..8cb3b760ee4bd1690332e9a494325d3b972c6524 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsProcess.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsProcess.h
@@ -116,6 +116,8 @@ private:
     /// Solutions of the previous time step
     std::array<std::unique_ptr<GlobalVector>, 2> _xs_previous_timestep;
 
+    void computeSecondaryVariableConcrete(const double t,
+                                          GlobalVector const& x) override;
     /**
      * @copydoc ProcessLib::Process::getDOFTableForExtrapolatorData()
      */
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsProcessData.h b/ProcessLib/HydroMechanics/HydroMechanicsProcessData.h
index 88022050f9d96d4a4c1cf4b95a005f34393f329f..59a61253424354d53ff4c764516d4c9035835a53 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsProcessData.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsProcessData.h
@@ -104,6 +104,8 @@ struct HydroMechanicsProcessData
 
     double const reference_temperature;
 
+    MeshLib::PropertyVector<double>* pressure_interpolated = nullptr;
+
     EIGEN_MAKE_ALIGNED_OPERATOR_NEW;
 };
 
diff --git a/ProcessLib/HydroMechanics/Tests.cmake b/ProcessLib/HydroMechanics/Tests.cmake
index 593c60553563e1438597f279651408542d0a7059..6e51efe9da9e0eed68f60db2ecd8fcb4514de210 100644
--- a/ProcessLib/HydroMechanics/Tests.cmake
+++ b/ProcessLib/HydroMechanics/Tests.cmake
@@ -12,6 +12,7 @@ AddTest(
     DIFF_DATA
     GLOB square_1e2_pcs_0_ts_*.vtu displacement displacement 1e-15 1e-15
     GLOB square_1e2_pcs_0_ts_*.vtu pressure pressure 1e-15 1e-15
+    GLOB square_1e2_pcs_0_ts_*.vtu pressure_interpolated pressure_interpolated 1e-15 1e-15
     GLOB square_1e2_pcs_0_ts_*.vtu velocity velocity 1e-15 1e-15
     GLOB square_1e2_pcs_0_ts_*.vtu HydraulicFlow HydraulicFlow 1e-15 1e-15
     GLOB square_1e2_pcs_0_ts_*.vtu NodalForces NodalForces 1e-15 1e-15
diff --git a/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerMatrix-impl.h b/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerMatrix-impl.h
index 5d6bcd9e6e4398bb2a69f035dbfd6bc79ef46caa..ac79ee973d06f0b35b68c4645f368c3a033112ce 100644
--- a/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerMatrix-impl.h
+++ b/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerMatrix-impl.h
@@ -15,7 +15,7 @@
 #include "MaterialLib/SolidModels/SelectSolidConstitutiveRelation.h"
 #include "MathLib/KelvinVector.h"
 #include "MeshLib/ElementStatus.h"
-#include "NumLib/Fem/CoordinatesMapping/NaturalNodeCoordinates.h"
+#include "NumLib/Function/Interpolation.h"
 #include "ProcessLib/Deformation/LinearBMatrix.h"
 #include "ProcessLib/Utils/InitShapeMatrices.h"
 
@@ -422,46 +422,10 @@ void HydroMechanicsLocalAssemblerMatrix<ShapeFunctionDisplacement,
         (*_process_data.mesh_prop_velocity)[element_id * 3 + i] =
             ele_velocity[i];
 
-    // For each higher order node evaluate the shape matrices for the lower
-    // order element (the base nodes)
-    // TODO (naumov) Extract this method to be useful for other processes.
-    auto interpolate_p = [&]() {
-        using FemType =
-            NumLib::TemplateIsoparametric<ShapeFunctionPressure,
-                                          ShapeMatricesTypePressure>;
-
-        FemType fe(
-            *static_cast<const typename ShapeFunctionPressure::MeshElement*>(
-                &_element));
-        int const number_base_nodes = _element.getNumberOfBaseNodes();
-        int const number_all_nodes = _element.getNumberOfNodes();
-
-        for (int n = 0; n < number_base_nodes; ++n)
-        {
-            std::size_t const global_index = _element.getNodeIndex(n);
-            (*_process_data.mesh_prop_nodal_p)[global_index] = p[n];
-        }
-
-        for (int n = number_base_nodes; n < number_all_nodes; ++n)
-        {
-            // Evaluated at higher order nodes' coordinates.
-            typename ShapeMatricesTypePressure::ShapeMatrices shape_matrices_p{
-                ShapeFunctionPressure::DIM, GlobalDim,
-                ShapeFunctionPressure::NPOINTS};
-
-            fe.computeShapeFunctions(
-                NumLib::NaturalCoordinates<typename ShapeFunctionDisplacement::
-                                               MeshElement>::coordinates[n]
-                    .data(),
-                shape_matrices_p, GlobalDim, _is_axially_symmetric);
-
-            auto const& N_p = shape_matrices_p.N;
-
-            std::size_t const global_index = _element.getNodeIndex(n);
-            (*_process_data.mesh_prop_nodal_p)[global_index] = N_p * p;
-        }
-    };
-    interpolate_p();
+    NumLib::interpolateToHigherOrderNodes<
+        ShapeFunctionPressure, typename ShapeFunctionDisplacement::MeshElement,
+        GlobalDim>(_element, _is_axially_symmetric, p,
+                   *_process_data.mesh_prop_nodal_p);
 }
 
 template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
index c7d2fb5f04e0667c5a2b4ca9b7bea8c60b2c5eae..f386dc0ddd57634c10b151ee9e1e543d91707e1e 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
@@ -832,6 +832,11 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
     saturation_avg /= n_integration_points;
 
     (*_process_data.element_saturation)[_element.getID()] = saturation_avg;
+
+    NumLib::interpolateToHigherOrderNodes<
+        ShapeFunctionPressure, typename ShapeFunctionDisplacement::MeshElement,
+        DisplacementDim>(_element, _is_axially_symmetric, p_L,
+                         *_process_data.pressure_interpolated);
 }
 
 }  // namespace RichardsMechanics
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp
index d1ef1a3cd66cb61a4845b9773352bd416071d9a1..2765c51b04dc5bf4f5cac2ebbb8c5cda7b1186ae 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp
@@ -195,6 +195,11 @@ void RichardsMechanicsProcess<DisplacementDim>::initializeConcreteProcess(
     _process_data.element_saturation = MeshLib::getOrCreateMeshProperty<double>(
         const_cast<MeshLib::Mesh&>(mesh), "saturation_avg",
         MeshLib::MeshItemType::Cell, 1);
+
+    _process_data.pressure_interpolated =
+        MeshLib::getOrCreateMeshProperty<double>(
+            const_cast<MeshLib::Mesh&>(mesh), "pressure_interpolated",
+            MeshLib::MeshItemType::Node, 1);
 }
 
 template <int DisplacementDim>
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsProcessData.h b/ProcessLib/RichardsMechanics/RichardsMechanicsProcessData.h
index f737f7ce28ada2b4bbc29dcfc7e26832ddf99a16..c788c645b81f205237cf31d402f6c21a74248b87 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsProcessData.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsProcessData.h
@@ -104,6 +104,7 @@ struct RichardsMechanicsProcessData
     double t = 0.0;
 
     MeshLib::PropertyVector<double>* element_saturation = nullptr;
+    MeshLib::PropertyVector<double>* pressure_interpolated = nullptr;
 
     EIGEN_MAKE_ALIGNED_OPERATOR_NEW;
 };
diff --git a/ProcessLib/RichardsMechanics/Tests.cmake b/ProcessLib/RichardsMechanics/Tests.cmake
index 0058a30a4475e1979176d20e12c6b097ff49cf20..bae6c283efe059f5757fd98eb5bc1150c26d81a4 100644
--- a/ProcessLib/RichardsMechanics/Tests.cmake
+++ b/ProcessLib/RichardsMechanics/Tests.cmake
@@ -27,8 +27,8 @@ AddTest(
     REQUIREMENTS NOT OGS_USE_MPI
     DIFF_DATA
     GLOB mechanics_linear_pcs_0_ts_*.vtu displacement displacement 1e-15 0
-    GLOB mechanics_linear_pcs_0_ts_*.vtu sigma sigma 1e-15 0
-    GLOB mechanics_linear_pcs_0_ts_*.vtu epsilon epsilon 1e-15 0
+    GLOB mechanics_linear_pcs_0_ts_*.vtu sigma sigma 5e-15 0
+    GLOB mechanics_linear_pcs_0_ts_*.vtu epsilon epsilon 5e-15 0
     GLOB mechanics_linear_pcs_0_ts_*.vtu pressure pressure 1e-15 1e-15
     GLOB mechanics_linear_pcs_0_ts_*.vtu velocity velocity 1e-15 1e-15
     GLOB mechanics_linear_pcs_0_ts_*.vtu HydraulicFlow HydraulicFlow 1e-15 0
@@ -84,6 +84,7 @@ AddTest(
     GLOB RichardsFlow_2d_small_pcs_0_ts_*.vtu sigma sigma 1e-8 0
     GLOB RichardsFlow_2d_small_pcs_0_ts_*.vtu epsilon epsilon 1e-15 0
     GLOB RichardsFlow_2d_small_pcs_0_ts_*.vtu pressure pressure 1e-7 1e-15
+    GLOB RichardsFlow_2d_small_pcs_0_ts_*.vtu pressure_interpolated pressure_interpolated 1e-7 1e-15
     GLOB RichardsFlow_2d_small_pcs_0_ts_*.vtu saturation saturation 1e-11 1e-15
     GLOB RichardsFlow_2d_small_pcs_0_ts_*.vtu saturation_avg saturation_avg 1e-11 1e-15
     GLOB RichardsFlow_2d_small_pcs_0_ts_*.vtu velocity velocity 1e-15 1e-15
@@ -118,6 +119,6 @@ AddTest(
     TESTER vtkdiff
     REQUIREMENTS NOT OGS_USE_MPI
     DIFF_DATA
-    GLOB RichardsFlow_2d_richardsflow_pcs_0_ts_*.vtu pressure pressure 1e-11 1e-15
+    GLOB RichardsFlow_2d_richardsflow_pcs_0_ts_*.vtu pressure pressure 5e-11 1e-15
     GLOB RichardsFlow_2d_richardsflow_pcs_0_ts_*.vtu saturation saturation 1e-14 1e-15
 )
diff --git a/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_120_t_1000.000000.vtu b/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_120_t_1000.000000.vtu
deleted file mode 100644
index 4a647f4e86f3d85fc0923bfe4f6e16f37f9b5893..0000000000000000000000000000000000000000
--- a/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_120_t_1000.000000.vtu
+++ /dev/null
@@ -1,3 +0,0 @@
-version https://git-lfs.github.com/spec/v1
-oid sha256:b576bd022a64ae72d8112bc37f00c3ed496fd7b4b9ddb04672e2f6fc0a2bd021
-size 33801
diff --git a/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_1_t_5.000000.vtu b/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_1_t_5.000000.vtu
deleted file mode 100644
index 4a559babba3bff5ceecaa3b9634a4e29a7a8b30c..0000000000000000000000000000000000000000
--- a/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_1_t_5.000000.vtu
+++ /dev/null
@@ -1,3 +0,0 @@
-version https://git-lfs.github.com/spec/v1
-oid sha256:5f68d60f488aa6b2e370e0ae199649faaf01144b9d261d9be90185b867a6cc1a
-size 35031
diff --git a/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_20_t_100.000000.vtu b/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_20_t_100.000000.vtu
deleted file mode 100644
index 321557402f20fc574769ea46cc0282b63e84821f..0000000000000000000000000000000000000000
--- a/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_20_t_100.000000.vtu
+++ /dev/null
@@ -1,3 +0,0 @@
-version https://git-lfs.github.com/spec/v1
-oid sha256:a78d1ef315921e45241b37bd554ae60525a08aa049fa092ddc54092feed12460
-size 32115
diff --git a/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_420_t_4000.000000.vtu b/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_420_t_4000.000000.vtu
deleted file mode 100644
index 35b001e2157ba707861a26e5c0a78a5052a58958..0000000000000000000000000000000000000000
--- a/Tests/Data/HydroMechanics/Linear/Confined_Compression/expected_square_1e2_pcs_0_ts_420_t_4000.000000.vtu
+++ /dev/null
@@ -1,3 +0,0 @@
-version https://git-lfs.github.com/spec/v1
-oid sha256:221b681e63a7c7e098dd0f139365d0c4aa5809c25fa03e105168b557d569a3a8
-size 36871
diff --git a/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_120_t_1000.000000.vtu b/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_120_t_1000.000000.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..537e94192547ea6a17fb09f407d62fc7f7fa00b4
--- /dev/null
+++ b/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_120_t_1000.000000.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:0b3a1cc0ba1b4e1403c73abd1bbf4e59a5d02b8f910312c7ff1169497708763c
+size 36381
diff --git a/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_1_t_5.000000.vtu b/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_1_t_5.000000.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..53198af58f6f6a229c259dbf9cd6f02d84516ade
--- /dev/null
+++ b/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_1_t_5.000000.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:e68238d8160c4b5bc280d49c3856fa139bc720c94538af6b5c74908ae57f7d4d
+size 35972
diff --git a/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_20_t_100.000000.vtu b/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_20_t_100.000000.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..31c46b8f0f13009bfaa6e77c0de7bc947a596c6a
--- /dev/null
+++ b/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_20_t_100.000000.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:6b2a55261b638241e2cedc923779fd8be9f09ce676e3713068139e21d5c0d696
+size 33155
diff --git a/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_420_t_4000.000000.vtu b/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_420_t_4000.000000.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..6fceceb02bc9c300c1aca808af9b5736373e2cda
--- /dev/null
+++ b/Tests/Data/HydroMechanics/Linear/Confined_Compression/square_1e2_pcs_0_ts_420_t_4000.000000.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:735cdbe865eeb4320bbd913fd0abb2c6f6c00b39cc5d786ddaffadadcde62423
+size 41019
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_0_t_0.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_0_t_0.000000.vtu
index 080b872b8f572a200ee62cecfa1f01ff4170437e..828d37e33717a6312654d9baa7beee013e41b36b 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_0_t_0.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_0_t_0.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:9b777be31f90ab6ab765aa0f75d0ffb288bf905e6a22cd6dab4e9601881e0d2d
-size 7907
+oid sha256:b1cd831a64de78b22ff3b11760937780898c46633d07fb3ac612dfaa3abcb1f7
+size 12641
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_104_t_2000.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_104_t_2000.000000.vtu
index 1bc0d62c1caa49b02f3e105fca4bce6f65c2b5f9..642e5b345515b30bc1f555e6d43331dc793e9554 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_104_t_2000.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_104_t_2000.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:6712e56a933409dab30238242d7ea73f7303012c31fb2f6cef5880dd03da51ab
-size 40439
+oid sha256:68fc9c5aba72e898fcc54bdeee27aa870337fa6a10bd9e0e2468ebbddb98df82
+size 47016
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_19_t_318.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_19_t_318.000000.vtu
index a3ea5cbceb8383018b78c016bbfc59ca67a17a45..ae7f2e71e1c0c28c3de4a221c7d0390ab7259104 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_19_t_318.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_19_t_318.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:f8ad15919d42474e466bb56bf298a0294044c6f8f9168554699f5a32abebdff7
-size 40771
+oid sha256:f3f31960e8194950a522354270d3460aff9fd5ebc7b88af9a7dd75c4eb105c59
+size 47068
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_1_t_1.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_1_t_1.000000.vtu
index 236b96d902fcd3f4707222519d57a005f0788a20..b136564e35b21d6c6fdb78a095c73dd4f33a2310 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_1_t_1.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_1_t_1.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:d945de267a2957763e3593bf0207ee5c7086d58d181d6df23b3e8200bad4115b
-size 40139
+oid sha256:bb58b40038941564364dc829b491d7bddb0051e1e06a4bad2dd0115495eea5dc
+size 47003
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_29_t_518.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_29_t_518.000000.vtu
index f14803a9dc0b28da7a07365ce4998565bbde8d66..d63a5641e2911bc86b141a5da59a72064398988c 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_29_t_518.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_29_t_518.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:ada9d4d5b1d7f7d84b4a9d244e28cfb0a1623edaf59d385e088a2396c5fb8f2b
-size 40282
+oid sha256:8af25d729b0de8daf90de81a9a0ac3e3ca96a87b6632de45bec8050feb50d056
+size 46782
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_2_t_3.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_2_t_3.000000.vtu
index b7472fd857e4fb2e34b0679d896305326d728173..b7e9290e55d40ad2a89f718f85ebb1d1897830cd 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_2_t_3.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_2_t_3.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:ee70c0dcf7cf252965f40cee6d3604d8ee7fc47bf4ce072cef66f663ef68cea9
-size 39308
+oid sha256:e28773f4449ff14cc9e2b2b0e9c242d85886d3a77a8bdb781f615a902336e20f
+size 46505
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_39_t_718.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_39_t_718.000000.vtu
index f9bda67dbdbbcf74788a8e6c7fddce134d8b1e27..fa6e9b704ddb9504c5917dc476d5af05e45be589 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_39_t_718.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_39_t_718.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:1b7b592b267fd583ac528fb7321edba8f26e1fd14e58942d202d753911cac9f9
-size 40467
+oid sha256:c5d9a27edb806493fbd1b977066c11fb7b81e900b121703e47fb09dbaa6ad239
+size 47036
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_3_t_8.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_3_t_8.000000.vtu
index 1e621e899aa90268e6bd1a4eb7dc4129fe853a36..d387dce2fd5f1818cdc7e511f657e710d5b20222 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_3_t_8.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_3_t_8.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:8388d9daefc0e0edd69fd992201bb32c7afe70f5b815e56e59559b9d5d297751
-size 39796
+oid sha256:98127ecf088b0badd5828b1aa8d988f0e53cd8825ede7d7059fee26ac562fd87
+size 46599
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_49_t_918.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_49_t_918.000000.vtu
index fec7eba0ccd1ae0466ba8c48e0d0d2fbd1ee7b29..46dc6e080009ee77e595e48d5287c7adadc844c8 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_49_t_918.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_49_t_918.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:eb72701046aa78f9719d1c2209d7b2cc39eac2409c37ed50445455b3bb8cf8ab
-size 40463
+oid sha256:3fd5bca114609468ffa46224a9aae42552d8ef1506c811c5beb1d768dc8db8e5
+size 46825
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_4_t_18.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_4_t_18.000000.vtu
index a587648e6d615bb77d82fdbbe707fb88b572a6c6..b2df7d746127562becb494f8f984a9284e73df8f 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_4_t_18.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_4_t_18.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:0bb02dc84476c959f044fbfe136a637c6ebdbf10bff61a48023e6ee040b47177
-size 39528
+oid sha256:a1ce7ea656a7e8a45e9e81f0bc19f564bef9a8a9c3ab50d4d3506e4cd786e1f9
+size 46532
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_59_t_1118.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_59_t_1118.000000.vtu
index 2296ca6296868ba63b6fcb9823688d26bdace41b..8f3bec67fc717975d25d2fb1a643497fc2a789b1 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_59_t_1118.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_59_t_1118.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:21ff9bf6479f701c03474482f047d4ff189298bf769f2a8e6fe3409d4278e4b0
-size 41155
+oid sha256:4a9571736fc38bc762e3478cafd38331f20568091666e454b84184ee46bc4c55
+size 47472
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_5_t_38.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_5_t_38.000000.vtu
index 5aed7fa73c3933bad2b5c07810dfcd5cb6e63fee..1ac9fae2bcefdc61843dc13976b0d3497ed579af 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_5_t_38.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_5_t_38.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:3e17d83929c78459462a6282e9c1fa04accf1faca232dc17e4b5631dd4006986
-size 39900
+oid sha256:b65ef652209f71e0039168bb03fa5e3c470c7123641128fe43e8f6c56044db7e
+size 46593
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_69_t_1318.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_69_t_1318.000000.vtu
index b6f5d111bb9d8a3fc7e8bc34a1a1c69bf60acaca..83a26e2b7ddc31c950697195749a171ecb7335ee 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_69_t_1318.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_69_t_1318.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:39f65cc2000506d12da621a6cf8ea4c3b6fd79f88c6ef9fd37ee208b4e8b53d6
-size 40575
+oid sha256:977f973eb5b267915923e1513f0c8246b9217cf695684a0f1f31f91ed35ae476
+size 46744
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_6_t_58.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_6_t_58.000000.vtu
index b0a4eb5a7e3e99de34d59c59d887e84cf8bdb357..15c2ba6a02e4abefa5121684621dd10cede615c5 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_6_t_58.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_6_t_58.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:09abff3b241912ca0adc9b1a8b027f8c4fc019496f4b35e7caa2e4cec431e3f8
-size 40335
+oid sha256:abc582474cf8949a86d4168055d8e6da225156464f04abdcdcfcea7df571603b
+size 46753
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_79_t_1518.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_79_t_1518.000000.vtu
index a392faab134dbf9e179647244bfb317f12337e55..15bf6387497757f9381174db077fcf7c7ab4c6c2 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_79_t_1518.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_79_t_1518.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:323b6efa89483400b93c48ab0b723cb5399cebd26be648b3c471c83584d1e04e
-size 40806
+oid sha256:c64b6999f512ebff36226be3ce3c91d15723407a4448fcba6d5d3c1883736dd0
+size 47274
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_7_t_78.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_7_t_78.000000.vtu
index 1a8afbedc08ffd1255f56d9df06b2591b05dc5be..c7542e4b34f329e7396b82797174e834b51af033 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_7_t_78.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_7_t_78.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:7602b7dbb8ea84b51b2e86ce3688931bad0a02f7e2862117451dc5593344ddb5
-size 40622
+oid sha256:057baf03f7376b8fbb4e08526ab22408c339fd1bc591308f1d4dcda761de2831
+size 47081
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_89_t_1718.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_89_t_1718.000000.vtu
index 87b9a181831248acea03be1111eed1e1e7b06787..94fa73b48c2add0b737af2b883ed8199c0b32666 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_89_t_1718.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_89_t_1718.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:d31d9d6f8481319bad0b4ee829b73c721f6f3eb3911ca969310ca2d6ad03e655
-size 40620
+oid sha256:17f79fba986c17c35179aaeaafc3215187a2e5be8740f96c8e2c540a2fc1b634
+size 47294
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_8_t_98.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_8_t_98.000000.vtu
index 0b91f3fbf420674825b4a378063b440106931ddb..c84aa000d6640115c69aeb561339bb182e1c9ce0 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_8_t_98.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_8_t_98.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:47249c58b1811a5aa7c2afc62094ea980ac52732f0b694452f78955fb4af3db0
-size 40616
+oid sha256:e60454423f661f44f5ab9ea5da10236eeb1f26a43ec657be34b27aca19fe70c8
+size 46967
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_99_t_1918.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_99_t_1918.000000.vtu
index 51f5dc89f5bf87dd25139afee07c9625e1782c2a..3698086098af26b64e160921f886168f75f9886c 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_99_t_1918.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_99_t_1918.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:d94176b5141381071f3d498f270c6b42990ce7e8060783f8bbbac3e54c1ac804
-size 41069
+oid sha256:17ecd0b70493765269f729c69207559e960eb70ed7ce69cbeefc38084fa017c5
+size 47399
diff --git a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_9_t_118.000000.vtu b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_9_t_118.000000.vtu
index 73486a6b3c086ed7a8fa8f63958eba035c37ec7f..36b958d74998e8ffb79ca19276bc52080d8b103e 100644
--- a/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_9_t_118.000000.vtu
+++ b/Tests/Data/RichardsMechanics/RichardsFlow_2d_small_pcs_0_ts_9_t_118.000000.vtu
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:d3d53e70c0e6c59b71124ffe1e473577afe734ae3d21547615a31f3e0eb7e2f4
-size 39900
+oid sha256:736e7e4353509ca2340313359ffb4978c21e0e4327b6adf703ea0487180a539f
+size 46515