3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
interfacialarea.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_MATERIAL_FLUIDMATRIX_TWO_P_INTERFACIAL_AREA_HH
25#define DUMUX_MATERIAL_FLUIDMATRIX_TWO_P_INTERFACIAL_AREA_HH
26
31
32namespace Dumux::FluidMatrix {
33
50template<class ScalarType,
51 class BaseLaw,
52 template<class> class InterfaceType,
53 class Regularization = NoRegularization,
54 class EffToAbsPolicy = TwoPEffToAbsDefaultPolicy>
55class InterfacialArea : public Adapter<InterfacialArea<ScalarType, BaseLaw, InterfaceType, Regularization, EffToAbsPolicy>, InterfaceType>
56{
57public:
58
59 using Scalar = ScalarType;
60
61 using BasicParams = typename BaseLaw::template Params<Scalar>;
62 using EffToAbsParams = typename EffToAbsPolicy::template Params<Scalar>;
63 using RegularizationParams = typename Regularization::template Params<Scalar>;
64
65 using EffToAbs = EffToAbsPolicy;
66
70 static constexpr bool isRegularized()
71 { return !std::is_same<Regularization, NoRegularization>::value; }
72
77 InterfacialArea() = delete;
78
83 explicit InterfacialArea(const std::string& paramGroup)
84 : basicParams_(makeBasicParams(paramGroup))
85 , effToAbsParams_(makeEffToAbsParams(paramGroup))
86 {
87 if constexpr (isRegularized())
88 regularization_.init(this, paramGroup);
89 }
90
95 InterfacialArea(const BasicParams& baseParams,
97 const RegularizationParams& regParams = {})
98 : basicParams_(baseParams)
99 , effToAbsParams_(effToAbsParams)
100 {
101 if constexpr (isRegularized())
102 regularization_.init(this, baseParams, effToAbsParams, regParams);
103 }
104
108 template<bool enableRegularization = isRegularized()>
109 Scalar area(const Scalar sw, const Scalar pc) const
110 {
111 const auto swe = EffToAbs::swToSwe(sw, effToAbsParams_);
112 if constexpr (enableRegularization)
113 {
114 const auto regularized = regularization_.area(swe, pc);
115 if (regularized)
116 return regularized.value();
117 }
118
119 return BaseLaw::area(swe, pc, basicParams_);
120 }
121
125 template<bool enableRegularization = isRegularized()>
126 Scalar darea_dpc(const Scalar sw, const Scalar pc) const
127 {
128 if constexpr (enableRegularization)
129 {
130 const auto regularized = regularization_.darea_dpc(sw, pc);
131 if (regularized)
132 return regularized.value();
133 }
134
135 return BaseLaw::darea_dpc(sw, pc, basicParams_);
136 }
137
141 template<bool enableRegularization = isRegularized()>
142 Scalar darea_dsw(const Scalar sw, const Scalar pc) const
143 {
144 if constexpr (enableRegularization)
145 {
146 const auto regularized = regularization_.darea_dsw(sw, pc);
147 if (regularized)
148 return regularized.value();
149 }
150
151 return BaseLaw::darea_dsw(sw, pc, basicParams_);
152 }
153
157 bool operator== (const InterfacialArea& o) const
158 {
159 return basicParams_ == o.basicParams_
160 && effToAbsParams_ == o.effToAbsParams_
161 && regularization_ == o.regularization_;
162 }
163
168 static BasicParams makeBasicParams(const std::string& paramGroup)
169 {
170 return BaseLaw::template makeParams<Scalar>(paramGroup);
171 }
172
177 { return basicParams_; }
178
183 static EffToAbsParams makeEffToAbsParams(const std::string& paramGroup)
184 {
185 return EffToAbs::template makeParams<Scalar>(paramGroup);
186 }
187
192 { return effToAbsParams_; }
193
194private:
195 BasicParams basicParams_;
196 EffToAbsParams effToAbsParams_;
197 Regularization regularization_;
198};
199
200} // end namespace Dumux::FluidMatrix
201
202#endif
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
This is a policy for 2p material laws how to convert absolute to relative saturations and vice versa.
A tag to turn off regularization and it's overhead.
Wrapper type to combine an arbitrary number of different laws for fluid-matrix interaction (e....
Definition: brookscorey.hh:35
Wrapper class to implement regularized laws (pc-sw-a) with a conversion policy between absolution and...
Definition: interfacialarea.hh:56
EffToAbsPolicy EffToAbs
Definition: interfacialarea.hh:65
const BasicParams & basicParams() const
Return the base law's parameters.
Definition: interfacialarea.hh:176
bool operator==(const InterfacialArea &o) const
Equality comparison with another instance.
Definition: interfacialarea.hh:157
ScalarType Scalar
Definition: interfacialarea.hh:59
Scalar darea_dpc(const Scalar sw, const Scalar pc) const
The partial derivative of the capillary pressure w.r.t. the saturation.
Definition: interfacialarea.hh:126
const EffToAbsParams & effToAbsParams() const
Return the parameters of the EffToAbs policy.
Definition: interfacialarea.hh:191
typename EffToAbsPolicy::template Params< Scalar > EffToAbsParams
Definition: interfacialarea.hh:62
InterfacialArea()=delete
Deleted default constructor (so we are never in an undefined state)
Scalar area(const Scalar sw, const Scalar pc) const
The capillary pressure-saturation curve.
Definition: interfacialarea.hh:109
Scalar darea_dsw(const Scalar sw, const Scalar pc) const
The partial derivative of the saturation to the capillary pressure.
Definition: interfacialarea.hh:142
static EffToAbsParams makeEffToAbsParams(const std::string &paramGroup)
Create the parameters of the EffToAbs policy using input file parameters.
Definition: interfacialarea.hh:183
InterfacialArea(const BasicParams &baseParams, const EffToAbsParams &effToAbsParams={}, const RegularizationParams &regParams={})
Construct from parameter structs.
Definition: interfacialarea.hh:95
static constexpr bool isRegularized()
Return whether this law is regularized.
Definition: interfacialarea.hh:70
static BasicParams makeBasicParams(const std::string &paramGroup)
Create the base law's parameters using input file parameters.
Definition: interfacialarea.hh:168
InterfacialArea(const std::string &paramGroup)
Construct from a subgroup from the global parameter tree.
Definition: interfacialarea.hh:83
typename Regularization::template Params< Scalar > RegularizationParams
Definition: interfacialarea.hh:63
typename BaseLaw::template Params< Scalar > BasicParams
Definition: interfacialarea.hh:61
Adapter to inherit from, allowing the inheriting class to be wrapped by the makeFluidMatrixInteractio...
Definition: fluidmatrixinteraction.hh:67