13#ifndef DUMUX_MATERIAL_FLUIDMATRIX_TWOP_MATERIAL_LAW_HH 
   14#define DUMUX_MATERIAL_FLUIDMATRIX_TWOP_MATERIAL_LAW_HH 
   40template<
class ScalarType,
 
   50    using BasicParams = 
typename BaseLaw::template Params<Scalar>;
 
   66    { 
return !std::is_same<Regularization, NoRegularization>::value; }
 
 
   83            regularization_.init(
this, paramGroup);
 
 
   93    : basicParams_(baseParams)
 
 
  103    template<
bool enableRegularization = isRegularized()>
 
  107        if constexpr (enableRegularization)
 
  109            const auto regularized = regularization_.pc(swe);
 
  111                return regularized.value();
 
  114        return BaseLaw::pc(swe, basicParams_);
 
 
  120    template<
bool enableRegularization = isRegularized()>
 
  124        if constexpr (enableRegularization)
 
  126            const auto regularized = regularization_.dpc_dswe(swe);
 
 
  139        return BaseLaw::endPointPc(basicParams_);
 
 
  145    template<
bool enableRegularization = isRegularized()>
 
  148        if constexpr (enableRegularization)
 
  150            const auto regularized = regularization_.swe(
pc);
 
 
  161    template<
bool enableRegularization = isRegularized()>
 
  164        if constexpr (enableRegularization)
 
  166            const auto regularized = regularization_.dswe_dpc(
pc);
 
 
  177    template<
bool enableRegularization = isRegularized()>
 
  181        if constexpr (enableRegularization)
 
  183            const auto regularized = regularization_.krw(swe);
 
  185                return regularized.value();
 
  188        return BaseLaw::krw(swe, basicParams_);
 
 
  194    template<
bool enableRegularization = isRegularized()>
 
  198        if constexpr (enableRegularization)
 
  200            const auto regularized = regularization_.dkrw_dswe(swe);
 
 
  211    template<
bool enableRegularization = isRegularized()>
 
  215        if constexpr (enableRegularization)
 
  217            const auto regularized = regularization_.krn(swe);
 
  219                return regularized.value();
 
  222        return BaseLaw::krn(swe, basicParams_);
 
 
  228    template<
bool enableRegularization = isRegularized()>
 
  232        if constexpr (enableRegularization)
 
  234            const auto regularized = regularization_.dkrn_dswe(swe);
 
 
  247        return basicParams_ == o.basicParams_
 
  248               && effToAbsParams_ == o.effToAbsParams_
 
  249               && regularization_ == o.regularization_;
 
 
  258        return BaseLaw::template makeParams<Scalar>(paramGroup);
 
 
  265    { 
return basicParams_; }
 
 
  273        return EffToAbs::template makeParams<Scalar>(paramGroup);
 
 
  280    { 
return effToAbsParams_; }
 
 
  285    Regularization regularization_;
 
 
This is a policy for 2p material laws how to convert absolute to relative saturations and vice versa.
Definition efftoabsdefaultpolicy.hh:34
static Scalar swToSwe(const Scalar sw, const Params< Scalar > ¶ms)
Convert an absolute wetting saturation to an effective one.
Definition efftoabsdefaultpolicy.hh:106
static Scalar sweToSw(const Scalar swe, const Params< Scalar > ¶ms)
Convert an effective wetting saturation to an absolute one.
Definition efftoabsdefaultpolicy.hh:121
static Scalar dswe_dsw(const Params< Scalar > ¶ms)
Derivative of the effective saturation w.r.t. the absolute saturation.
Definition efftoabsdefaultpolicy.hh:135
static Scalar dsw_dswe(const Params< Scalar > ¶ms)
Derivative of the absolute saturation w.r.t. the effective saturation.
Definition efftoabsdefaultpolicy.hh:149
TwoPMaterialLaw()=delete
Deleted default constructor (so we are never in an undefined state)
static constexpr int numFluidPhases()
Return the number of fluid phases.
Definition materiallaw.hh:59
Scalar endPointPc() const
The capillary pressure at Swe = 1.0 also called end point capillary pressure.
Definition materiallaw.hh:137
Scalar dsw_dpc(const Scalar pc) const
The partial derivative of the saturation to the capillary pressure.
Definition materiallaw.hh:162
TwoPMaterialLaw(const std::string ¶mGroup)
Construct from a subgroup from the global parameter tree.
Definition materiallaw.hh:78
TwoPMaterialLaw(const BasicParams &baseParams, const EffToAbsParams &effToAbsParams={}, const RegularizationParams ®Params={})
Construct from parameter structs.
Definition materiallaw.hh:90
static constexpr bool isRegularized()
Return whether this law is regularized.
Definition materiallaw.hh:65
Scalar dkrn_dsw(const Scalar sw) const
The derivative of the relative permeability for the non-wetting phase w.r.t. saturation.
Definition materiallaw.hh:229
static EffToAbsParams makeEffToAbsParams(const std::string ¶mGroup)
Definition materiallaw.hh:271
Scalar Scalar
Definition materiallaw.hh:48
Scalar krn(const Scalar sw) const
The relative permeability for the non-wetting phase.
Definition materiallaw.hh:212
const BasicParams & basicParams() const
Return the base law's parameters.
Definition materiallaw.hh:264
static BasicParams makeBasicParams(const std::string ¶mGroup)
Definition materiallaw.hh:256
Scalar dkrw_dsw(const Scalar sw) const
The derivative of the relative permeability for the wetting phase w.r.t. saturation.
Definition materiallaw.hh:195
typename BrooksCoreyRegularization< Scalar >::template Params< Scalar > RegularizationParams
Definition materiallaw.hh:52
const EffToAbsParams & effToAbsParams() const
Definition materiallaw.hh:279
bool operator==(const TwoPMaterialLaw &o) const
Equality comparison with another instance.
Definition materiallaw.hh:245
Scalar dpc_dsw(const Scalar sw) const
The partial derivative of the capillary pressure w.r.t. the saturation.
Definition materiallaw.hh:121
TwoPEffToAbsDefaultPolicy EffToAbs
Definition materiallaw.hh:54
Scalar sw(const Scalar pc) const
Definition materiallaw.hh:146
typename TwoPEffToAbsDefaultPolicy::template Params< Scalar > EffToAbsParams
Definition materiallaw.hh:51
typename BrooksCorey::template Params< Scalar > BasicParams
Definition materiallaw.hh:50
Scalar pc(const Scalar sw) const
The capillary pressure-saturation curve.
Definition materiallaw.hh:104
Scalar krw(const Scalar sw) const
The relative permeability for the wetting phase.
Definition materiallaw.hh:178
This is a policy for 2p material laws how to convert absolute to relative saturations and vice versa.
Wrapper type to combine an arbitrary number of different laws for fluid-matrix interaction (e....
Definition brookscorey.hh:23
A tag to turn off regularization and it's overhead.
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
Adapter to inherit from, allowing the inheriting class to be wrapped by the makeFluidMatrixInteractio...
Definition fluidmatrixinteraction.hh:56
A tag to turn off regularization and it's overhead.
Definition noregularization.hh:22