13#ifndef DUMUX_MATERIAL_FLUIDMATRIX_TWO_P_INTERFACIAL_AREA_POLYNOMIAL_SECOND_ORDER_HH 
   14#define DUMUX_MATERIAL_FLUIDMATRIX_TWO_P_INTERFACIAL_AREA_POLYNOMIAL_SECOND_ORDER_HH 
   17#include <dune/common/exceptions.hh> 
   18#include <dune/common/float_cmp.hh> 
   32    template<
class Scalar>
 
   38        Scalar 
a00()
 const { 
return a00_; }
 
   41        Scalar 
a01()
 const { 
return a01_; }
 
   44        Scalar 
a02()
 const { 
return a02_; }
 
   47        Scalar 
a10()
 const { 
return a10_; }
 
   50        Scalar 
a11()
 const { 
return a11_; }
 
   53        Scalar 
a20()
 const { 
return a20_; }
 
   58            return Dune::FloatCmp::eq(
a00(), p.
a00(), 1e-6)
 
   59                   && Dune::FloatCmp::eq(
a01(), p.
a01(), 1e-6)
 
   60                   && Dune::FloatCmp::eq(
a02(), p.
a02(), 1e-6)
 
   61                   && Dune::FloatCmp::eq(
a10(), p.
a10(), 1e-6)
 
   62                   && Dune::FloatCmp::eq(
a11(), p.
a11(), 1e-6)
 
   63                   && Dune::FloatCmp::eq(
a20(), p.
a20(), 1e-6);
 
 
   67        Scalar a00_, a01_, a02_, a10_, a11_, a20_;
 
 
   74    template<
class Scalar = 
double>
 
   83        return {a00, a01, a02, a10, a11, a20};
 
 
   97    template<
class Scalar>
 
  100        const Scalar a00 = params.
a00();
 
  101        const Scalar a10 = params.
a10();
 
  102        const Scalar a20 = params.
a20();
 
  103        const Scalar a11 = params.
a11();
 
  104        const Scalar a01 = params.
a01();
 
  105        const Scalar a02 = params.
a02();
 
  107        return a00 + a10 * swe + a20 * swe*swe + a11*swe*pc +  a01*pc + a02*pc*pc;
 
 
  117    template<
class Scalar>
 
  120        return params.
a11()*swe + params.
a01() + 2.0*params.
a02() * pc;
 
 
  130    template<
class Scalar>
 
  133        return params.
a11()*pc + params.
a10() + 2.0*params.
a20()*swe;
 
 
 
Implementation of the polynomial of second order relating specific interfacial area to wetting phase ...
Definition polynomial2ndorder.hh:29
static Scalar darea_dpc(const Scalar swe, const Scalar pc, const Params< Scalar > ¶ms)
the derivative of specific interfacial area function w.r.t. capillary pressure
Definition polynomial2ndorder.hh:118
static Params< Scalar > makeParams(const std::string ¶mGroup)
Construct from a subgroup from the global parameter tree.
Definition polynomial2ndorder.hh:75
static Scalar darea_dsw(const Scalar swe, const Scalar pc, const Params< Scalar > ¶ms)
the derivative of specific interfacial area function w.r.t. saturation
Definition polynomial2ndorder.hh:131
static Scalar area(const Scalar swe, const Scalar pc, const Params< Scalar > ¶ms)
The interfacial area.
Definition polynomial2ndorder.hh:98
T getParamFromGroup(Args &&... args)
A free function to get a parameter from the parameter tree singleton with a model group.
Definition parameters.hh:149
Definition brookscorey.hh:23
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
Definition polynomial2ndorder.hh:34
Params(Scalar a00=0, Scalar a01=0, Scalar a02=0, Scalar a10=0, Scalar a11=0, Scalar a20=0)
Definition polynomial2ndorder.hh:35
void setA20(Scalar a20)
Definition polynomial2ndorder.hh:54
void setA02(Scalar a02)
Definition polynomial2ndorder.hh:45
bool operator==(const Params &p) const
Definition polynomial2ndorder.hh:56
Scalar a20() const
Definition polynomial2ndorder.hh:53
Scalar a00() const
Definition polynomial2ndorder.hh:38
Scalar a11() const
Definition polynomial2ndorder.hh:50
void setA01(Scalar a01)
Definition polynomial2ndorder.hh:42
Scalar a10() const
Definition polynomial2ndorder.hh:47
Scalar a02() const
Definition polynomial2ndorder.hh:44
void setA00(Scalar a00)
Definition polynomial2ndorder.hh:39
void setA11(Scalar a11)
Definition polynomial2ndorder.hh:51
Scalar a01() const
Definition polynomial2ndorder.hh:41
void setA10(Scalar a10)
Definition polynomial2ndorder.hh:48