27 #ifndef OPM_ECL_MULTIPLEXER_MATERIAL_HPP
28 #define OPM_ECL_MULTIPLEXER_MATERIAL_HPP
47 template <
class TraitsT,
48 class GasOilMaterialLawT,
49 class OilWaterMaterialLawT,
50 class GasWaterMaterialLawT,
51 class ParamsT = EclMultiplexerMaterialParams<TraitsT,
54 GasWaterMaterialLawT> >
58 using GasOilMaterialLaw = GasOilMaterialLawT;
59 using OilWaterMaterialLaw = OilWaterMaterialLawT;
60 using GasWaterMaterialLaw = GasWaterMaterialLawT;
68 static_assert(TraitsT::numPhases == 3,
69 "The number of phases considered by this capillary pressure "
70 "law is always three!");
71 static_assert(GasOilMaterialLaw::numPhases == 2,
72 "The number of phases considered by the gas-oil capillary "
73 "pressure law must be two!");
74 static_assert(OilWaterMaterialLaw::numPhases == 2,
75 "The number of phases considered by the oil-water capillary "
76 "pressure law must be two!");
77 static_assert(GasWaterMaterialLaw::numPhases == 2,
78 "The number of phases considered by the gas-water capillary "
79 "pressure law must be two!");
80 static_assert(std::is_same<
typename GasOilMaterialLaw::Scalar,
81 typename OilWaterMaterialLaw::Scalar>::value,
82 "The two two-phase capillary pressure laws must use the same "
83 "type of floating point values.");
85 using Traits = TraitsT;
86 using Params = ParamsT;
87 using Scalar =
typename Traits::Scalar;
89 static constexpr
int numPhases = 3;
90 static constexpr
int waterPhaseIdx = Traits::wettingPhaseIdx;
91 static constexpr
int oilPhaseIdx = Traits::nonWettingPhaseIdx;
92 static constexpr
int gasPhaseIdx = Traits::gasPhaseIdx;
132 template <
class ContainerT,
class Flu
idState>
134 const Params& params,
135 const FluidState& fluidState)
137 switch (params.approach()) {
138 case EclMultiplexerApproach::EclStone1Approach:
140 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
144 case EclMultiplexerApproach::EclStone2Approach:
146 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
150 case EclMultiplexerApproach::EclDefaultApproach:
152 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
156 case EclMultiplexerApproach::EclTwoPhaseApproach:
158 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
162 case EclMultiplexerApproach::EclOnePhaseApproach:
174 static void oilWaterHysteresisParams(Scalar& pcSwMdc,
176 const Params& params)
178 switch (params.approach()) {
179 case EclMultiplexerApproach::EclStone1Approach:
180 Stone1Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
181 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
184 case EclMultiplexerApproach::EclStone2Approach:
185 Stone2Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
186 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
189 case EclMultiplexerApproach::EclDefaultApproach:
190 DefaultMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
191 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
194 case EclMultiplexerApproach::EclTwoPhaseApproach:
195 TwoPhaseMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
196 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
199 case EclMultiplexerApproach::EclOnePhaseApproach:
211 static void setOilWaterHysteresisParams(
const Scalar& pcSwMdc,
212 const Scalar& krnSwMdc,
215 switch (params.approach()) {
216 case EclMultiplexerApproach::EclStone1Approach:
217 Stone1Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
218 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
221 case EclMultiplexerApproach::EclStone2Approach:
222 Stone2Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
223 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
226 case EclMultiplexerApproach::EclDefaultApproach:
227 DefaultMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
228 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
231 case EclMultiplexerApproach::EclTwoPhaseApproach:
232 TwoPhaseMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
233 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
236 case EclMultiplexerApproach::EclOnePhaseApproach:
248 static void gasOilHysteresisParams(Scalar& pcSwMdc,
250 const Params& params)
252 switch (params.approach()) {
253 case EclMultiplexerApproach::EclStone1Approach:
254 Stone1Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
255 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
258 case EclMultiplexerApproach::EclStone2Approach:
259 Stone2Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
260 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
263 case EclMultiplexerApproach::EclDefaultApproach:
264 DefaultMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
265 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
268 case EclMultiplexerApproach::EclTwoPhaseApproach:
269 TwoPhaseMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
270 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
273 case EclMultiplexerApproach::EclOnePhaseApproach:
285 static void setGasOilHysteresisParams(
const Scalar& pcSwMdc,
286 const Scalar& krnSwMdc,
289 switch (params.approach()) {
290 case EclMultiplexerApproach::EclStone1Approach:
291 Stone1Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
292 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
295 case EclMultiplexerApproach::EclStone2Approach:
296 Stone2Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
297 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
300 case EclMultiplexerApproach::EclDefaultApproach:
301 DefaultMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
302 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
305 case EclMultiplexerApproach::EclTwoPhaseApproach:
306 TwoPhaseMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
307 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
310 case EclMultiplexerApproach::EclOnePhaseApproach:
325 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
326 static Evaluation
pcgn(
const Params& ,
329 throw std::logic_error(
"Not implemented: pcgn()");
341 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
342 static Evaluation
pcnw(
const Params& ,
345 throw std::logic_error(
"Not implemented: pcnw()");
351 template <
class ContainerT,
class Flu
idState>
356 throw std::logic_error(
"Not implemented: saturations()");
362 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
363 static Evaluation
Sg(
const Params& ,
366 throw std::logic_error(
"Not implemented: Sg()");
372 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
373 static Evaluation
Sn(
const Params& ,
376 throw std::logic_error(
"Not implemented: Sn()");
382 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
383 static Evaluation
Sw(
const Params& ,
386 throw std::logic_error(
"Not implemented: Sw()");
404 template <
class ContainerT,
class Flu
idState>
406 const Params& params,
407 const FluidState& fluidState)
409 switch (params.approach()) {
410 case EclMultiplexerApproach::EclStone1Approach:
412 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
416 case EclMultiplexerApproach::EclStone2Approach:
418 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
422 case EclMultiplexerApproach::EclDefaultApproach:
424 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
428 case EclMultiplexerApproach::EclTwoPhaseApproach:
430 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
434 case EclMultiplexerApproach::EclOnePhaseApproach:
439 throw std::logic_error(
"Not implemented: relativePermeabilities() option for unknown EclMultiplexerApproach (="
440 + std::to_string(
static_cast<int>(params.approach())) +
")");
447 template <
class Evaluation,
class Flu
idState>
449 const FluidState& fluidState)
451 switch (params.approach()) {
452 case EclMultiplexerApproach::EclStone1Approach:
453 return Stone1Material::template relpermOilInOilGasSystem<Evaluation>
454 (params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
457 case EclMultiplexerApproach::EclStone2Approach:
458 return Stone2Material::template relpermOilInOilGasSystem<Evaluation>
459 (params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
462 case EclMultiplexerApproach::EclDefaultApproach:
463 return DefaultMaterial::template relpermOilInOilGasSystem<Evaluation>
464 (params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
468 throw std::logic_error {
469 "relpermOilInOilGasSystem() is specific to three phases"
477 template <
class Evaluation,
class Flu
idState>
479 const FluidState& fluidState)
481 switch (params.approach()) {
482 case EclMultiplexerApproach::EclStone1Approach:
483 return Stone1Material::template relpermOilInOilWaterSystem<Evaluation>
484 (params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
487 case EclMultiplexerApproach::EclStone2Approach:
488 return Stone2Material::template relpermOilInOilWaterSystem<Evaluation>
489 (params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
492 case EclMultiplexerApproach::EclDefaultApproach:
493 return DefaultMaterial::template relpermOilInOilWaterSystem<Evaluation>
494 (params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
498 throw std::logic_error {
499 "relpermOilInOilWaterSystem() is specific to three phases"
507 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
508 static Evaluation
krg(
const Params& ,
511 throw std::logic_error(
"Not implemented: krg()");
517 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
518 static Evaluation
krw(
const Params& ,
521 throw std::logic_error(
"Not implemented: krw()");
527 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
528 static Evaluation
krn(
const Params& ,
531 throw std::logic_error(
"Not implemented: krn()");
542 template <
class Flu
idState>
545 switch (params.approach()) {
546 case EclMultiplexerApproach::EclStone1Approach:
551 case EclMultiplexerApproach::EclStone2Approach:
556 case EclMultiplexerApproach::EclDefaultApproach:
561 case EclMultiplexerApproach::EclTwoPhaseApproach:
565 case EclMultiplexerApproach::EclOnePhaseApproach:
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Multiplexer implementation for the parameters required by the multiplexed three-phase material law.
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclDefaultMaterial.hpp:61
static void capillaryPressures(ContainerT &values, const Params ¶ms, const FluidState &state)
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclDefaultMaterial.hpp:136
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclDefaultMaterial.hpp:310
static void updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclDefaultMaterial.hpp:424
Implements a multiplexer class that provides all three phase capillary pressure laws used by the ECLi...
Definition: EclMultiplexerMaterial.hpp:56
static void capillaryPressures(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
Implements the multiplexer three phase capillary pressure law used by the ECLipse simulator.
Definition: EclMultiplexerMaterial.hpp:133
static Evaluation relpermOilInOilWaterSystem(const Params ¶ms, const FluidState &fluidState)
The relative permeability of oil in oil/water system.
Definition: EclMultiplexerMaterial.hpp:478
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition: EclMultiplexerMaterial.hpp:518
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition: EclMultiplexerMaterial.hpp:363
static Evaluation relpermOilInOilGasSystem(const Params ¶ms, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition: EclMultiplexerMaterial.hpp:448
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:326
static Evaluation pcnw(const Params &, const FluidState &)
Capillary pressure between the non-wetting liquid (i.e., oil) and the wetting liquid (i....
Definition: EclMultiplexerMaterial.hpp:342
static void updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclMultiplexerMaterial.hpp:543
static constexpr bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: EclMultiplexerMaterial.hpp:108
static constexpr bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclMultiplexerMaterial.hpp:100
static constexpr bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclMultiplexerMaterial.hpp:104
static constexpr bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclMultiplexerMaterial.hpp:112
static constexpr bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclMultiplexerMaterial.hpp:96
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclMultiplexerMaterial.hpp:405
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:373
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition: EclMultiplexerMaterial.hpp:508
static constexpr bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: EclMultiplexerMaterial.hpp:116
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: EclMultiplexerMaterial.hpp:352
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:528
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition: EclMultiplexerMaterial.hpp:383
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition: EclStone1Material.hpp:60
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclStone1Material.hpp:313
static void capillaryPressures(ContainerT &values, const Params ¶ms, const FluidState &state)
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclStone1Material.hpp:135
static void updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclStone1Material.hpp:420
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition: EclStone2Material.hpp:61
static void updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclStone2Material.hpp:404
static void capillaryPressures(ContainerT &values, const Params ¶ms, const FluidState &state)
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclStone2Material.hpp:136
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclStone2Material.hpp:314
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Definition: EclTwoPhaseMaterial.hpp:57
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclTwoPhaseMaterial.hpp:317
static void updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclTwoPhaseMaterial.hpp:393
static void capillaryPressures(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
Implements the multiplexer three phase capillary pressure law used by the ECLipse simulator.
Definition: EclTwoPhaseMaterial.hpp:129