27 #ifndef OPM_REGULARIZED_VAN_GENUCHTEN_HPP
28 #define OPM_REGULARIZED_VAN_GENUCHTEN_HPP
70 template <
class TraitsT,
class ParamsT = RegularizedVanGenuchtenParams<TraitsT> >
79 typedef typename Traits::Scalar
Scalar;
83 static_assert(numPhases == 2,
84 "The regularized van Genuchten capillary pressure law only "
85 "applies to the case of two fluid phases");
89 static const bool implementsTwoPhaseApi =
true;
93 static const bool implementsTwoPhaseSatApi =
true;
97 static const bool isSaturationDependent =
true;
101 static const bool isPressureDependent =
false;
105 static const bool isTemperatureDependent =
false;
109 static const bool isCompositionDependent =
false;
115 template <
class Container,
class Flu
idState>
116 static void capillaryPressures(Container &values,
const Params ¶ms,
const FluidState &fs)
118 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
120 values[Traits::wettingPhaseIdx] = 0.0;
121 values[Traits::nonWettingPhaseIdx] = pcnw<FluidState, Evaluation>(params, fs);
128 template <
class Container,
class Flu
idState>
129 static void saturations(Container &values,
const Params ¶ms,
const FluidState &fs)
131 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
133 values[Traits::wettingPhaseIdx] = Sw<FluidState, Evaluation>(params, fs);
134 values[Traits::nonWettingPhaseIdx] = 1 - values[Traits::wettingPhaseIdx];
141 template <
class Container,
class Flu
idState>
142 static void relativePermeabilities(Container &values,
const Params ¶ms,
const FluidState &fs)
144 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
146 values[Traits::wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fs);
147 values[Traits::nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fs);
162 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
163 static Evaluation pcnw(
const Params ¶ms,
const FluidState &fs)
165 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
167 const auto& Sw = FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
168 return twoPhaseSatPcnw(params, Sw);
171 template <
class Evaluation>
172 static Evaluation twoPhaseSatPcnw(
const Params ¶ms,
const Evaluation& Sw)
176 const Scalar SwThLow = params.pcnwLowSw();
177 const Scalar SwThHigh = params.pcnwHighSw();
185 return params.pcnwLow() + params.pcnwSlopeLow()*(Sw - SwThLow);
187 else if (Sw > SwThHigh)
189 Scalar yTh = params.pcnwHigh();
190 Scalar m1 = (0.0 - yTh)/(1.0 - SwThHigh)*2;
194 const Spline<Scalar> &sp = params.pcnwHighSpline();
200 return m1*(Sw - 1.0) + 0.0;
206 return VanGenuchten::twoPhaseSatPcnw(params, Sw);
223 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
224 static Evaluation Sw(
const Params ¶ms,
const FluidState &fs)
226 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
228 const Evaluation& pC =
229 FsToolbox::template toLhs<Evaluation>(fs.pressure(Traits::nonWettingPhaseIdx))
230 - FsToolbox::template toLhs<Evaluation>(fs.pressure(Traits::wettingPhaseIdx));
231 return twoPhaseSatSw(params, pC);
234 template <
class Evaluation>
235 static Evaluation twoPhaseSatSw(
const Params ¶ms,
const Evaluation& pC)
239 const Scalar SwThLow = params.pcnwLowSw();
240 const Scalar SwThHigh = params.pcnwHighSw();
247 Scalar m1 = params.pcnwSlopeHigh();
251 Evaluation Sw = VanGenuchten::twoPhaseSatSw(params, pC);
256 Scalar pC_SwLow = VanGenuchten::twoPhaseSatPcnw(params, SwThLow);
257 return (pC - pC_SwLow)/params.pcnwSlopeLow() + SwThLow;
259 else if (SwThHigh < Sw )
262 const Spline<Scalar>& spline = params.pcnwHighSpline();
264 return spline.template intersectInterval<Evaluation>(SwThHigh, 1.0,
275 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
276 static Evaluation Sn(
const Params ¶ms,
const FluidState &fs)
277 {
return 1 - Sw<FluidState, Evaluation>(params, fs); }
279 template <
class Evaluation>
280 static Evaluation twoPhaseSatSn(
const Params ¶ms,
const Evaluation& pc)
281 {
return 1 - twoPhaseSatSw(params, pc); }
297 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
298 static Evaluation krw(
const Params ¶ms,
const FluidState &fs)
300 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
302 const auto& Sw = FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
303 return twoPhaseSatKrw(params, Sw);
306 template <
class Evaluation>
307 static Evaluation twoPhaseSatKrw(
const Params ¶ms,
const Evaluation& Sw)
309 typedef MathToolbox<Evaluation> Toolbox;
313 return Toolbox::createConstant(0);
315 return Toolbox::createConstant(1);
317 return VanGenuchten::twoPhaseSatKrw(params, Sw);
334 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
335 static Evaluation krn(
const Params ¶ms,
const FluidState &fs)
337 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
339 const Evaluation& Sw =
340 1.0 - FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
341 return twoPhaseSatKrn(params, Sw);
344 template <
class Evaluation>
345 static Evaluation twoPhaseSatKrn(
const Params ¶ms,
const Evaluation& Sw)
347 typedef MathToolbox<Evaluation> Toolbox;
351 return Toolbox::createConstant(1);
353 return Toolbox::createConstant(0);
355 return VanGenuchten::twoPhaseSatKrn(params, Sw);
Definition: Air_Mesitylene.hpp:31
Implementation of the van Genuchten capillary pressure - saturation relation.
Definition: VanGenuchten.hpp:54
Implementation of the van Genuchten capillary pressure - saturation relation.
Class implementing cubic splines.
Traits::Scalar Scalar
Definition: RegularizedVanGenuchten.hpp:79
ParamsT Params
Definition: RegularizedVanGenuchten.hpp:77
TraitsT Traits
Definition: RegularizedVanGenuchten.hpp:76
static const int numPhases
The number of fluid phases.
Definition: RegularizedVanGenuchten.hpp:82
Parameters that are necessary for the regularization of VanGenuchten "material law".
Implementation of the regularized van Genuchten's capillary pressure / relative permeability <-> satu...
Definition: RegularizedVanGenuchten.hpp:71