27 #ifndef OPM_SAT_CURVE_MULTIPLEXER_HPP
28 #define OPM_SAT_CURVE_MULTIPLEXER_HPP
42 template <
class TraitsT,
class ParamsT = SatCurveMultiplexerParams<TraitsT> >
46 using Traits = TraitsT;
47 using Params = ParamsT;
48 using Scalar =
typename Traits::Scalar;
56 "The Brooks-Corey capillary pressure law only applies "
57 "to the case of two fluid phases");
83 static_assert(Traits::numPhases == 2,
84 "The number of fluid phases must be two if you want to use "
85 "this material law!");
90 template <
class Container,
class Flu
idState>
91 static void capillaryPressures(Container& values,
const Params& params,
const FluidState& fluidState)
93 switch (params.approach()) {
94 case SatCurveMultiplexerApproach::LETApproach:
96 params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
100 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
102 params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
112 template <
class Container,
class Flu
idState>
113 static void saturations(Container& values,
const Params& params,
const FluidState& fluidState)
115 switch (params.approach()) {
116 case SatCurveMultiplexerApproach::LETApproach:
118 params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
122 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
124 params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
140 template <
class Container,
class Flu
idState>
143 switch (params.approach()) {
144 case SatCurveMultiplexerApproach::LETApproach:
146 params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
150 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
152 params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
161 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
162 static Evaluation
pcnw(
const Params& params,
const FluidState& fluidState)
164 switch (params.approach()) {
165 case SatCurveMultiplexerApproach::LETApproach:
166 return LETTwoPhaseLaw::pcnw(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
170 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
171 return PLTwoPhaseLaw::pcnw(params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
179 template <
class Evaluation>
180 static Evaluation twoPhaseSatPcnw(
const Params& params,
const Evaluation&
Sw)
182 switch (params.approach()) {
183 case SatCurveMultiplexerApproach::LETApproach:
184 return LETTwoPhaseLaw::twoPhaseSatPcnw(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
188 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
197 template <
class Evaluation>
198 static Evaluation twoPhaseSatPcnwInv(
const Params&,
const Evaluation&)
200 throw std::logic_error(
"SatCurveMultiplexer::twoPhaseSatPcnwInv"
201 " not implemented!");
207 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
208 static Evaluation
Sw(
const Params& params,
const FluidState& fluidstate)
210 switch (params.approach()) {
211 case SatCurveMultiplexerApproach::LETApproach:
212 return LETTwoPhaseLaw::Sw(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
216 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
217 return PLTwoPhaseLaw::Sw(params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
225 template <
class Evaluation>
226 static Evaluation twoPhaseSatSw(
const Params&,
const Evaluation&)
228 throw std::logic_error(
"SatCurveMultiplexer::twoPhaseSatSw"
229 " not implemented!");
236 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
237 static Evaluation
Sn(
const Params& params,
const FluidState& fluidstate)
239 switch (params.approach()) {
240 case SatCurveMultiplexerApproach::LETApproach:
241 return LETTwoPhaseLaw::Sn(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
245 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
246 return PLTwoPhaseLaw::Sn(params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
254 template <
class Evaluation>
255 static Evaluation twoPhaseSatSn(
const Params& params,
const Evaluation& pc)
257 switch (params.approach()) {
258 case SatCurveMultiplexerApproach::LETApproach:
259 return LETTwoPhaseLaw::twoPhaseSatSn(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
263 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
264 return PLTwoPhaseLaw::twoPhaseSatSn(params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
276 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
277 static Evaluation
krw(
const Params& params,
const FluidState& fluidstate)
279 switch (params.approach()) {
280 case SatCurveMultiplexerApproach::LETApproach:
281 return LETTwoPhaseLaw::krw(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
285 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
286 return PLTwoPhaseLaw::krw(params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
294 template <
class Evaluation>
295 static Evaluation twoPhaseSatKrw(
const Params& params,
const Evaluation&
Sw)
297 switch (params.approach()) {
298 case SatCurveMultiplexerApproach::LETApproach:
299 return LETTwoPhaseLaw::twoPhaseSatKrw(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
303 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
304 return PLTwoPhaseLaw::twoPhaseSatKrw(params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
312 template <
class Evaluation>
313 static Evaluation twoPhaseSatKrwInv(
const Params&,
const Evaluation&)
315 throw std::logic_error(
"Not implemented: twoPhaseSatKrwInv()");
322 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
323 static Evaluation
krn(
const Params& params,
const FluidState& fluidstate)
325 switch (params.approach()) {
326 case SatCurveMultiplexerApproach::LETApproach:
327 return LETTwoPhaseLaw::krn(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
331 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
332 return PLTwoPhaseLaw::krn(params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
340 template <
class Evaluation>
341 static Evaluation twoPhaseSatKrn(
const Params& params,
const Evaluation&
Sw)
343 switch (params.approach()) {
344 case SatCurveMultiplexerApproach::LETApproach:
345 return LETTwoPhaseLaw::twoPhaseSatKrn(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
349 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
350 return PLTwoPhaseLaw::twoPhaseSatKrn(params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
358 template <
class Evaluation>
359 static Evaluation twoPhaseSatKrnInv(
const Params& params,
const Evaluation&
krn)
361 switch (params.approach()) {
362 case SatCurveMultiplexerApproach::LETApproach:
363 return LETTwoPhaseLaw::twoPhaseSatKrnInv(params.template getRealParams<SatCurveMultiplexerApproach::LETApproach>(),
367 case SatCurveMultiplexerApproach::PiecewiseLinearApproach:
368 return PLTwoPhaseLaw::twoPhaseSatKrnInv(params.template getRealParams<SatCurveMultiplexerApproach::PiecewiseLinearApproach>(),
Specification of the material parameters for the saturation function multiplexer.
Implementation of a tabulated, piecewise linear capillary pressure law.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:50
static Evaluation krn(const Params ¶ms, const FluidState &fs)
The relative permeability for the non-wetting phase of the porous medium.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:197
static Evaluation twoPhaseSatPcnw(const Params ¶ms, const Evaluation &Sw)
The saturation-capillary pressure curve.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:141
static Evaluation krw(const Params ¶ms, const FluidState &fs)
The relative permeability for the wetting phase of the porous medium.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:176
static void relativePermeabilities(Container &values, const Params ¶ms, const FluidState &fs)
The relative permeabilities.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:117
static Evaluation Sn(const Params ¶ms, const FluidState &fs)
Calculate the non-wetting phase saturations depending on the phase pressures.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:164
static void saturations(Container &, const Params &, const FluidState &)
The saturations of the fluid phases starting from their pressure differences.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:110
static Evaluation pcnw(const Params ¶ms, const FluidState &fs)
The capillary pressure-saturation curve.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:129
static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fs)
The capillary pressure-saturation curve.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:97
static Evaluation Sw(const Params &, const FluidState &)
The saturation-capillary pressure curve.
Definition: PiecewiseLinearTwoPhaseMaterial.hpp:152
Implements a multiplexer class that provides LET curves and piecewise linear saturation functions.
Definition: SatCurveMultiplexer.hpp:44
static Evaluation Sw(const Params ¶ms, const FluidState &fluidstate)
The saturation-capillary pressure curve.
Definition: SatCurveMultiplexer.hpp:208
static Evaluation krn(const Params ¶ms, const FluidState &fluidstate)
The relative permeability for the non-wetting phase of the medium.
Definition: SatCurveMultiplexer.hpp:323
static Evaluation krw(const Params ¶ms, const FluidState &fluidstate)
The relative permeability for the wetting phase of the medium.
Definition: SatCurveMultiplexer.hpp:277
static constexpr bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: SatCurveMultiplexer.hpp:73
static void saturations(Container &values, const Params ¶ms, const FluidState &fluidState)
Calculate the saturations of the phases starting from their pressure differences.
Definition: SatCurveMultiplexer.hpp:113
static constexpr bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: SatCurveMultiplexer.hpp:65
static constexpr bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: SatCurveMultiplexer.hpp:77
static constexpr bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: SatCurveMultiplexer.hpp:81
static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fluidState)
The capillary pressure-saturation curves.
Definition: SatCurveMultiplexer.hpp:91
static constexpr int numPhases
The number of fluid phases to which this material law applies.
Definition: SatCurveMultiplexer.hpp:54
static Evaluation Sn(const Params ¶ms, const FluidState &fluidstate)
Calculate the non-wetting phase saturations depending on the phase pressures.
Definition: SatCurveMultiplexer.hpp:237
static void relativePermeabilities(Container &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability-saturation curves.
Definition: SatCurveMultiplexer.hpp:141
static constexpr bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: SatCurveMultiplexer.hpp:61
static constexpr bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: SatCurveMultiplexer.hpp:69
static Evaluation pcnw(const Params ¶ms, const FluidState &fluidState)
The capillary pressure-saturation curve.
Definition: SatCurveMultiplexer.hpp:162
Implementation of the LET curve saturation functions.
Definition: TwoPhaseLETCurves.hpp:49
static Evaluation krn(const Params &, const FluidState &)
The relative permeability for the non-wetting phase of the medium as implied by the LET parameterizat...
Definition: TwoPhaseLETCurves.hpp:236
static Evaluation pcnw(const Params &, const FluidState &)
The capillary pressure-saturation curve.
Definition: TwoPhaseLETCurves.hpp:129
static Evaluation krw(const Params &, const FluidState &)
The relative permeability for the wetting phase of the medium implied by the LET parameterization.
Definition: TwoPhaseLETCurves.hpp:192
static void saturations(Container &, const Params &, const FluidState &)
Calculate the saturations of the phases starting from their pressure differences.
Definition: TwoPhaseLETCurves.hpp:102
static void capillaryPressures(Container &, const Params &, const FluidState &)
The capillary pressure-saturation curves.
Definition: TwoPhaseLETCurves.hpp:92
static void relativePermeabilities(Container &, const Params &, const FluidState &)
The relative permeability-saturation curves.
Definition: TwoPhaseLETCurves.hpp:118