27#ifndef OPM_REGULARIZED_VAN_GENUCHTEN_HPP
28#define OPM_REGULARIZED_VAN_GENUCHTEN_HPP
70template <
class TraitsT,
class ParamsT = RegularizedVanGenuchtenParams<TraitsT> >
73 typedef ::Opm::VanGenuchten<TraitsT, ParamsT>
VanGenuchten;
76 typedef TraitsT Traits;
77 typedef ParamsT Params;
79 typedef typename Traits::Scalar Scalar;
84 "The regularized van Genuchten capillary pressure law only "
85 "applies to the case of two fluid phases");
115 template <
class Container,
class Flu
idState>
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& params,
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>
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& params,
const FluidState& fs)
165 const auto&
Sw = decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
166 return twoPhaseSatPcnw(params,
Sw);
169 template <
class Evaluation>
170 static Evaluation twoPhaseSatPcnw(
const Params& params,
const Evaluation&
Sw)
174 const Scalar SwThLow = params.pcnwLowSw();
175 const Scalar SwThHigh = params.pcnwHighSw();
183 return params.pcnwLow() + params.pcnwSlopeLow()*(
Sw - SwThLow);
185 else if (
Sw > SwThHigh)
187 Scalar yTh = params.pcnwHigh();
188 Scalar m1 = (0.0 - yTh)/(1.0 - SwThHigh)*2;
192 const Spline<Scalar>& sp = params.pcnwHighSpline();
198 return m1*(
Sw - 1.0) + 0.0;
221 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
222 static Evaluation
Sw(
const Params& params,
const FluidState& fs)
224 const Evaluation& pC =
225 decay<Evaluation>(fs.pressure(Traits::nonWettingPhaseIdx))
226 - decay<Evaluation>(fs.pressure(Traits::wettingPhaseIdx));
227 return twoPhaseSatSw(params, pC);
230 template <
class Evaluation>
231 static Evaluation twoPhaseSatSw(
const Params& params,
const Evaluation& pC)
235 const Scalar SwThLow = params.pcnwLowSw();
236 const Scalar SwThHigh = params.pcnwHighSw();
243 Scalar m1 = params.pcnwSlopeHigh();
247 Evaluation
Sw = VanGenuchten::twoPhaseSatSw(params, pC);
253 return (pC - pC_SwLow)/params.pcnwSlopeLow() + SwThLow;
255 else if (SwThHigh <
Sw )
258 const Spline<Scalar>& spline = params.pcnwHighSpline();
260 return spline.template intersectInterval<Evaluation>(SwThHigh, 1.0,
271 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
272 static Evaluation
Sn(
const Params& params,
const FluidState& fs)
273 {
return 1 - Sw<FluidState, Evaluation>(params, fs); }
275 template <
class Evaluation>
276 static Evaluation twoPhaseSatSn(
const Params& params,
const Evaluation& pc)
277 {
return 1 - twoPhaseSatSw(params, pc); }
293 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
294 static Evaluation
krw(
const Params& params,
const FluidState& fs)
296 const auto&
Sw = decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
297 return twoPhaseSatKrw(params,
Sw);
300 template <
class Evaluation>
301 static Evaluation twoPhaseSatKrw(
const Params& params,
const Evaluation&
Sw)
309 return VanGenuchten::twoPhaseSatKrw(params,
Sw);
326 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
327 static Evaluation
krn(
const Params& params,
const FluidState& fs)
329 const Evaluation&
Sw =
330 1.0 - decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
331 return twoPhaseSatKrn(params,
Sw);
334 template <
class Evaluation>
335 static Evaluation twoPhaseSatKrn(
const Params& params,
const Evaluation&
Sw)
343 return VanGenuchten::twoPhaseSatKrn(params,
Sw);
Parameters that are necessary for the regularization of VanGenuchten "material law".
Class implementing cubic splines.
Implementation of the van Genuchten capillary pressure - saturation relation.
Implementation of the regularized van Genuchten's capillary pressure / relative permeability <-> satu...
Definition RegularizedVanGenuchten.hpp:72
static const int numPhases
The number of fluid phases.
Definition RegularizedVanGenuchten.hpp:82
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition RegularizedVanGenuchten.hpp:109
static Evaluation krn(const Params ¶ms, const FluidState &fs)
Regularized version of the relative permeability for the non-wetting phase of the medium implied by t...
Definition RegularizedVanGenuchten.hpp:327
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition RegularizedVanGenuchten.hpp:89
static void saturations(Container &values, const Params ¶ms, const FluidState &fs)
Calculate the saturations of the phases starting from their pressure differences.
Definition RegularizedVanGenuchten.hpp:129
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition RegularizedVanGenuchten.hpp:97
static Evaluation Sn(const Params ¶ms, const FluidState &fs)
Calculate the non-wetting phase saturations depending on the phase pressures.
Definition RegularizedVanGenuchten.hpp:272
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition RegularizedVanGenuchten.hpp:101
static Evaluation pcnw(const Params ¶ms, const FluidState &fs)
A regularized van Genuchten capillary pressure-saturation curve.
Definition RegularizedVanGenuchten.hpp:163
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition RegularizedVanGenuchten.hpp:105
static Evaluation Sw(const Params ¶ms, const FluidState &fs)
A regularized van Genuchten saturation-capillary pressure curve.
Definition RegularizedVanGenuchten.hpp:222
static void relativePermeabilities(Container &values, const Params ¶ms, const FluidState &fs)
Returns the relative permeabilities of the phases dependening on the phase saturations.
Definition RegularizedVanGenuchten.hpp:142
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition RegularizedVanGenuchten.hpp:93
static Evaluation krw(const Params ¶ms, const FluidState &fs)
Regularized version of the relative permeability for the wetting phase of the medium implied by the v...
Definition RegularizedVanGenuchten.hpp:294
static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fs)
Calculate the pressure difference of the phases in the most generic way.
Definition RegularizedVanGenuchten.hpp:116
Implementation of the van Genuchten capillary pressure - saturation relation.
Definition VanGenuchten.hpp:56
static Evaluation twoPhaseSatPcnw(const Params ¶ms, const Evaluation &Sw)
The saturation-capillary pressure curve according to van Genuchten using a material law specific API.
Definition VanGenuchten.hpp:194
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30