12#ifndef DUMUX_MATERIAL_FLUIDMATRIX_TWOP_SPLINE_MATERIAL_LAW_HH 
   13#define DUMUX_MATERIAL_FLUIDMATRIX_TWOP_SPLINE_MATERIAL_LAW_HH 
   34template<
class TwoPMaterialLaw, 
bool approximatePcSwInverse = false>
 
   37, 
public Adapter<SplineTwoPMaterialLaw<TwoPMaterialLaw, approximatePcSwInverse>, PcKrSw>
 
   72        const std::array<Scalar, 2> defaultInterval{{ 0.01, 1.0 }};
 
   74        swInterval_ = {{ TwoPMaterialLaw::EffToAbs::sweToSw(sweInterval[0], this->
effToAbsParams()),
 
   75                         TwoPMaterialLaw::EffToAbs::sweToSw(sweInterval[1], this->
effToAbsParams()) }};
 
   80        pcSpline_ = makeSweSpline_(
 
   82            approximatePcSwInverse
 
   85        krwSpline_ = makeSweSpline_(
 
   89        krnSpline_ = makeSweSpline_(
 
 
   99                          std::size_t numSwSamples,
 
  102    , numSwSamples_(numSwSamples)
 
  104        swInterval_ = {{ TwoPMaterialLaw::EffToAbs::sweToSw(sweInterval[0], this->
effToAbsParams()),
 
  105                         TwoPMaterialLaw::EffToAbs::sweToSw(sweInterval[1], this->
effToAbsParams()) }};
 
 
  115        if (
sw > swInterval_[0] && 
sw < swInterval_[1])
 
  117            if constexpr (approximatePcSwInverse)
 
  118                return pcSpline_->evalInverse(
sw);
 
  120                return pcSpline_->eval(
sw);
 
 
  131        if (
sw > swInterval_[0] && 
sw < swInterval_[1])
 
  133            if constexpr (approximatePcSwInverse)
 
  134                return 1.0/pcSpline_->evalDerivative(pcSpline_->evalInverse(
sw));
 
  136                return pcSpline_->evalDerivative(
sw);
 
 
  147        if (
pc > swIntervalPc_[0] && 
pc < swIntervalPc_[1])
 
  149            if constexpr (approximatePcSwInverse)
 
  150                return pcSpline_->eval(
pc);
 
  152                return pcSpline_->evalInverse(
pc);
 
 
  163        if (
pc > swIntervalPc_[0] && 
pc < swIntervalPc_[1])
 
  165            if constexpr (approximatePcSwInverse)
 
  166                return pcSpline_->evalDerivative(
pc);
 
  168                return 1.0/pcSpline_->evalDerivative(pcSpline_->evalInverse(
pc));
 
 
  179        if (
sw > swInterval_[0] && 
sw < swInterval_[1])
 
  180            return krwSpline_->eval(
sw);
 
 
  190        if (
sw > swInterval_[0] && 
sw < swInterval_[1])
 
  191            return krwSpline_->evalDerivative(
sw);
 
 
  201        if (
sw > swInterval_[0] && 
sw < swInterval_[1])
 
  202            return krnSpline_->eval(
sw);
 
 
  212        if (
sw > swInterval_[0] && 
sw < swInterval_[1])
 
  213            return krnSpline_->evalDerivative(
sw);
 
 
  219    template<
class Function>
 
  220    std::unique_ptr<MonotoneCubicSpline<Scalar>>
 
  221    makeSweSpline_(
const Function& f, 
bool invert = 
false)
 const 
  223        const auto sw = 
linspace(swInterval_[0], swInterval_[1], numSwSamples_);
 
  226        std::transform(
sw.begin(), 
sw.end(), values.begin(), f);
 
  229            return std::make_unique<MonotoneCubicSpline<Scalar>>(values, 
sw);
 
  231            return std::make_unique<MonotoneCubicSpline<Scalar>>(
sw, values);
 
  234    std::unique_ptr<MonotoneCubicSpline<Scalar>> pcSpline_;
 
  235    std::unique_ptr<MonotoneCubicSpline<Scalar>> krwSpline_;
 
  236    std::unique_ptr<MonotoneCubicSpline<Scalar>> krnSpline_;
 
  238    std::array<Scalar, 2> swInterval_;
 
  239    std::array<Scalar, 2> swIntervalPc_;
 
  240    std::size_t numSwSamples_;
 
 
typename TwoPMaterialLaw::RegularizationParams RegularizationParams
Definition splinemateriallaw.hh:44
static constexpr bool isRegularized()
We are always regularized in the sense that we replace the original curve by a cubic spline.
Definition splinemateriallaw.hh:56
Scalar dkrn_dsw(const Scalar sw) const
The derivative of the relative permeability for the non-wetting phase w.r.t. saturation.
Definition splinemateriallaw.hh:210
typename TwoPMaterialLaw::BasicParams BasicParams
Definition splinemateriallaw.hh:42
Scalar pc(const Scalar sw) const
The capillary pressure-saturation curve.
Definition splinemateriallaw.hh:113
static constexpr int numFluidPhases()
Return the number of fluid phases.
Definition splinemateriallaw.hh:49
typename TwoPMaterialLaw::Scalar Scalar
Definition splinemateriallaw.hh:40
SplineTwoPMaterialLaw(const std::string ¶mGroup)
Construct from a subgroup from the global parameter tree.
Definition splinemateriallaw.hh:69
SplineTwoPMaterialLaw(const std::array< Scalar, 2 > &sweInterval, std::size_t numSwSamples, TwoPMaterialLaw &&twoP)
Construct from parameter structs.
Definition splinemateriallaw.hh:98
Scalar dsw_dpc(const Scalar pc) const
The partial derivative of the saturation to the capillary pressure.
Definition splinemateriallaw.hh:161
Scalar krn(const Scalar sw) const
The relative permeability for the non-wetting phase.
Definition splinemateriallaw.hh:199
Scalar krw(const Scalar sw) const
The relative permeability for the wetting phase.
Definition splinemateriallaw.hh:177
typename TwoPMaterialLaw::BasicParams EffToAbsParams
Definition splinemateriallaw.hh:43
Scalar dkrw_dsw(const Scalar sw) const
The derivative of the relative permeability for the wetting phase w.r.t. saturation.
Definition splinemateriallaw.hh:188
SplineTwoPMaterialLaw()=delete
Deleted default constructor (so we are never in an undefined state)
Scalar dpc_dsw(const Scalar sw) const
The partial derivative of the capillary pressure w.r.t. the saturation.
Definition splinemateriallaw.hh:129
Scalar sw(const Scalar pc) const
The saturation-capillary pressure curve.
Definition splinemateriallaw.hh:145
TwoPMaterialLaw()=delete
Deleted default constructor (so we are never in an undefined state)
Scalar dsw_dpc(const Scalar pc) const
The partial derivative of the saturation to the capillary pressure.
Definition materiallaw.hh:162
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
ScalarType Scalar
Definition materiallaw.hh:48
Scalar krn(const Scalar sw) const
The relative permeability for the non-wetting phase.
Definition materiallaw.hh:212
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 Regularization::template Params< Scalar > RegularizationParams
Definition materiallaw.hh:52
const EffToAbsParams & effToAbsParams() const
Return the parameters of the EffToAbs policy.
Definition materiallaw.hh:279
Scalar dpc_dsw(const Scalar sw) const
The partial derivative of the capillary pressure w.r.t. the saturation.
Definition materiallaw.hh:121
Scalar sw(const Scalar pc) const
The saturation-capillary pressure curve.
Definition materiallaw.hh:146
typename BaseLaw::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
Wrapper type to combine an arbitrary number of different laws for fluid-matrix interaction (e....
std::vector< Scalar > linspace(const Scalar begin, const Scalar end, std::size_t samples, bool endPoint=true)
Generates linearly spaced vectors.
Definition math.hh:611
T getParamFromGroup(Args &&... args)
A free function to get a parameter from the parameter tree singleton with a model group.
Definition parameters.hh:149
Implementation helper for capillary pressure and relative permeability <-> saturation relations for t...
Definition brookscorey.hh:23
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