13#ifndef DUMUX_MATERIAL_FLUIDMATRIX_TWO_P_INTERFACIAL_AREA_EXPONENTIAL_CUBIC 
   14#define DUMUX_MATERIAL_FLUIDMATRIX_TWO_P_INTERFACIAL_AREA_EXPONENTIAL_CUBIC 
   17#include <dune/common/exceptions.hh> 
   18#include <dune/common/float_cmp.hh> 
   32    template<
class Scalar>
 
   38        Scalar 
swr()
 const { 
return swr_; }
 
   41        Scalar 
a1()
 const { 
return a1_; }
 
   44        Scalar 
a2()
 const { 
return a2_; }
 
   47        Scalar 
a3()
 const { 
return a3_; }
 
   52            return Dune::FloatCmp::eq(
swr(), p.
swr(), 1e-6)
 
   53                   && Dune::FloatCmp::eq(
a1(), p.
a1(), 1e-6)
 
   54                   && Dune::FloatCmp::eq(
a2(), p.
a2(), 1e-6)
 
   55                   && Dune::FloatCmp::eq(
a3(), p.
a3(), 1e-6);
 
 
   59        Scalar swr_, a1_, a2_, a3_;
 
 
   66    template<
class Scalar = 
double>
 
   73        return {swr, a1, a2, a3};
 
 
   86    template<
class Scalar>
 
   89        const Scalar a1 = params.
a1();
 
   90        const Scalar a2 = params.
a2();
 
   91        const Scalar a3 = params.
a3();
 
   92        const Scalar swr = params.
swr();
 
   93        const Scalar factor = (swr-sw)*(1.0-sw) ;
 
   94        return a1*factor + a2*factor*pc + a3*factor*pc*pc;
 
 
  104    template<
class Scalar>
 
  107        const Scalar swr = params.
swr();
 
  108        const Scalar a1 = params.
a1();
 
  109        const Scalar a2 = params.
a2();
 
  110        const Scalar a3 = params.
a3();
 
  111        return a2*(swr-sw)*(1-sw) + a3*(swr-sw)*(1-sw)*pc;
 
 
  121    template<
class Scalar>
 
  124        const Scalar swr = params.
swr();
 
  125        const Scalar a1 = params.
a1();
 
  126        const Scalar a2 = params.
a2();
 
  127        const Scalar a3 = params.
a3();
 
  128        const Scalar derivativeFactor = (sw-1.0)+(sw-swr);
 
  129        return a1 * derivativeFactor + a2 * derivativeFactor * pc + a3 * derivativeFactor * pc*pc ;
 
 
 
Implementation of the polynomial of second order relating specific interfacial area to wetting phase ...
Definition polynomialedgezero2ndorder.hh:29
static Params< Scalar > makeParams(const std::string ¶mGroup)
Construct from a subgroup from the global parameter tree.
Definition polynomialedgezero2ndorder.hh:67
static Scalar darea_dpc(const Scalar sw, const Scalar pc, const Params< Scalar > ¶ms)
the derivative of specific interfacial area function w.r.t. capillary pressure
Definition polynomialedgezero2ndorder.hh:105
static Scalar darea_dsw(const Scalar sw, const Scalar pc, const Params< Scalar > ¶ms)
the derivative of specific interfacial area function w.r.t. saturation
Definition polynomialedgezero2ndorder.hh:122
static Scalar area(const Scalar sw, const Scalar pc, const Params< Scalar > ¶ms)
The interfacial area the suggested (as estimated from pore network models) awn surface:   ].
Definition polynomialedgezero2ndorder.hh:87
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 polynomialedgezero2ndorder.hh:34
void setA2(Scalar a2)
Definition polynomialedgezero2ndorder.hh:45
Scalar a1() const
Definition polynomialedgezero2ndorder.hh:41
void setSwr(Scalar swr)
Definition polynomialedgezero2ndorder.hh:39
Scalar a3() const
Definition polynomialedgezero2ndorder.hh:47
Params(Scalar swr=0, Scalar a1=0, Scalar a2=0, Scalar a3=0)
Definition polynomialedgezero2ndorder.hh:35
bool operator==(const Params &p) const
Definition polynomialedgezero2ndorder.hh:50
void setA3(Scalar a3)
Definition polynomialedgezero2ndorder.hh:48
void setA1(Scalar a1)
Definition polynomialedgezero2ndorder.hh:42
Scalar swr() const
Definition polynomialedgezero2ndorder.hh:38
Scalar a2() const
Definition polynomialedgezero2ndorder.hh:44