27 #ifndef OPM_THREE_PHASE_PARKER_VAN_GENUCHTEN_HPP
28 #define OPM_THREE_PHASE_PARKER_VAN_GENUCHTEN_HPP
48 template <
class TraitsT,
49 class ParamsT = ThreePhaseParkerVanGenuchtenParams<TraitsT> >
53 static_assert(TraitsT::numPhases == 3,
54 "The number of phases considered by this capillary pressure "
55 "law is always three!");
57 typedef TraitsT Traits;
58 typedef ParamsT Params;
59 typedef typename Traits::Scalar Scalar;
61 static const int numPhases = 3;
62 static const int wettingPhaseIdx = Traits::wettingPhaseIdx;
63 static const int nonWettingPhaseIdx = Traits::nonWettingPhaseIdx;
64 static const int gasPhaseIdx = Traits::gasPhaseIdx;
68 static const bool implementsTwoPhaseApi =
false;
72 static const bool implementsTwoPhaseSatApi =
false;
76 static const bool isSaturationDependent =
true;
80 static const bool isPressureDependent =
false;
84 static const bool isTemperatureDependent =
false;
88 static const bool isCompositionDependent =
false;
101 template <
class ContainerT,
class Flu
idState>
102 static void capillaryPressures(ContainerT &values,
103 const Params ¶ms,
104 const FluidState &fluidState)
106 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
108 values[gasPhaseIdx] = pcgn<FluidState, Evaluation>(params, fluidState);
109 values[nonWettingPhaseIdx] = 0;
110 values[wettingPhaseIdx] = - pcnw<FluidState, Evaluation>(params, fluidState);
122 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
123 static Evaluation pcgn(
const Params ¶ms,
const FluidState &fluidState)
128 Scalar PC_VG_REG = 0.01;
132 FsToolbox::template toLhs<Evaluation>(fluidState.saturation(wettingPhaseIdx))
133 + FsToolbox::template toLhs<Evaluation>(fluidState.saturation(nonWettingPhaseIdx));
135 Evaluation Se = (St - params.Swrx())/(1. - params.Swrx());
143 if (Se>PC_VG_REG && Se<1-PC_VG_REG)
145 const Evaluation& x =
Toolbox::pow(Se,-1/params.vgM()) - 1;
146 return Toolbox::pow(x, 1.0 - params.vgM())/params.vgAlpha();
154 Se_regu = 1-PC_VG_REG;
155 const Evaluation& x =
std::pow(Se_regu,-1/params.vgM())-1;
156 const Evaluation& pc =
Toolbox::pow(x, 1.0/params.vgN())/params.vgAlpha();
157 const Evaluation& pc_prime =
159 *
std::pow(Se_regu,-1/params.vgM()-1)
162 / (1 - params.Sgr() - params.Swrx())
166 return ((Se-Se_regu)*pc_prime + pc)/params.betaGN();
178 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
179 static Evaluation pcnw(
const Params ¶ms,
const FluidState &fluidState)
184 const Evaluation& Sw =
185 FsToolbox::template toLhs<Evaluation>(fluidState.saturation(wettingPhaseIdx));
186 Evaluation Se = (Sw-params.Swr())/(1.-params.Snr());
188 Scalar PC_VG_REG = 0.01;
196 if (Se>PC_VG_REG && Se<1-PC_VG_REG) {
199 return x/params.vgAlpha();
207 Se_regu = 1.0 - PC_VG_REG;
209 const Evaluation& x =
std::pow(Se_regu,-1/params.vgM())-1;
210 const Evaluation& pc =
Toolbox::pow(x, 1/params.vgN())/params.vgAlpha();
211 const Evaluation& pc_prime =
213 *
std::pow(Se_regu, -1.0/params.vgM() - 1)
216 / (1-params.Snr()-params.Swr())
220 return ((Se-Se_regu)*pc_prime + pc)/params.betaNW();
227 template <
class ContainerT,
class Flu
idState>
228 static void saturations(ContainerT &,
231 { OPM_THROW(std::logic_error,
"Not implemented: inverse capillary pressures"); }
236 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
237 static Evaluation Sg(
const Params &,
const FluidState &)
238 { OPM_THROW(std::logic_error,
"Not implemented: Sg()"); }
243 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
244 static Evaluation Sn(
const Params &,
const FluidState &)
245 { OPM_THROW(std::logic_error,
"Not implemented: Sn()"); }
250 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
251 static Evaluation Sw(
const Params &,
const FluidState &)
252 { OPM_THROW(std::logic_error,
"Not implemented: Sw()"); }
257 template <
class ContainerT,
class Flu
idState>
258 static void relativePermeabilities(ContainerT &values,
259 const Params ¶ms,
260 const FluidState &fluidState)
262 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
264 values[wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fluidState);
265 values[nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fluidState);
266 values[gasPhaseIdx] = krg<FluidState, Evaluation>(params, fluidState);
277 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
278 static Evaluation krw(
const Params ¶ms,
const FluidState &fluidState)
283 const Evaluation& Sw =
284 FsToolbox::template toLhs<Evaluation>(fluidState.saturation(wettingPhaseIdx));
286 const Evaluation& Se = (Sw - params.Swr()) / (1-params.Swr());
289 if(Se > 1.0)
return 1.;
290 if(Se < 0.0)
return 0.;
308 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
309 static Evaluation krn(
const Params ¶ms,
const FluidState &fluidState)
314 const Evaluation& Sn =
315 FsToolbox::template toLhs<Evaluation>(fluidState.saturation(nonWettingPhaseIdx));
316 const Evaluation& Sw =
317 FsToolbox::template toLhs<Evaluation>(fluidState.saturation(wettingPhaseIdx));
318 Evaluation Swe =
Toolbox::min((Sw - params.Swr()) / (1 - params.Swr()), 1.);
319 Evaluation Ste =
Toolbox::min((Sw + Sn - params.Swr()) / (1 - params.Swr()), 1.);
322 if(Swe <= 0.0) Swe = 0.;
323 if(Ste <= 0.0) Ste = 0.;
324 if(Ste - Swe <= 0.0)
return 0.;
331 if (params.krRegardsSnr())
335 const Evaluation& resIncluded =
356 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
357 static Evaluation krg(
const Params ¶ms,
const FluidState &fluidState)
362 const Evaluation& Sg =
363 FsToolbox::template toLhs<Evaluation>(fluidState.saturation(gasPhaseIdx));
364 const Evaluation& Se =
Toolbox::min(((1-Sg) - params.Sgr()) / (1 - params.Sgr()), 1.);
372 Evaluation scaleFactor = 1.;
374 scaleFactor = (Sg - params.Sgr())/(0.1 - params.Sgr());
375 if (scaleFactor < 0.)
Definition: Air_Mesitylene.hpp:31
Evaluation< Scalar, VarSetTag, numVars > max(const Evaluation< Scalar, VarSetTag, numVars > &x1, const Evaluation< Scalar, VarSetTag, numVars > &x2)
Definition: Math.hpp:114
Evaluation< Scalar, VarSetTag, numVars > sqrt(const Evaluation< Scalar, VarSetTag, numVars > &x)
Definition: Math.hpp:278
Specification of the material params for the three-phase van Genuchten capillary pressure model...
Evaluation< Scalar, VarSetTag, numVars > min(const Evaluation< Scalar, VarSetTag, numVars > &x1, const Evaluation< Scalar, VarSetTag, numVars > &x2)
Definition: Math.hpp:61
Implementation of three-phase capillary pressure and relative permeability relations proposed by Park...
Definition: ThreePhaseParkerVanGenuchten.hpp:50
Evaluation< Scalar, VarSetTag, numVars > pow(const Evaluation< Scalar, VarSetTag, numVars > &base, Scalar exp)
Definition: Math.hpp:312