26 #ifndef OPM_VAN_GENUCHTEN_HPP
27 #define OPM_VAN_GENUCHTEN_HPP
53 template <
class TraitsT,
class ParamsT = VanGenuchtenParams<TraitsT> >
64 typedef typename Traits::Scalar
Scalar;
68 static_assert(numPhases == 2,
69 "The van Genuchten capillary pressure law only "
70 "applies to the case of two fluid phases");
74 static const bool implementsTwoPhaseApi =
true;
78 static const bool implementsTwoPhaseSatApi =
true;
82 static const bool isSaturationDependent =
true;
86 static const bool isPressureDependent =
false;
90 static const bool isTemperatureDependent =
false;
94 static const bool isCompositionDependent =
false;
112 template <
class Container,
class Flu
idState>
113 static void capillaryPressures(Container &values,
const Params ¶ms,
const FluidState &fs)
115 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
117 values[Traits::wettingPhaseIdx] = 0.0;
118 values[Traits::nonWettingPhaseIdx] = pcnw<FluidState, Evaluation>(params, fs);
125 template <
class Container,
class Flu
idState>
126 static void saturations(Container &values,
const Params ¶ms,
const FluidState &fs)
128 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
130 values[Traits::wettingPhaseIdx] = Sw<FluidState, Evaluation>(params, fs);
131 values[Traits::nonWettingPhaseIdx] = 1 - values[Traits::wettingPhaseIdx];
144 template <
class Container,
class Flu
idState>
145 static void relativePermeabilities(Container &values,
const Params ¶ms,
const FluidState &fs)
147 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
149 values[Traits::wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fs);
150 values[Traits::nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fs);
167 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
168 static Evaluation pcnw(
const Params ¶ms,
const FluidState &fs)
170 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
172 const Evaluation& Sw =
173 FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
175 assert(0 <= Sw && Sw <= 1);
177 return twoPhaseSatPcnw(params, Sw);
194 template <
class Evaluation>
195 static Evaluation twoPhaseSatPcnw(
const Params ¶ms,
const Evaluation& Sw)
197 typedef MathToolbox<Evaluation> Toolbox;
214 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
215 static Evaluation Sw(
const Params ¶ms,
const FluidState &fs)
217 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
220 FsToolbox::template toLhs<Evaluation>(fs.pressure(Traits::nonWettingPhaseIdx))
221 - FsToolbox::template toLhs<Evaluation>(fs.pressure(Traits::wettingPhaseIdx));
222 return twoPhaseSatSw(params, pC);
225 template <
class Evaluation>
226 static Evaluation twoPhaseSatSw(
const Params ¶ms,
const Evaluation& pC)
228 typedef MathToolbox<Evaluation> Toolbox;
239 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
240 static Evaluation Sn(
const Params ¶ms,
const FluidState &fs)
241 {
return 1 - Sw<FluidState, Evaluation>(params, fs); }
243 template <
class Evaluation>
244 static Evaluation twoPhaseSatSn(
const Params ¶ms,
const Evaluation& pC)
245 {
return 1 - twoPhaseSatSw(params, pC); }
257 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
258 static Evaluation krw(
const Params ¶ms,
const FluidState &fs)
260 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
262 const Evaluation& Sw =
263 FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
265 return twoPhaseSatKrw(params, Sw);
268 template <
class Evaluation>
269 static Evaluation twoPhaseSatKrw(
const Params ¶ms,
const Evaluation& Sw)
271 typedef MathToolbox<Evaluation> Toolbox;
273 assert(0.0 <= Sw && Sw <= 1.0);
288 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
289 static Evaluation krn(
const Params ¶ms,
const FluidState &fs)
291 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
293 const Evaluation& Sw =
294 1.0 - FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
296 return twoPhaseSatKrn(params, Sw);
299 template <
class Evaluation>
300 static Evaluation twoPhaseSatKrn(
const Params ¶ms, Evaluation Sw)
302 typedef MathToolbox<Evaluation> Toolbox;
304 assert(0 <= Sw && Sw <= 1);
Definition: Air_Mesitylene.hpp:31
Implementation of the van Genuchten capillary pressure - saturation relation.
Definition: VanGenuchten.hpp:54
Evaluation< Scalar, VarSetTag, numVars > sqrt(const Evaluation< Scalar, VarSetTag, numVars > &x)
Definition: Math.hpp:278
ParamsT Params
The type of the parameter objects for this law.
Definition: VanGenuchten.hpp:61
static const int numPhases
The number of fluid phases.
Definition: VanGenuchten.hpp:67
Evaluation< Scalar, VarSetTag, numVars > pow(const Evaluation< Scalar, VarSetTag, numVars > &base, Scalar exp)
Definition: Math.hpp:312
Traits::Scalar Scalar
The type of the scalar values for this law.
Definition: VanGenuchten.hpp:64
Specification of the material parameters for the van Genuchten constitutive relations.
TraitsT Traits
The traits class for this material law.
Definition: VanGenuchten.hpp:58