12#ifndef DUMUX_MATERIAL_FLUIDMATRIX_DATA_SPLINE_MATERIAL_HH 
   13#define DUMUX_MATERIAL_FLUIDMATRIX_DATA_SPLINE_MATERIAL_HH 
   33template <
class S, 
bool approximatePcSwInverse = false>
 
   35: 
public Adapter<DataSplineTwoPMaterialLaw<S, approximatePcSwInverse>, PcKrSw>
 
   61        using V = std::vector<Scalar>;
 
   62        const std::string swPcGroup = paramGroup.empty() ? 
"Pc" : paramGroup + 
".Pc";
 
   64        const std::string krwPcGroup = paramGroup.empty() ? 
"Krw" : paramGroup + 
".Krw";
 
   66        const std::string krnPcGroup = paramGroup.empty() ? 
"Krn" : paramGroup + 
".Krn";
 
   73        updateData_(swPc, 
pc, swKrw, 
krw, swKrn, 
krn);
 
 
   80                              const std::vector<Scalar>& 
pc,
 
   81                              const std::vector<Scalar>& swKrw,
 
   82                              const std::vector<Scalar>& 
krw,
 
   83                              const std::vector<Scalar>& swKrn,
 
   84                              const std::vector<Scalar>& 
krn)
 
   86        updateData_(swPc, 
pc, swKrw, 
krw, swKrn, 
krn);
 
 
   94        if constexpr (approximatePcSwInverse)
 
   95            return pcSpline_.evalInverse(
sw);
 
   97            return pcSpline_.eval(
sw);
 
 
  105        if constexpr (approximatePcSwInverse)
 
  106            return 1.0/pcSpline_.evalDerivative(pcSpline_.evalInverse(
sw));
 
  108            return pcSpline_.evalDerivative(
sw);
 
 
  116        if constexpr (approximatePcSwInverse)
 
  117            return pcSpline_.eval(
pc);
 
  119            return pcSpline_.evalInverse(
pc);
 
 
  127        if constexpr (approximatePcSwInverse)
 
  128            return pcSpline_.evalDerivative(
pc);
 
  130            return 1.0/pcSpline_.evalDerivative(pcSpline_.evalInverse(
pc));
 
 
  138        return krwSpline_.eval(
sw);
 
 
  146        return krwSpline_.evalDerivative(
sw);
 
 
  154        return krnSpline_.eval(
sw);
 
 
  162        return krnSpline_.evalDerivative(
sw);
 
 
  166    void updateData_(
const std::vector<Scalar>& swPc,
 
  167                     const std::vector<Scalar>& 
pc,
 
  168                     const std::vector<Scalar>& swKrw,
 
  169                     const std::vector<Scalar>& 
krw,
 
  170                     const std::vector<Scalar>& swKrn,
 
  171                     const std::vector<Scalar>& 
krn)
 
  173        if constexpr (approximatePcSwInverse)
 
 
Scalar pc(const Scalar sw) const
The capillary pressure-saturation curve.
Definition datasplinemateriallaw.hh:92
DataSplineTwoPMaterialLaw()=delete
Deleted default constructor (so we are never in an undefined state)
static constexpr int numFluidPhases()
Return the number of fluid phases.
Definition datasplinemateriallaw.hh:44
Scalar dkrw_dsw(const Scalar sw) const
The derivative of the relative permeability for the wetting phase w.r.t. saturation.
Definition datasplinemateriallaw.hh:144
Scalar dkrn_dsw(const Scalar sw) const
The derivative of the relative permeability for the non-wetting phase w.r.t. saturation.
Definition datasplinemateriallaw.hh:160
S Scalar
Definition datasplinemateriallaw.hh:39
DataSplineTwoPMaterialLaw(const std::vector< Scalar > &swPc, const std::vector< Scalar > &pc, const std::vector< Scalar > &swKrw, const std::vector< Scalar > &krw, const std::vector< Scalar > &swKrn, const std::vector< Scalar > &krn)
Construct from user data.
Definition datasplinemateriallaw.hh:79
static constexpr bool isRegularized()
Definition datasplinemateriallaw.hh:47
Scalar dpc_dsw(const Scalar sw) const
The partial derivative of the capillary pressure w.r.t. the saturation.
Definition datasplinemateriallaw.hh:103
Scalar krw(const Scalar sw) const
The relative permeability for the wetting phase.
Definition datasplinemateriallaw.hh:136
DataSplineTwoPMaterialLaw(const std::string ¶mGroup)
Construct from a subgroup from the global parameter tree.
Definition datasplinemateriallaw.hh:59
Scalar sw(const Scalar pc) const
The saturation-capillary pressure curve.
Definition datasplinemateriallaw.hh:114
Scalar krn(const Scalar sw) const
The relative permeability for the non-wetting phase.
Definition datasplinemateriallaw.hh:152
Scalar dsw_dpc(const Scalar pc) const
The partial derivative of the saturation to the capillary pressure.
Definition datasplinemateriallaw.hh:125
A monotone cubic spline.
Definition monotonecubicspline.hh:39
void updatePoints(const std::vector< Scalar > &x, const std::vector< Scalar > &y)
Create a monotone cubic spline from the control points (x[i], y[i])
Definition monotonecubicspline.hh:63
Wrapper type to combine an arbitrary number of different laws for fluid-matrix interaction (e....
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