diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/CreateRelativePermeabilityModel.cpp b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/CreateRelativePermeabilityModel.cpp index 5fd34bc7eb693997f809eebf4fbd3b7c2ea9afd8..71aa2a46c21eff6608615807db64cefb0fbd2134 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/CreateRelativePermeabilityModel.cpp +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/CreateRelativePermeabilityModel.cpp @@ -34,97 +34,167 @@ namespace MaterialLib { namespace PorousMedium { -enum class BrookCoreyOrVanGenuchtenModelType +/** + \param config ConfigTree object which contains the input data + including <type>WettingPhaseVanGenuchten</type> + and it has a tag of <relative_permeability> +*/ +std::unique_ptr<RelativePermeability> createWettingPhaseVanGenuchten( + BaseLib::ConfigTree const& config) { - WettingPhaseVanGenuchten, - NonWettingPhaseVanGenuchten, - WettingPhaseBrookCoreyOilGas, - NonWettingPhaseBrookCoreyOilGas -}; + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type} + config.checkConfigParameter("type", "WettingPhaseVanGenuchten"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__WettingPhaseVanGenuchten__sr} + const double Sr = config.getConfigParameter<double>("sr"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__WettingPhaseVanGenuchten__smax} + const double Smax = config.getConfigParameter<double>("smax"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__WettingPhaseVanGenuchten__m} + const double m = config.getConfigParameter<double>("m"); + if (m > 1.0) // m <= 1 + { + OGS_FATAL( + "The exponent parameter of WettingPhaseVanGenuchten relative\n" + " permeability model, m, must be in an interval of [0, 1]"); + } + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__WettingPhaseVanGenuchten__m} + const double krel_max = config.getConfigParameter<double>("krel_min"); + + return std::unique_ptr<RelativePermeability>( + new WettingPhaseVanGenuchten(Sr, Smax, m, krel_max)); +} /** - * \param config ConfigTree object which contains the input data - * including <type>WettingPhaseBrookCoreyOilGas</type> - * or <type>NonWettingPhaseBrookCoreyOilGas</type> - * or <type>WettingPhaseVanGenuchten</type> - * or <type>NonWettingPhaseVanGenuchten</type> - * and it has a tag of <relative_permeability> + \param config ConfigTree object which contains the input data + including <type>NonWettingPhaseVanGenuchten</type> + and it has a tag of <relative_permeability> */ -static std::unique_ptr<RelativePermeability> createBrookCoreyOrVanGenuchten( - BaseLib::ConfigTree const& config, - const BrookCoreyOrVanGenuchtenModelType mode_type) +std::unique_ptr<RelativePermeability> createNonWettingPhaseVanGenuchten( + BaseLib::ConfigTree const& config) { - std::array<double, 4> parameters = { - {//! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type__sr} - config.getConfigParameter<double>("sr"), - //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type__smax} - config.getConfigParameter<double>("smax"), - //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type__m} - config.getConfigParameter<double>("m"), - //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type__krel_min} - config.getConfigParameter<double>("krel_min")}}; - - switch (mode_type) + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type} + config.checkConfigParameter("type", "NonWettingPhaseVanGenuchten"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__NonWettingPhaseVanGenuchten__sr} + const double Sr = config.getConfigParameter<double>("sr"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__NonWettingPhaseVanGenuchten__smax} + const double Smax = config.getConfigParameter<double>("smax"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__NonWettingPhaseVanGenuchten__m} + const double m = config.getConfigParameter<double>("m"); + if (m > 1.0) // m <= 1 { - case BrookCoreyOrVanGenuchtenModelType::WettingPhaseVanGenuchten: - assert(parameters[2] <= 1.0); // m <= 1 - return std::unique_ptr<RelativePermeability>( - new WettingPhaseVanGenuchten(parameters)); - case BrookCoreyOrVanGenuchtenModelType::NonWettingPhaseVanGenuchten: - assert(parameters[2] <= 1.0); // m <= 1 - return std::unique_ptr<RelativePermeability>( - new NonWettingPhaseVanGenuchten(parameters)); - case BrookCoreyOrVanGenuchtenModelType::WettingPhaseBrookCoreyOilGas: - assert(parameters[2] >= 1.0); // m >= 1 - return std::unique_ptr<RelativePermeability>( - new WettingPhaseBrookCoreyOilGas(parameters)); - case BrookCoreyOrVanGenuchtenModelType::NonWettingPhaseBrookCoreyOilGas: - assert(parameters[2] >= 1.0); // m >= 1 - return std::unique_ptr<RelativePermeability>( - new NonWettingPhaseBrookCoreyOilGas(parameters)); + OGS_FATAL( + "The exponent parameter of NonWettingPhaseVanGenuchten relative\n" + " permeability model, m, must be in an interval of [0, 1]"); } + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__NonWettingPhaseVanGenuchten__m} + const double krel_max = config.getConfigParameter<double>("krel_min"); - return nullptr; + return std::unique_ptr<RelativePermeability>( + new NonWettingPhaseVanGenuchten(Sr, Smax, m, krel_max)); +} + +/** + \param config ConfigTree object which contains the input data + including <type>WettingPhaseBrookCoreyOilGas</type> + and it has a tag of <relative_permeability> +*/ +std::unique_ptr<RelativePermeability> createWettingPhaseBrookCoreyOilGas( + BaseLib::ConfigTree const& config) +{ + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type} + config.checkConfigParameter("type", "WettingPhaseBrookCoreyOilGas"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__WettingPhaseBrookCoreyOilGas__sr} + const double Sr = config.getConfigParameter<double>("sr"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__WettingPhaseBrookCoreyOilGas__smax} + const double Smax = config.getConfigParameter<double>("smax"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__WettingPhaseBrookCoreyOilGas__m} + const double m = config.getConfigParameter<double>("m"); + if (m < 1.0) // m >= 1 + { + OGS_FATAL( + "The exponent parameter of WettingPhaseBrookCoreyOilGas\n" + "relative permeability model, m, must not be smaller than 1"); + } + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__WettingPhaseBrookCoreyOilGas__m} + const double krel_max = config.getConfigParameter<double>("krel_min"); + + return std::unique_ptr<RelativePermeability>( + new WettingPhaseBrookCoreyOilGas(Sr, Smax, m, krel_max)); +} + +/** + \param config ConfigTree object which contains the input data + including <type>NonWettingPhaseBrookCoreyOilGas</type> + and it has a tag of <relative_permeability> +*/ +std::unique_ptr<RelativePermeability> createNonWettingPhaseBrookCoreyOilGas( + BaseLib::ConfigTree const& config) +{ + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type} + config.checkConfigParameter("type", "NonWettingPhaseBrookCoreyOilGas"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__NonWettingPhaseBrookCoreyOilGas__sr} + const double Sr = config.getConfigParameter<double>("sr"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__NonWettingPhaseBrookCoreyOilGas__smax} + const double Smax = config.getConfigParameter<double>("smax"); + + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__NonWettingPhaseBrookCoreyOilGas__m} + const double m = config.getConfigParameter<double>("m"); + if (m < 1.0) // m >= 1 + { + OGS_FATAL( + "The exponent parameter of NonWettingPhaseBrookCoreyOilGas\n" + "relative permeability model, m, must not be smaller than 1"); + } + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__NonWettingPhaseBrookCoreyOilGas__m} + const double krel_max = config.getConfigParameter<double>("krel_min"); + + return std::unique_ptr<RelativePermeability>( + new NonWettingPhaseBrookCoreyOilGas(Sr, Smax, m, krel_max)); } std::unique_ptr<RelativePermeability> createRelativePermeabilityModel( BaseLib::ConfigTree const& config) { //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type} - auto const type = config.getConfigParameter<std::string>("type"); + auto const type = config.peekConfigParameter<std::string>("type"); if (type == "WettingPhaseVanGenuchten") { - return createBrookCoreyOrVanGenuchten( - config, - BrookCoreyOrVanGenuchtenModelType::WettingPhaseVanGenuchten); + return createWettingPhaseVanGenuchten(config); } else if (type == "NonWettingPhaseVanGenuchten") { - return createBrookCoreyOrVanGenuchten( - config, - BrookCoreyOrVanGenuchtenModelType::NonWettingPhaseVanGenuchten); + return createNonWettingPhaseVanGenuchten(config); } else if (type == "WettingPhaseBrookCoreyOilGas") { - return createBrookCoreyOrVanGenuchten( - config, - BrookCoreyOrVanGenuchtenModelType::WettingPhaseBrookCoreyOilGas); + return createWettingPhaseBrookCoreyOilGas(config); } else if (type == "NonWettingPhaseBrookCoreyOilGas") { - return createBrookCoreyOrVanGenuchten( - config, - BrookCoreyOrVanGenuchtenModelType::NonWettingPhaseBrookCoreyOilGas); + return createNonWettingPhaseBrookCoreyOilGas(config); } else if (type == "Curve") { + config.checkConfigParameter("type", "Curve"); + std::vector<double> variables, values; - //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type__curve} + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__Curve__curve} auto const& curve_config = config.getConfigSubtree("curve"); - //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__type__curve__data} - for (auto const& data_string : - curve_config.getConfigParameterList<std::string>("data")) + for ( + auto const& data_string : + //! \ogs_file_param{material_property__porous_medium__porous_medium__relative_permeability__Curve_curve__data} + curve_config.getConfigParameterList<std::string>("data")) { std::stringstream ss(data_string); double var, val; @@ -142,12 +212,13 @@ std::unique_ptr<RelativePermeability> createRelativePermeabilityModel( else { OGS_FATAL( - "The capillary pressure model %s is unavailable.\n" + "The relative permeability model %s is unavailable.\n" "The available models are:" "\n\tWettingPhaseVanGenuchten," "\n\tNonWettingPhaseVanGenuchten," "\n\tWettingPhaseBrookCoreyOilGas," - "\n\tNonWettingPhaseBrookCoreyOilGas.\n", + "\n\tNonWettingPhaseBrookCoreyOilGas,", + "\n\tCurve.\n", type.data()); } } diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseBrookCoreyOilGas.cpp b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseBrookCoreyOilGas.cpp index d0ed9396204755be745c51702c3c131603a950a7..67c0377484a20126b89a7cb9a5a242f382aedf94 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseBrookCoreyOilGas.cpp +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseBrookCoreyOilGas.cpp @@ -26,9 +26,9 @@ double NonWettingPhaseBrookCoreyOilGas::getValue( const double S = MathLib::limitValueInInterval( saturation_w, _Sr + _minor_offset, _Smax - _minor_offset); const double Se = (S - _Sr) / (_Smax - _Sr); - const double Krel = + const double krel = std::pow(1.0 - Se, 2) * (1.0 - std::pow(Se, 1.0 + 2.0 / _mm)); - return Krel < _Krel_min ? _Krel_min : Krel; + return krel < _krel_min ? _krel_min : krel; } } // end namespace diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseBrookCoreyOilGas.h b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseBrookCoreyOilGas.h index 78b0f2a69c098aafa814d09b9f917a71ee003a9d..fbbc8afc1fe84dc71e5eaaa33ccab41489cb2501 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseBrookCoreyOilGas.h +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseBrookCoreyOilGas.h @@ -13,9 +13,6 @@ #ifndef OGS_NON_WETTING_PHASE_BROOK_COREY_OIL_GAS_H #define OGS_NON_WETTING_PHASE_BROOK_COREY_OIL_GAS_H -#include <array> -#include <limits> // std::numeric_limits - #include "RelativePermeability.h" namespace MaterialLib @@ -27,28 +24,29 @@ namespace PorousMedium * * \f[k{rel}= (1-S_e)^2 (1 - S_e^{1 + 2 /m})\f] * with - * \f[S_e=\dfrac{S-S_r}{S_{\mbox{max}}-S_r}\f] + * \f[S_e=\dfrac{S^w-S_r}{S^w_{\mbox{max}}-S_r}\f] * where * \f{eqnarray*}{ - * &S_r& \mbox{ residual saturation,}\\ - * &S_{\mbox{max}}& \mbox{ maximum saturation,}\\ + * &S^w_r& \mbox{residual saturation of wetting phase,}\\ + * &S^w_{\mbox{max}}& \mbox{maximum saturation of wetting phase,}\\ * &m(>=1) & \mbox{ exponent.}\\ * \f} */ class NonWettingPhaseBrookCoreyOilGas final : public RelativePermeability { public: - /** \param parameters An array contains the five parameters: - * [0] \f$ S_{nr} \f$ - * [1] \f$ S_{n\mbox{max}} \f$ - * [2] \f$ m \f$ - * [3] \f$ K_{rel}^{\mbox{min}}\f$ + /** + * @param Snr Residual saturation of the non-wetting phase, + * \f$ S^n_r \f$ + * @param Snmax Maximum saturation of the non-wetting phase, + * \f$ S^n_{\mbox{max}} \f$ + * @param m Exponent, \f$ m \f$ + * @param krel_min Minimum relative permeability, + * \f$ k_{rel}^{\mbox{min}}\f$ */ - NonWettingPhaseBrookCoreyOilGas(std::array<double, 4> const& parameters) - : _Sr(1. - parameters[1]), - _Smax(1. - parameters[0]), - _mm(parameters[2]), - _Krel_min(parameters[3]) + NonWettingPhaseBrookCoreyOilGas(const double Snr, const double Snmax, + const double m, const double krel_min) + : _Sr(1. - Snmax), _Smax(1. - Snr), _mm(m), _krel_min(krel_min) { } @@ -63,10 +61,10 @@ public: double getValue(const double saturation_w) const override; private: - const double _Sr; ///< Residual saturation of wetting phase, 1-Snr. - const double _Smax; ///< Maximum saturation of wetting phase., 1-Sn_max + const double _Sr; ///< Residual saturation of wetting phase, 1-Snmax. + const double _Smax; ///< Maximum saturation of wetting phase, 1-Snr. const double _mm; ///< Exponent (>=1.0), n=1/(1-mm). - const double _Krel_min; ///< Minimum relative permeability + const double _krel_min; ///< Minimum relative permeability }; } // end namespace diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseVanGenuchten.cpp b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseVanGenuchten.cpp index 8fe157012660c2715a9638c9b52c7ca317819188..fc4a3df74aebc074e5a13454dca451d44153fc60 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseVanGenuchten.cpp +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseVanGenuchten.cpp @@ -25,9 +25,9 @@ double NonWettingPhaseVanGenuchten::getValue(const double saturation_w) const const double S = MathLib::limitValueInInterval( saturation_w, _Sr + _minor_offset, _Smax - _minor_offset); const double Se = (S - _Sr) / (_Smax - _Sr); - const double Krel = std::pow(1.0 - Se, 1.0 / 3.0) * + const double krel = std::pow(1.0 - Se, 1.0 / 3.0) * std::pow(1.0 - std::pow(Se, 1.0 / _mm), 2.0 * _mm); - return Krel < _Krel_min ? _Krel_min : Krel; + return krel < _krel_min ? _krel_min : krel; } } // end namespace } // end namespace diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseVanGenuchten.h b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseVanGenuchten.h index 0314fb22350c258aa8024d4ac8768a3048607c4a..91e1c3dfa2ddc167e53af9907ece2a80f3627c4a 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseVanGenuchten.h +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/NonWettingPhaseVanGenuchten.h @@ -12,9 +12,6 @@ #ifndef OGS_NON_WETTING_PHASE_VAN_GENUCHTEN_H #define OGS_NON_WETTING_PHASE_VAN_GENUCHTEN_H -#include <array> -#include <limits> // std::numeric_limits - #include "RelativePermeability.h" namespace MaterialLib @@ -26,12 +23,12 @@ namespace PorousMedium * * \f[k{rel}= (1 - S_e)^{1/3} (1 - S_e^{1/m})^{2m}\f] * with - * \f[S_e=\dfrac{S-S_r}{S_{\mbox{max}}-S_r}\f] + * \f[S_e=\dfrac{S^w-S_r}{S^w_{\mbox{max}}-S^w_r}\f] * where * \f{eqnarray*}{ - * &S_r& \mbox{ residual saturation,}\\ - * &S_{\mbox{max}}& \mbox{ maximum saturation,}\\ - * &m(<=1) & \mbox{ exponent.}\\ + * &S^w_r& \mbox{residual saturation of wetting phase,}\\ + * &S^w_{\mbox{max}}& \mbox{maximum saturation of wetting phase,}\\ + * &m(<=1) & \mbox{ exponent.}\\ * \f} * * Note: @@ -40,24 +37,25 @@ namespace PorousMedium class NonWettingPhaseVanGenuchten final : public RelativePermeability { public: - /** \param parameters An array contains the five parameters: - * [0] \f$ S_{nr} \f$ - * [1] \f$ S_{n\mbox{max}} \f$ - * [2] \f$ m \f$ - * [3] \f$ K_{rel}^{\mbox{min}}\f$ + /** + * @param Snr Residual saturation of the non-wetting phase, + * \f$ S^n_r \f$ + * @param Snmax Maximum saturation of the non-wetting phase, + * \f$ S^n_{\mbox{max}} \f$ + * @param m Exponent, \f$ m \in [0,1]\f$ + * @param krel_min Minimum relative permeability, + * \f$ k_{rel}^{\mbox{min}}\f$ */ - NonWettingPhaseVanGenuchten(std::array<double, 4> const& parameters) - : _Sr(1. - parameters[1]), - _Smax(1. - parameters[0]), - _mm(parameters[2]), - _Krel_min(parameters[3]) + NonWettingPhaseVanGenuchten(const double Snr, const double Snmax, + const double m, const double krel_min) + : _Sr(1. - Snmax), _Smax(1. - Snr), _mm(m), _krel_min(krel_min) { } /// Get model name. std::string getName() const override { - return "Non-wetting phase van Genuchten model."; + return "Non-wetting phase van Genuchten relative permeability model."; } /// Get relative permeability value. @@ -65,12 +63,10 @@ public: double getValue(const double saturation_w) const override; private: - const double _Sr; ///< Residual saturation of wetting phase, 1-Snr. - const double _Smax; ///< Maximum saturation of wetting phase., 1-Sn_max + const double _Sr; ///< Residual saturation of wetting phase, 1-Snmax. + const double _Smax; ///< Maximum saturation of wetting phase, 1-Snr. const double _mm; ///< Exponent (<=1.0), n=1/(1-mm). - const double _Krel_min; ///< Minimum relative permeability - - const double _perturbation = std::numeric_limits<double>::epsilon(); + const double _krel_min; ///< Minimum relative permeability }; } // end namespace diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/RelativePermeability.h b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/RelativePermeability.h index 73dd98264149ded1faa0ceea001f2f719754a021..ace2c329dc87bef48187d14de088af460242304f 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/RelativePermeability.h +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/RelativePermeability.h @@ -13,6 +13,7 @@ #define OGS_RELATIVE_PERMEABILITY_H #include <string> +#include <limits> namespace MaterialLib { diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseBrookCoreyOilGas.cpp b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseBrookCoreyOilGas.cpp index 1b8703c6775ad796645c74d4376dd5306bb2a1ea..938bf57d2505b8e040b25ee297bddc473b595a57 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseBrookCoreyOilGas.cpp +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseBrookCoreyOilGas.cpp @@ -25,8 +25,8 @@ double WettingPhaseBrookCoreyOilGas::getValue(const double saturation) const const double S = MathLib::limitValueInInterval( saturation, _Sr + _minor_offset, _Smax - _minor_offset); const double Se = (S - _Sr) / (_Smax - _Sr); - const double Krel = std::pow(Se, 3.0 + 2.0 / _mm); - return Krel < _Krel_min ? _Krel_min : Krel; + const double krel = std::pow(Se, 3.0 + 2.0 / _mm); + return krel < _krel_min ? _krel_min : krel; } } // end namespace diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseBrookCoreyOilGas.h b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseBrookCoreyOilGas.h index f9732c43bd084526c98487977f6b063e1215055e..28b095d9982e1670e7bd985352e0807f768df88a 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseBrookCoreyOilGas.h +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseBrookCoreyOilGas.h @@ -13,9 +13,6 @@ #ifndef OGS_WETTING_PHASE_BROOK_COREY_OIL_GAS_H #define OGS_WETTING_PHASE_BROOK_COREY_OIL_GAS_H -#include <array> -#include <limits> // std::numeric_limits - #include "RelativePermeability.h" namespace MaterialLib @@ -38,17 +35,15 @@ namespace PorousMedium class WettingPhaseBrookCoreyOilGas final : public RelativePermeability { public: - /** \param parameters An array contains the five parameters: - * [0] \f$ S_{r} \f$ - * [1] \f$ S_{\mbox{max}} \f$ - * [2] \f$ m \f$ - * [3] \f$ K_{rel}^{\mbox{min}}\f$ + /** + * @param Sr Residual saturation, \f$ S_r \f$ + * @param Smax Maximum saturation, \f$ S_{\mbox{max}} \f$ + * @param m Exponent, \f$ m \f$ + * @param krel_min Minimum relative permeability,\f$ k_{rel}^{\mbox{min}}\f$ */ - WettingPhaseBrookCoreyOilGas(std::array<double, 4> const& parameters) - : _Sr(parameters[0]), - _Smax(parameters[1]), - _mm(parameters[2]), - _Krel_min(parameters[3]) + WettingPhaseBrookCoreyOilGas(const double Sr, const double Smax, + const double m, const double krel_min) + : _Sr(Sr), _Smax(Smax), _mm(m), _krel_min(krel_min) { } @@ -65,7 +60,7 @@ private: const double _Sr; ///< Residual saturation. const double _Smax; ///< Maximum saturation. const double _mm; ///< Exponent (>=1.0), n=1/(1-mm). - const double _Krel_min; ///< Minimum relative permeability + const double _krel_min; ///< Minimum relative permeability }; } // end namespace diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseVanGenuchten.cpp b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseVanGenuchten.cpp index 56532c924d8b00583760e8be7174d2477f8d399e..f57b66716b8e3fce42369eb1736a5c8a797acedd 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseVanGenuchten.cpp +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseVanGenuchten.cpp @@ -25,10 +25,10 @@ double WettingPhaseVanGenuchten::getValue(const double saturation) const const double S = MathLib::limitValueInInterval( saturation, _Sr + _minor_offset, _Smax - _minor_offset); const double Se = (S - _Sr) / (_Smax - _Sr); - const double Krel = + const double krel = std::sqrt(Se) * std::pow(1.0 - std::pow(1.0 - std::pow(Se, 1.0 / _mm), _mm), 2); - return Krel < _Krel_min ? _Krel_min : Krel; + return krel < _krel_min ? _krel_min : krel; } } // end namespace diff --git a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseVanGenuchten.h b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseVanGenuchten.h index db0a495845aa8e4252c0fbb70a919c685989eb50..8fd580f378e0bfa0c03dd284f5aa630de67e7b65 100644 --- a/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseVanGenuchten.h +++ b/MaterialLib/PorousMedium/UnsaturatedProperty/RelativePermeability/WettingPhaseVanGenuchten.h @@ -13,9 +13,6 @@ #ifndef OGS_WETTING_PHASE_VAN_GENUCHTEN_H #define OGS_WETTING_PHASE_VAN_GENUCHTEN_H -#include <array> -#include <limits> // std::numeric_limits - #include "RelativePermeability.h" namespace MaterialLib @@ -41,24 +38,22 @@ namespace PorousMedium class WettingPhaseVanGenuchten final : public RelativePermeability { public: - /** \param parameters An array contains the five parameters: - * [0] \f$ S_{r} \f$ - * [1] \f$ S_{\mbox{max}} \f$ - * [2] \f$ m \f$ - * [3] \f$ K_{rel}^{\mbox{min}}\f$ + /** + * @param Sr Residual saturation, \f$ S_r \f$ + * @param Smax Maximum saturation, \f$ S_{\mbox{max}} \f$ + * @param m Exponent, \f$ m \f$ + * @param krel_min Minimum relative permeability,\f$ k_{rel}^{\mbox{min}}\f$ */ - WettingPhaseVanGenuchten(std::array<double, 4> const& parameters) - : _Sr(parameters[0]), - _Smax(parameters[1]), - _mm(parameters[2]), - _Krel_min(parameters[3]) + WettingPhaseVanGenuchten(const double Sr, const double Smax, const double m, + const double krel_min) + : _Sr(Sr), _Smax(Smax), _mm(m), _krel_min(krel_min) { } /// Get model name. std::string getName() const override { - return "Wetting phase van Genuchten model."; + return "Wetting phase van Genuchten relative permeability model."; } /// Get relative permeability value. @@ -68,7 +63,7 @@ private: const double _Sr; ///< Residual saturation. const double _Smax; ///< Maximum saturation. const double _mm; ///< Exponent (<=1.0), n=1/(1-mm). - const double _Krel_min; ///< Minimum relative permeability + const double _krel_min; ///< Minimum relative permeability }; } // end namespace diff --git a/Tests/MaterialLib/TestRelativePermeabilityModel.cpp b/Tests/MaterialLib/TestRelativePermeabilityModel.cpp index e5cae9a96bc5dcb789a1c089fa5edb02fc2b9167..5207c8486cbb47fa3cc4183fc31f2cf6c8735613 100644 --- a/Tests/MaterialLib/TestRelativePermeabilityModel.cpp +++ b/Tests/MaterialLib/TestRelativePermeabilityModel.cpp @@ -52,7 +52,7 @@ TEST(MaterialPorousMedium, checkWettingPhaseVanGenuchten) auto const perm_model = createRelativePermeabilityModel(xml); std::vector<double> S = {0.2, 0.33, 0.45, 0.52, 0.6, 0.85}; - std::vector<double> Krel = {1.e-9, + std::vector<double> krel = {1.e-9, 1.9291770102827742e-006, 0.00041114113180510661, 0.0019569840357319015, @@ -61,7 +61,7 @@ TEST(MaterialPorousMedium, checkWettingPhaseVanGenuchten) for (std::size_t i = 0; i < S.size(); i++) { - ASSERT_NEAR(Krel[i], perm_model->getValue(S[i]), 1.e-5); + ASSERT_NEAR(krel[i], perm_model->getValue(S[i]), 1.e-9); } } @@ -78,13 +78,13 @@ TEST(MaterialPorousMedium, checkNonWettingPhaseVanGenuchten) auto const perm_model = createRelativePermeabilityModel(xml); std::vector<double> S = {0.2, 0.33, 0.45, 0.52, 0.6, 0.85}; - std::vector<double> Krel = {0.87422700239237161, 0.74331414436457388, + std::vector<double> krel = {0.87422700239237161, 0.74331414436457388, 0.59527539448807487, 0.49976666464188485, 0.38520070797257489, 0.041219134248319585}; for (std::size_t i = 0; i < S.size(); i++) { - ASSERT_NEAR(Krel[i], perm_model->getValue(S[i]), 1.e-5); + ASSERT_NEAR(krel[i], perm_model->getValue(S[i]), 1.e-9); } } @@ -101,7 +101,7 @@ TEST(MaterialPorousMedium, checkWettingPhaseBrookCoreyOilGas) auto const perm_model = createRelativePermeabilityModel(xml); std::vector<double> S = {0.2, 0.33, 0.45, 0.52, 0.6, 0.85}; - std::vector<double> Krel = {0., + std::vector<double> krel = {0., 0.0022037808641975302, 0.030140817901234556, 0.080908641975308573, @@ -110,7 +110,7 @@ TEST(MaterialPorousMedium, checkWettingPhaseBrookCoreyOilGas) for (std::size_t i = 0; i < S.size(); i++) { - ASSERT_NEAR(Krel[i], perm_model->getValue(S[i]), 1.e-5); + ASSERT_NEAR(krel[i], perm_model->getValue(S[i]), 1.e-9); } } @@ -127,7 +127,7 @@ TEST(MaterialPorousMedium, checkNonWettingPhaseBrookCoreyOilGas) auto const perm_model = createRelativePermeabilityModel(xml); std::vector<double> S = {0.2, 0.33, 0.45, 0.52, 0.6, 0.85}; - std::vector<double> Krel = {1.0, + std::vector<double> krel = {1.0, 0.58480547839506147, 0.28120177469135793, 0.15583209876543211, @@ -136,7 +136,7 @@ TEST(MaterialPorousMedium, checkNonWettingPhaseBrookCoreyOilGas) for (std::size_t i = 0; i < S.size(); i++) { - ASSERT_NEAR(Krel[i], perm_model->getValue(S[i]), 1.e-5); + ASSERT_NEAR(krel[i], perm_model->getValue(S[i]), 1.e-9); } } @@ -154,10 +154,10 @@ TEST(MaterialPorousMedium, checkReletivePermeabilityCurve) auto const perm_model = createRelativePermeabilityModel(xml); std::vector<double> S = {0.2, 0.33, 0.45, 0.52, 0.6, 0.85}; - std::vector<double> Krel = {0.7, 0.57, 0.451, 0.3824, 0.304, 0.059}; + std::vector<double> krel = {0.7, 0.57, 0.451, 0.3824, 0.304, 0.059}; for (std::size_t i = 0; i < S.size(); i++) { - ASSERT_NEAR(Krel[i], perm_model->getValue(S[i]), 1.e-5); + ASSERT_NEAR(krel[i], perm_model->getValue(S[i]), 1.e-9); } }