27#ifndef OPM_THREE_PHASE_PARKER_VAN_GENUCHTEN_HPP
28#define OPM_THREE_PHASE_PARKER_VAN_GENUCHTEN_HPP
48template <
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!");
59 typedef typename Traits::Scalar
Scalar;
101 template <
class ContainerT,
class Flu
idState>
104 const FluidState& fluidState)
106 typedef typename std::remove_reference<
decltype(values[0])>::type Evaluation;
108 values[
gasPhaseIdx] = pcgn<FluidState, Evaluation>(params, fluidState);
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& params,
const FluidState& fluidState)
132 Evaluation Se = (St - params.Swrx())/(1. - params.Swrx());
140 if (Se>PC_VG_REG && Se<1-PC_VG_REG)
142 const Evaluation& x =
pow(Se,-1/params.vgM()) - 1;
143 return pow(x, 1.0 - params.vgM())/params.vgAlpha();
151 Se_regu = 1-PC_VG_REG;
152 const Evaluation& x =
std::pow(Se_regu,-1/params.vgM())-1;
153 const Evaluation& pc =
pow(x, 1.0/params.vgN())/params.vgAlpha();
154 const Evaluation& pc_prime =
155 pow(x, 1/params.vgN()-1)
156 *
std::pow(Se_regu,-1/params.vgM()-1)
159 / (1 - params.Sgr() - params.Swrx())
163 return ((Se-Se_regu)*pc_prime + pc)/params.betaGN();
175 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
176 static Evaluation
pcnw(
const Params& params,
const FluidState& fluidState)
178 const Evaluation&
Sw =
180 Evaluation Se = (
Sw-params.Swr())/(1.-params.Snr());
190 if (Se>PC_VG_REG && Se<1-PC_VG_REG) {
191 Evaluation x =
pow(Se,-1/params.vgM()) - 1.0;
192 x =
pow(x, 1 - params.vgM());
193 return x/params.vgAlpha();
201 Se_regu = 1.0 - PC_VG_REG;
203 const Evaluation& x =
std::pow(Se_regu,-1/params.vgM())-1;
204 const Evaluation& pc =
pow(x, 1/params.vgN())/params.vgAlpha();
205 const Evaluation& pc_prime =
206 pow(x,1/params.vgN()-1)
207 *
std::pow(Se_regu, -1.0/params.vgM() - 1)
210 / (1-params.Snr()-params.Swr())
214 return ((Se-Se_regu)*pc_prime + pc)/params.betaNW();
221 template <
class ContainerT,
class Flu
idState>
225 {
throw std::logic_error(
"Not implemented: inverse capillary pressures"); }
230 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
231 static Evaluation
Sg(
const Params& ,
const FluidState& )
232 {
throw std::logic_error(
"Not implemented: Sg()"); }
237 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
238 static Evaluation
Sn(
const Params& ,
const FluidState& )
239 {
throw std::logic_error(
"Not implemented: Sn()"); }
244 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
245 static Evaluation
Sw(
const Params& ,
const FluidState& )
246 {
throw std::logic_error(
"Not implemented: Sw()"); }
251 template <
class ContainerT,
class Flu
idState>
254 const FluidState& fluidState)
256 typedef typename std::remove_reference<
decltype(values[0])>::type Evaluation;
258 values[
wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fluidState);
260 values[
gasPhaseIdx] = krg<FluidState, Evaluation>(params, fluidState);
271 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
272 static Evaluation
krw(
const Params& params,
const FluidState& fluidState)
274 const Evaluation&
Sw =
277 const Evaluation& Se = (
Sw - params.Swr()) / (1-params.Swr());
280 if(Se > 1.0)
return 1.;
281 if(Se < 0.0)
return 0.;
283 const Evaluation& r = 1. -
pow(1 -
pow(Se, 1/params.vgM()), params.vgM());
299 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
300 static Evaluation
krn(
const Params& params,
const FluidState& fluidState)
302 const Evaluation&
Sn =
304 const Evaluation&
Sw =
306 Evaluation Swe =
min((
Sw - params.Swr()) / (1 - params.Swr()), 1.);
307 Evaluation Ste =
min((
Sw +
Sn - params.Swr()) / (1 - params.Swr()), 1.);
310 if(Swe <= 0.0) Swe = 0.;
311 if(Ste <= 0.0) Ste = 0.;
312 if(Ste - Swe <= 0.0)
return 0.;
315 krn_ =
pow(1 -
pow(Swe, 1/params.vgM()), params.vgM());
316 krn_ -=
pow(1 -
pow(Ste, 1/params.vgM()), params.vgM());
319 if (params.krRegardsSnr())
323 const Evaluation& resIncluded =
324 max(
min(
Sw - params.Snr() / (1-params.Swr()), 1.0), 0.0);
325 krn_ *=
sqrt(resIncluded );
328 krn_ *=
sqrt(
Sn / (1 - params.Swr()));
344 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
345 static Evaluation
krg(
const Params& params,
const FluidState& fluidState)
347 const Evaluation&
Sg =
348 decay<Evaluation>(fluidState.saturation(
gasPhaseIdx));
349 const Evaluation& Se =
min(((1-
Sg) - params.Sgr()) / (1 - params.Sgr()), 1.);
357 Evaluation scaleFactor = 1.;
359 scaleFactor = (
Sg - params.Sgr())/(0.1 - params.Sgr());
360 if (scaleFactor < 0.)
365 *
pow(1 - Se, 1.0/3.)
366 *
pow(1 -
pow(Se, 1/params.vgM()), 2*params.vgM());
Implementation of three-phase capillary pressure and relative permeability relations proposed by Park...
Definition: ThreePhaseParkerVanGenuchten.hpp:51
static const int gasPhaseIdx
Definition: ThreePhaseParkerVanGenuchten.hpp:64
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition: ThreePhaseParkerVanGenuchten.hpp:238
static void capillaryPressures(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
Implements the three phase capillary pressure law proposed by Parker and van Genuchten.
Definition: ThreePhaseParkerVanGenuchten.hpp:102
static const int nonWettingPhaseIdx
Definition: ThreePhaseParkerVanGenuchten.hpp:63
static const bool isCompositionDependent
Definition: ThreePhaseParkerVanGenuchten.hpp:88
static const int numPhases
Definition: ThreePhaseParkerVanGenuchten.hpp:61
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition: ThreePhaseParkerVanGenuchten.hpp:252
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition: ThreePhaseParkerVanGenuchten.hpp:231
static const bool isPressureDependent
Definition: ThreePhaseParkerVanGenuchten.hpp:80
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition: ThreePhaseParkerVanGenuchten.hpp:245
ParamsT Params
Definition: ThreePhaseParkerVanGenuchten.hpp:58
static const bool implementsTwoPhaseSatApi
Definition: ThreePhaseParkerVanGenuchten.hpp:72
static void saturations(ContainerT &, const Params &, const FluidState &)
The saturation-capillary pressure curve.
Definition: ThreePhaseParkerVanGenuchten.hpp:222
Traits::Scalar Scalar
Definition: ThreePhaseParkerVanGenuchten.hpp:59
static const bool isSaturationDependent
Definition: ThreePhaseParkerVanGenuchten.hpp:76
static Evaluation pcnw(const Params ¶ms, const FluidState &fluidState)
Capillary pressure between the non-wetting liquid (i.e., oil) and the wetting liquid (i....
Definition: ThreePhaseParkerVanGenuchten.hpp:176
static const int wettingPhaseIdx
Definition: ThreePhaseParkerVanGenuchten.hpp:62
static Evaluation krg(const Params ¶ms, const FluidState &fluidState)
The relative permeability for the non-wetting phase of the medium implied by van Genuchten's paramete...
Definition: ThreePhaseParkerVanGenuchten.hpp:345
static Evaluation pcgn(const Params ¶ms, const FluidState &fluidState)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition: ThreePhaseParkerVanGenuchten.hpp:123
static Evaluation krw(const Params ¶ms, const FluidState &fluidState)
The relative permeability for the wetting phase of the medium implied by van Genuchten's parameteriza...
Definition: ThreePhaseParkerVanGenuchten.hpp:272
static Evaluation krn(const Params ¶ms, const FluidState &fluidState)
The relative permeability for the non-wetting phase due to the model of Parker et al....
Definition: ThreePhaseParkerVanGenuchten.hpp:300
TraitsT Traits
Definition: ThreePhaseParkerVanGenuchten.hpp:57
static const bool isTemperatureDependent
Definition: ThreePhaseParkerVanGenuchten.hpp:84
static const bool implementsTwoPhaseApi
Definition: ThreePhaseParkerVanGenuchten.hpp:68
Definition: Air_Mesitylene.hpp:34
Evaluation sqrt(const Evaluation &value)
Definition: MathToolbox.hpp:399
ReturnEval_< Evaluation1, Evaluation2 >::type min(const Evaluation1 &arg1, const Evaluation2 &arg2)
Definition: MathToolbox.hpp:346
ReturnEval_< Evaluation1, Evaluation2 >::type max(const Evaluation1 &arg1, const Evaluation2 &arg2)
Definition: MathToolbox.hpp:341
ReturnEval_< Evaluation1, Evaluation2 >::type pow(const Evaluation1 &base, const Evaluation2 &exp)
Definition: MathToolbox.hpp:416