3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
exponentialcubic.hh
Go to the documentation of this file.
1/*****************************************************************************
2 * See the file COPYING for full copying permissions. *
3 * *
4 * This program is free software: you can redistribute it and/or modify *
5 * it under the terms of the GNU General Public License as published by *
6 * the Free Software Foundation, either version 3 of the License, or *
7 * (at your option) any later version. *
8 * *
9 * This program is distributed in the hope that it will be useful, *
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
12 * GNU General Public License for more details. *
13 * *
14 * You should have received a copy of the GNU General Public License *
15 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
16 *****************************************************************************/
26#ifndef DUMUX_MATERIAL_FLUIDMATRIX_TWO_P_INTERFACIAL_AREA_EXPONENTIAL_CUBIC
27#define DUMUX_MATERIAL_FLUIDMATRIX_TWO_P_INTERFACIAL_AREA_EXPONENTIAL_CUBIC
28
29#include <cmath>
30#include <dune/common/exceptions.hh>
31#include <dune/common/float_cmp.hh>
33
34namespace Dumux::FluidMatrix {
35
42{
43public:
44
45 template<class Scalar>
46 struct Params
47 {
48 Params(Scalar a1 = 0, Scalar a2 = 0, Scalar a3 = 0)
49 : a1_(a1), a2_(a2), a3_(a3) {}
50
51 Scalar a1() const { return a1_; }
52 void setA1(Scalar a1) { a1_ = a1; }
53
54 Scalar a2() const { return a2_; }
55 void setA2(Scalar a2) { a2_ = a2; }
56
57 Scalar a3() const { return a3_; }
58 void setA3(Scalar a3) { a3_ = a3; }
59
60 bool operator== (const Params& p) const
61 {
62 return Dune::FloatCmp::eq(a1(), p.a1(), 1e-6)
63 && Dune::FloatCmp::eq(a2(), p.a2(), 1e-6)
64 && Dune::FloatCmp::eq(a3(), p.a3(), 1e-6);
65 }
66
67 private:
68 Scalar a1_, a2_, a3_;
69 };
70
75 template<class Scalar = double>
76 static Params<Scalar> makeParams(const std::string& paramGroup)
77 {
78 const auto a1 = getParamFromGroup<Scalar>(paramGroup, "A1");
79 const auto a2 = getParamFromGroup<Scalar>(paramGroup, "A2");
80 const auto a3 = getParamFromGroup<Scalar>(paramGroup, "A3");
81 return {a1, a2, a3};
82 }
83
95 template<class Scalar>
96 static Scalar area(const Scalar swe, const Scalar pc, const Params<Scalar>& params)
97 {
98 // TODO think about awn surface for relative saturation
99 const Scalar a1 = params.a1();
100 const Scalar a2 = params.a2();
101 const Scalar a3 = params.a3();
102
103 using std::exp;
104 return a1 * exp( a2 * swe) + a3 * pc * pc * pc ;
105 }
106};
107
108} // namespace Dumux::FluidMatrix
109
110#endif
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
Definition: brookscorey.hh:35
Implementation of a exponential function relating specific interfacial area to wetting phase saturati...
Definition: exponentialcubic.hh:42
static Params< Scalar > makeParams(const std::string &paramGroup)
Construct from a subgroup from the global parameter tree.
Definition: exponentialcubic.hh:76
static Scalar area(const Scalar swe, const Scalar pc, const Params< Scalar > &params)
The interfacial area.
Definition: exponentialcubic.hh:96
void setA2(Scalar a2)
Definition: exponentialcubic.hh:55
Params(Scalar a1=0, Scalar a2=0, Scalar a3=0)
Definition: exponentialcubic.hh:48
Scalar a1() const
Definition: exponentialcubic.hh:51
void setA1(Scalar a1)
Definition: exponentialcubic.hh:52
bool operator==(const Params &p) const
Definition: exponentialcubic.hh:60
Scalar a2() const
Definition: exponentialcubic.hh:54
Scalar a3() const
Definition: exponentialcubic.hh:57
void setA3(Scalar a3)
Definition: exponentialcubic.hh:58