diff --git a/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/c_SaturationVanGenuchten.md b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/c_SaturationVanGenuchten.md
new file mode 100644
index 0000000000000000000000000000000000000000..52a5da84e631b26049d2acc01c416417653728ba
--- /dev/null
+++ b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/c_SaturationVanGenuchten.md
@@ -0,0 +1 @@
+\copydoc MaterialPropertyLib::SaturationVanGenuchtenWithVolumetricStrain
diff --git a/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_a.md b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_a.md
new file mode 100644
index 0000000000000000000000000000000000000000..77ef249420d3bf133719c97000b11667ce43af1f
--- /dev/null
+++ b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_a.md
@@ -0,0 +1 @@
+Scaling parameter for the effect of the volumetric strain.
diff --git a/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_d_diff.md b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_d_diff.md
new file mode 100644
index 0000000000000000000000000000000000000000..a3a411c4dd46cdb0d8c5e2b7e1eb1666c5b8465d
--- /dev/null
+++ b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_d_diff.md
@@ -0,0 +1 @@
+Pore size density parameter to calculate \f$p_b\f$ values between micro and macropores.
diff --git a/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_e_0.md b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_e_0.md
new file mode 100644
index 0000000000000000000000000000000000000000..c1a181af997609c9806a0ce470f0d757a0e05d20
--- /dev/null
+++ b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_e_0.md
@@ -0,0 +1,4 @@
+Total void ratio.
+
+It can be calculated using the porosity (\f$\phi\f$) as follows:
+\f$e_0 =\phi/(1-\phi)\f$.
diff --git a/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_e_m.md b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_e_m.md
new file mode 100644
index 0000000000000000000000000000000000000000..1a72232ee3571aab284f68408d8501022ba01b32
--- /dev/null
+++ b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_e_m.md
@@ -0,0 +1,6 @@
+Void ratio of the micropores.
+
+It can be calculated as
+\f$e_{\text{m}} = \frac{V_{\text{m}}}{V_{\text{s}}_{\text{m}}}\f$, where
+\f$V_{\text{m}}\f$ is the volume of the micropores, and \f$V_{\text{s}}\f$ is
+the solid volume.
diff --git a/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_exponent.md b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_exponent.md
new file mode 100644
index 0000000000000000000000000000000000000000..260b1dd015e7c9cb8b04a65e4df7587dc81b9224
--- /dev/null
+++ b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_exponent.md
@@ -0,0 +1 @@
+The exponent of the van Genuchten saturation function.
diff --git a/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_p_b.md b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_p_b.md
new file mode 100644
index 0000000000000000000000000000000000000000..146d2d144354f812ae014c59557ed90cd7496242
--- /dev/null
+++ b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_p_b.md
@@ -0,0 +1,9 @@
+Pressure scaling parameter with units of pressure for the micropores.
+
+Corresponds to the van Genuchten parameter \f$\rho g / alpha\f$.
+The inverse of \f$\alpha\f$ "corresponds to the water contents lower then
+saturation" [1].
+
+[1] van Lier, Quirijn de Jong, and Everton Alves Rodrigues Pinheiro. "An alert
+regarding a common misinterpretation of the van Genuchten α parameter." Revista
+Brasileira de Ciência do Solo 42 (2018).
diff --git a/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_residual_gas_saturation.md b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_residual_gas_saturation.md
new file mode 100644
index 0000000000000000000000000000000000000000..89df995e0324903f98b9b6308ac7164848603adc
--- /dev/null
+++ b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_residual_gas_saturation.md
@@ -0,0 +1,2 @@
+The smallest degree of saturation of the gas (or non-wetting) phase in the
+medium.
diff --git a/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_residual_liquid_saturation.md b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_residual_liquid_saturation.md
new file mode 100644
index 0000000000000000000000000000000000000000..fd93070baedc37e28fa006e76180cddc21afa8e5
--- /dev/null
+++ b/Documentation/ProjectFile/properties/property/SaturationVanGenuchtenWithVolumetricStrain/t_residual_liquid_saturation.md
@@ -0,0 +1,2 @@
+The smallest degree of saturation of the liquid (or wetting) phase in the
+medium.
diff --git a/MaterialLib/MPL/CreateProperty.cpp b/MaterialLib/MPL/CreateProperty.cpp
index 2c4f5035361ee78d4a774e585f984391f54aea40..20279a58b356353bb570c013ad3df40dcbba799d 100644
--- a/MaterialLib/MPL/CreateProperty.cpp
+++ b/MaterialLib/MPL/CreateProperty.cpp
@@ -195,6 +195,12 @@ std::unique_ptr<MaterialPropertyLib::Property> createProperty(
return createSaturationVanGenuchten(config);
}
+ if (boost::iequals(property_type,
+ "SaturationVanGenuchtenWithVolumetricStrain"))
+ {
+ return createSaturationVanGenuchtenWithVolumetricStrain(config);
+ }
+
if (boost::iequals(property_type, "CapillaryPressureVanGenuchten"))
{
return createCapillaryPressureVanGenuchten(config);
diff --git a/MaterialLib/MPL/Properties/CapillaryPressureSaturation/CreateSaturationVanGenuchtenWithVolumetricStrain.cpp b/MaterialLib/MPL/Properties/CapillaryPressureSaturation/CreateSaturationVanGenuchtenWithVolumetricStrain.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..4f67eff06d9ad27bbb92e8370e9f6b6a2aa3a401
--- /dev/null
+++ b/MaterialLib/MPL/Properties/CapillaryPressureSaturation/CreateSaturationVanGenuchtenWithVolumetricStrain.cpp
@@ -0,0 +1,60 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2023, 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 "BaseLib/ConfigTree.h"
+#include "SaturationVanGenuchtenWithVolumetricStrain.h"
+
+namespace MaterialPropertyLib
+{
+std::unique_ptr<SaturationVanGenuchtenWithVolumetricStrain>
+createSaturationVanGenuchtenWithVolumetricStrain(
+ BaseLib::ConfigTree const& config)
+{
+ //! \ogs_file_param{properties__property__type}
+ config.checkConfigParameter("type",
+ "SaturationVanGenuchtenWithVolumetricStrain");
+
+ // Second access for storage.
+ //! \ogs_file_param{properties__property__name}
+ auto property_name = config.peekConfigParameter<std::string>("name");
+
+ DBUG(
+ "Create SaturationVanGenuchtenWithVolumetricStrain medium property "
+ "{:s}.",
+ property_name);
+
+ auto const residual_liquid_saturation =
+ //! \ogs_file_param{properties__property__SaturationVanGenuchtenWithVolumetricStrain__residual_liquid_saturation}
+ config.getConfigParameter<double>("residual_liquid_saturation");
+ auto const residual_gas_saturation =
+ //! \ogs_file_param{properties__property__SaturationVanGenuchtenWithVolumetricStrain__residual_gas_saturation}
+ config.getConfigParameter<double>("residual_gas_saturation");
+ auto const exponent =
+ //! \ogs_file_param{properties__property__SaturationVanGenuchtenWithVolumetricStrain__exponent}
+ config.getConfigParameter<double>("exponent");
+ //! \ogs_file_param{properties__property__SaturationVanGenuchtenWithVolumetricStrain__p_b}
+ auto const p_b = config.getConfigParameter<double>("p_b");
+ auto const e_0 =
+ //! \ogs_file_param{properties__property__SaturationVanGenuchtenWithVolumetricStrain__e_0}
+ config.getConfigParameter<double>("e_0");
+ auto const e_m =
+ //! \ogs_file_param{properties__property__SaturationVanGenuchtenWithVolumetricStrain__e_m}
+ config.getConfigParameter<double>("e_m");
+ auto const a =
+ //! \ogs_file_param{properties__property__SaturationVanGenuchtenWithVolumetricStrain__a}
+ config.getConfigParameter<double>("a");
+ auto const d_diff =
+ //! \ogs_file_param{properties__property__SaturationVanGenuchtenWithVolumetricStrain__d_diff}
+ config.getConfigParameter<double>("d_diff");
+
+ return std::make_unique<SaturationVanGenuchtenWithVolumetricStrain>(
+ std::move(property_name), residual_liquid_saturation,
+ residual_gas_saturation, exponent, p_b, e_0, e_m, a, d_diff);
+}
+} // namespace MaterialPropertyLib
diff --git a/MaterialLib/MPL/Properties/CapillaryPressureSaturation/CreateSaturationVanGenuchtenWithVolumetricStrain.h b/MaterialLib/MPL/Properties/CapillaryPressureSaturation/CreateSaturationVanGenuchtenWithVolumetricStrain.h
new file mode 100644
index 0000000000000000000000000000000000000000..78c9658f863c532eba76afaca002a70385c60dae
--- /dev/null
+++ b/MaterialLib/MPL/Properties/CapillaryPressureSaturation/CreateSaturationVanGenuchtenWithVolumetricStrain.h
@@ -0,0 +1,29 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2023, 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 <memory>
+
+namespace BaseLib
+{
+class ConfigTree;
+}
+
+namespace MaterialPropertyLib
+{
+class SaturationVanGenuchtenWithVolumetricStrain;
+}
+
+namespace MaterialPropertyLib
+{
+std::unique_ptr<SaturationVanGenuchtenWithVolumetricStrain>
+createSaturationVanGenuchtenWithVolumetricStrain(
+ BaseLib::ConfigTree const& config);
+} // namespace MaterialPropertyLib
diff --git a/MaterialLib/MPL/Properties/CapillaryPressureSaturation/SaturationVanGenuchtenWithVolumetricStrain.cpp b/MaterialLib/MPL/Properties/CapillaryPressureSaturation/SaturationVanGenuchtenWithVolumetricStrain.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..b4c08d67dc5c66e91150ff9318cdb52214677ed2
--- /dev/null
+++ b/MaterialLib/MPL/Properties/CapillaryPressureSaturation/SaturationVanGenuchtenWithVolumetricStrain.cpp
@@ -0,0 +1,129 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2023, 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 "SaturationVanGenuchtenWithVolumetricStrain.h"
+
+#include <algorithm>
+#include <cmath>
+#include <limits>
+
+#include "BaseLib/Error.h"
+#include "MaterialLib/MPL/Medium.h"
+#include "MathLib/MathTools.h"
+
+namespace MaterialPropertyLib
+{
+SaturationVanGenuchtenWithVolumetricStrain::
+ SaturationVanGenuchtenWithVolumetricStrain(
+ std::string name,
+ double const residual_liquid_saturation,
+ double const residual_gas_saturation,
+ double const exponent,
+ double const p_b,
+ double const e_0,
+ double const e_m,
+ double const a,
+ double const d_diff)
+ : S_L_res_(residual_liquid_saturation),
+ S_L_max_(1. - residual_gas_saturation),
+ m_(exponent),
+ p_b_(p_b),
+ e_0_(e_0),
+ e_m_(e_m),
+ a_(a),
+ d_diff_(d_diff)
+{
+ name_ = std::move(name);
+
+ if (!(m_ > 0 && m_ < 1))
+ {
+ OGS_FATAL("The exponent value m = {:g}, is out of its range of (0, 1)",
+ m_);
+ }
+}
+
+PropertyDataType SaturationVanGenuchtenWithVolumetricStrain::value(
+ VariableArray const& variable_array,
+ ParameterLib::SpatialPosition const& /*pos*/, double const /*t*/,
+ double const /*dt*/) const
+{
+ const double p_cap = variable_array.capillary_pressure;
+
+ if (p_cap <= 0)
+ {
+ return S_L_max_;
+ }
+
+ double const e_vol = variable_array.volumetric_strain;
+
+ double const n = 1. / (1. - m_);
+ double const d_e = -1 * (1 + e_0_) * e_vol / e_0_;
+ double const p_b_M = p_b_ * (1 / d_diff_);
+ double const p = p_cap / p_b_;
+ double const p_to_n = std::pow(p, n);
+ double const S_eff_mi = std::pow(p_to_n + 1., -m_);
+ double const p_M = p_cap / p_b_M;
+ double const p_to_n_M = std::pow(p_M, n);
+ double const S_eff_M = std::pow(p_to_n_M + 1., -m_);
+ double const S_eff = (((e_m_ + (a_ * d_e)) * (S_eff_mi)) +
+ ((e_0_ - e_m_ - (a_ * d_e)) * (S_eff_M))) /
+ e_0_;
+ double const S = S_eff * S_L_max_ - S_eff * S_L_res_ + S_L_res_;
+
+ return std::clamp(S, S_L_res_, S_L_max_);
+}
+
+PropertyDataType SaturationVanGenuchtenWithVolumetricStrain::dValue(
+ VariableArray const& variable_array, Variable const variable,
+ ParameterLib::SpatialPosition const& /*pos*/, double const /*t*/,
+ double const /*dt*/) const
+{
+ if (variable != Variable::capillary_pressure)
+ {
+ OGS_FATAL(
+ "SaturationVanGenuchtenWithVolumetricStrain::dValue is implemented "
+ "for derivatives with respect to capillary pressure only.");
+ }
+ const double p_cap = variable_array.capillary_pressure;
+
+ if (p_cap <= 0)
+ {
+ return 0.;
+ }
+
+ double const e_vol = variable_array.volumetric_strain;
+
+ double const n = 1. / (1. - m_);
+ double const d_e = -1 * (1 + e_0_) * e_vol / e_0_;
+ double const p_b_M = p_b_ * (1 / d_diff_);
+ double const p = p_cap / p_b_;
+ double const p_to_n = std::pow(p, n);
+ double const S_eff_mi = std::pow(p_to_n + 1., -m_);
+ double const p_M = p_cap / p_b_M;
+ double const p_to_n_M = std::pow(p_M, n);
+ double const S_eff_M = std::pow(p_to_n_M + 1., -m_);
+ double const S_eff = (((e_m_ + (a_ * d_e)) * (S_eff_mi)) +
+ ((e_0_ - e_m_ - (a_ * d_e)) * (S_eff_M))) /
+ e_0_;
+ double const S = S_eff * S_L_max_ - S_eff * S_L_res_ + S_L_res_;
+
+ if (S < S_L_res_ || S > S_L_max_)
+ {
+ return 0.;
+ }
+
+ double const dS_eff_dp_cap =
+ (-(e_m_ + (a_ * d_e)) * n * m_ * p_to_n *
+ std::pow(p_to_n + 1., -m_ - 1) -
+ (e_0_ - e_m_ - (a_ * d_e)) * n * m_ * p_to_n_M *
+ std::pow(p_to_n_M + 1., -m_ - 1)) /
+ (e_0_ * p_cap);
+ return dS_eff_dp_cap * (S_L_max_ - S_L_res_);
+}
+} // namespace MaterialPropertyLib
diff --git a/MaterialLib/MPL/Properties/CapillaryPressureSaturation/SaturationVanGenuchtenWithVolumetricStrain.h b/MaterialLib/MPL/Properties/CapillaryPressureSaturation/SaturationVanGenuchtenWithVolumetricStrain.h
new file mode 100644
index 0000000000000000000000000000000000000000..7cf978082b806c02162af6f488541cd3708f4871
--- /dev/null
+++ b/MaterialLib/MPL/Properties/CapillaryPressureSaturation/SaturationVanGenuchtenWithVolumetricStrain.h
@@ -0,0 +1,92 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2023, 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 "MaterialLib/MPL/Property.h"
+
+namespace MaterialPropertyLib
+{
+class Medium;
+class Phase;
+class Component;
+/**
+ * \brief A strain dependent bimodal water retention model.
+ *
+ * It is based on the van Genuchten model.
+ *
+ * The equation is as follows
+ * \f[S_{e} = \left((e_{m}+a \Delta e)
+ * \left[\frac{1}{((\frac{p_{c}}{p_{b}})^{n}+1)}\right]^{m}
+ * + (e_{M} - a \Delta e)
+ * \left[\frac{1}{((\frac{p_{c}}{p_{{b,M}}})^{n} +1)}\right]^{n}
+ * \right) \frac{1}{e}\f]
+ * with effective saturation defined as \f$S_{e}=\frac{S-S_r}{S_{max}-S_r}\f$,
+ * where \f$S_r\f$ and \f$S_{max}\f$ are the residual saturation and the maximum
+ * saturation, respectively.
+ * The exponent \f$m \in (0,1)\f$ is the same as in the van Genuchten equation (
+ * see SaturationVanGenuchten). In the original work another exponent \f$n\f$ is
+ * used, but usually set to \f$n = 1 / (1 - m)\f$, and also in this
+ * implementation.
+ * The pressure scaling parameter \f$p_{b}\f$ is added by the user and is the
+ * scaling parameter of the micropores. The scaling parameter of the macropores
+ * can be calculated as follows \f$p_{b,M} = p_{b} d_{diff}\f$ The total void
+ * ratio and the void ratio of the micropores are \f$e_0\f$ and \f$e_m\f$.
+ * Another scaling factor \f$a\f$ scales the effect of the strain
+ *
+ * The changing void ratio is calculated as
+ * \f$\Delta e = -\frac{(1-e)\epsilon_{vol}}{e}\f$,
+ * with \f$\epsilon_{vol}\f$ as the volumetric strain. The result is then
+ * clamped between the residual and maximum liquid saturations.
+ */
+class SaturationVanGenuchtenWithVolumetricStrain final : public Property
+{
+public:
+ SaturationVanGenuchtenWithVolumetricStrain(
+ std::string name,
+ double const residual_liquid_saturation,
+ double const residual_gas_saturation,
+ double const exponent,
+ double const p_b,
+ double const e_0,
+ double const e_m,
+ double const a,
+ double const d_diff);
+
+ void checkScale() const override
+ {
+ if (!std::holds_alternative<Medium*>(scale_))
+ {
+ OGS_FATAL(
+ "The property 'SaturationVanGenuchtenWithVolumetricStrain' is "
+ "implemented on the 'media' scale only.");
+ }
+ }
+
+ /// Those methods override the base class implementations and actually
+ /// compute and set the property values_ and dValues_.
+ PropertyDataType value(VariableArray const& variable_array,
+ ParameterLib::SpatialPosition const& pos,
+ double const t, double const dt) const override;
+ PropertyDataType dValue(VariableArray const& variable_array,
+ Variable const variable,
+ ParameterLib::SpatialPosition const& pos,
+ double const t, double const dt) const override;
+
+private:
+ double const S_L_res_;
+ double const S_L_max_;
+ double const m_;
+ double const p_b_;
+ double const e_0_;
+ double const e_m_;
+ double const a_;
+ double const d_diff_;
+};
+} // namespace MaterialPropertyLib
diff --git a/MaterialLib/MPL/Properties/CreateProperties.h b/MaterialLib/MPL/Properties/CreateProperties.h
index dfe4b2059fdc61321ff6ab4dd1157d38d754a465..5d16b8c14608c500df0aaf21000803885a2f3db8 100644
--- a/MaterialLib/MPL/Properties/CreateProperties.h
+++ b/MaterialLib/MPL/Properties/CreateProperties.h
@@ -17,6 +17,7 @@
#include "CapillaryPressureSaturation/CreateSaturationExponential.h"
#include "CapillaryPressureSaturation/CreateSaturationLiakopoulos.h"
#include "CapillaryPressureSaturation/CreateSaturationVanGenuchten.h"
+#include "CapillaryPressureSaturation/CreateSaturationVanGenuchtenWithVolumetricStrain.h"
#include "CreateAverageMolarMass.h"
#include "CreateBishopsPowerLaw.h"
#include "CreateBishopsSaturationCutoff.h"
diff --git a/MaterialLib/MPL/Properties/Properties.h b/MaterialLib/MPL/Properties/Properties.h
index 7324d4255d3eca480c3f1c0ae64a2d0a586a5bc3..dada0e0cfef4034cfefd906cf2bfa71632ac8134 100644
--- a/MaterialLib/MPL/Properties/Properties.h
+++ b/MaterialLib/MPL/Properties/Properties.h
@@ -18,6 +18,7 @@
#include "CapillaryPressureSaturation/SaturationExponential.h"
#include "CapillaryPressureSaturation/SaturationLiakopoulos.h"
#include "CapillaryPressureSaturation/SaturationVanGenuchten.h"
+#include "CapillaryPressureSaturation/SaturationVanGenuchtenWithVolumetricStrain.h"
#include "ClausiusClapeyron.h"
#include "Constant.h"
#include "Curve.h"
diff --git a/Tests/MaterialLib/TestMPLSaturationVanGenuchtenWithVolumetricStrain.cpp b/Tests/MaterialLib/TestMPLSaturationVanGenuchtenWithVolumetricStrain.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..489212b31b725def5421c110f212700f92860773
--- /dev/null
+++ b/Tests/MaterialLib/TestMPLSaturationVanGenuchtenWithVolumetricStrain.cpp
@@ -0,0 +1,77 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2023, 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 <gtest/gtest.h>
+
+#include "MaterialLib/MPL/Medium.h"
+#include "MaterialLib/MPL/Properties/CapillaryPressureSaturation/SaturationVanGenuchtenWithVolumetricStrain.h"
+#include "TestMPL.h"
+#include "Tests/TestTools.h"
+
+namespace MPL = MaterialPropertyLib;
+
+TEST(MaterialPropertyLib, SaturationVanGenuchtenWithVolumetricStrain)
+{
+ double const residual_liquid_saturation = 0.1;
+ double const residual_gas_saturation = 0.05;
+ double const exponent = 0.5;
+ double const p_b = 5000000;
+ double const e_0 = 0.7;
+ double const e_m = 0.5;
+ double const a = 8;
+ double const d_diff = 20;
+
+ MPL::Property const& pressure_saturation =
+ MPL::SaturationVanGenuchtenWithVolumetricStrain{
+ "saturation",
+ residual_liquid_saturation,
+ residual_gas_saturation,
+ exponent,
+ p_b,
+ e_0,
+ e_m,
+ a,
+ d_diff};
+
+ MPL::VariableArray variable_array;
+ ParameterLib::SpatialPosition const pos;
+ double const t = std::numeric_limits<double>::quiet_NaN();
+ double const dt = std::numeric_limits<double>::quiet_NaN();
+
+ /// No volumetric strain
+ {
+ double const vol_s = 0;
+ double const p_L = 1000000;
+
+ variable_array.capillary_pressure = p_L;
+ variable_array.volumetric_strain = vol_s;
+
+ double const S = pressure_saturation.template value<double>(
+ variable_array, pos, t, dt);
+
+ double const S_expected = 0.75426406;
+
+ ASSERT_NEAR(S, S_expected, 1e-5);
+ }
+ /// Volumetric strain
+ {
+ double const vol_s = 0.002;
+ double const p_L = 1000000;
+
+ variable_array.capillary_pressure = p_L;
+ variable_array.volumetric_strain = vol_s;
+
+ double const S = pressure_saturation.template value<double>(
+ variable_array, pos, t, dt);
+
+ double const S_expected = 0.71943039;
+
+ ASSERT_NEAR(S, S_expected, 1e-5);
+ }
+}