27#ifndef OPM_ECL_MULTIPLEXER_MATERIAL_HPP
28#define OPM_ECL_MULTIPLEXER_MATERIAL_HPP
36#include <opm/common/TimingMacros.hpp>
49template <
class TraitsT,
50 class GasOilMaterialLawT,
51 class OilWaterMaterialLawT,
52 class GasWaterMaterialLawT,
53 class ParamsT = EclMultiplexerMaterialParams<TraitsT,
56 GasWaterMaterialLawT> >
60 using GasOilMaterialLaw = GasOilMaterialLawT;
61 using OilWaterMaterialLaw = OilWaterMaterialLawT;
62 using GasWaterMaterialLaw = GasWaterMaterialLawT;
70 static_assert(TraitsT::numPhases == 3,
71 "The number of phases considered by this capillary pressure "
72 "law is always three!");
73 static_assert(GasOilMaterialLaw::numPhases == 2,
74 "The number of phases considered by the gas-oil capillary "
75 "pressure law must be two!");
76 static_assert(OilWaterMaterialLaw::numPhases == 2,
77 "The number of phases considered by the oil-water capillary "
78 "pressure law must be two!");
79 static_assert(GasWaterMaterialLaw::numPhases == 2,
80 "The number of phases considered by the gas-water capillary "
81 "pressure law must be two!");
82 static_assert(std::is_same<
typename GasOilMaterialLaw::Scalar,
83 typename OilWaterMaterialLaw::Scalar>::value,
84 "The two two-phase capillary pressure laws must use the same "
85 "type of floating point values.");
87 using Traits = TraitsT;
88 using Params = ParamsT;
89 using Scalar =
typename Traits::Scalar;
91 static constexpr int numPhases = 3;
92 static constexpr int waterPhaseIdx = Traits::wettingPhaseIdx;
93 static constexpr int oilPhaseIdx = Traits::nonWettingPhaseIdx;
94 static constexpr int gasPhaseIdx = Traits::gasPhaseIdx;
134 template <
class ContainerT,
class Flu
idState>
136 const Params& params,
137 const FluidState& fluidState)
139 OPM_TIMEFUNCTION_LOCAL();
140 switch (params.approach()) {
141 case EclMultiplexerApproach::Stone1:
143 params.template getRealParams<EclMultiplexerApproach::Stone1>(),
147 case EclMultiplexerApproach::Stone2:
149 params.template getRealParams<EclMultiplexerApproach::Stone2>(),
153 case EclMultiplexerApproach::Default:
155 params.template getRealParams<EclMultiplexerApproach::Default>(),
159 case EclMultiplexerApproach::TwoPhase:
161 params.template getRealParams<EclMultiplexerApproach::TwoPhase>(),
165 case EclMultiplexerApproach::OnePhase:
177 static void oilWaterHysteresisParams(Scalar& pcSwMdc,
179 const Params& params)
181 OPM_TIMEFUNCTION_LOCAL();
182 switch (params.approach()) {
183 case EclMultiplexerApproach::Stone1:
184 Stone1Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
185 params.template getRealParams<EclMultiplexerApproach::Stone1>());
188 case EclMultiplexerApproach::Stone2:
189 Stone2Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
190 params.template getRealParams<EclMultiplexerApproach::Stone2>());
193 case EclMultiplexerApproach::Default:
194 DefaultMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
195 params.template getRealParams<EclMultiplexerApproach::Default>());
198 case EclMultiplexerApproach::TwoPhase:
199 TwoPhaseMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
200 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
203 case EclMultiplexerApproach::OnePhase:
215 static void setOilWaterHysteresisParams(
const Scalar& pcSwMdc,
216 const Scalar& krnSwMdc,
219 OPM_TIMEFUNCTION_LOCAL();
220 switch (params.approach()) {
221 case EclMultiplexerApproach::Stone1:
222 Stone1Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
223 params.template getRealParams<EclMultiplexerApproach::Stone1>());
226 case EclMultiplexerApproach::Stone2:
227 Stone2Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
228 params.template getRealParams<EclMultiplexerApproach::Stone2>());
231 case EclMultiplexerApproach::Default:
232 DefaultMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
233 params.template getRealParams<EclMultiplexerApproach::Default>());
236 case EclMultiplexerApproach::TwoPhase:
237 TwoPhaseMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
238 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
241 case EclMultiplexerApproach::OnePhase:
253 static void gasOilHysteresisParams(Scalar& pcSwMdc,
255 const Params& params)
257 OPM_TIMEFUNCTION_LOCAL();
258 switch (params.approach()) {
259 case EclMultiplexerApproach::Stone1:
260 Stone1Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
261 params.template getRealParams<EclMultiplexerApproach::Stone1>());
264 case EclMultiplexerApproach::Stone2:
265 Stone2Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
266 params.template getRealParams<EclMultiplexerApproach::Stone2>());
269 case EclMultiplexerApproach::Default:
270 DefaultMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
271 params.template getRealParams<EclMultiplexerApproach::Default>());
274 case EclMultiplexerApproach::TwoPhase:
275 TwoPhaseMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
276 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
279 case EclMultiplexerApproach::OnePhase:
285 static Scalar trappedGasSaturation(
const Params& params)
287 OPM_TIMEFUNCTION_LOCAL();
288 switch (params.approach()) {
289 case EclMultiplexerApproach::Stone1:
290 return params.template getRealParams<EclMultiplexerApproach::Stone1>().gasOilParams().SnTrapped();
291 case EclMultiplexerApproach::Stone2:
292 return params.template getRealParams<EclMultiplexerApproach::Stone2>().gasOilParams().SnTrapped();
293 case EclMultiplexerApproach::Default:
294 return params.template getRealParams<EclMultiplexerApproach::Default>().gasOilParams().SnTrapped();
295 case EclMultiplexerApproach::TwoPhase:
296 if(params.template getRealParams<EclMultiplexerApproach::TwoPhase>().approach() == EclTwoPhaseApproach::GasOil)
297 return params.template getRealParams<EclMultiplexerApproach::TwoPhase>().gasOilParams().SnTrapped();
298 if(params.template getRealParams<EclMultiplexerApproach::TwoPhase>().approach() == EclTwoPhaseApproach::GasWater)
299 return params.template getRealParams<EclMultiplexerApproach::TwoPhase>().gasWaterParams().SnTrapped();
301 case EclMultiplexerApproach::OnePhase:
313 static void setGasOilHysteresisParams(
const Scalar& pcSwMdc,
314 const Scalar& krnSwMdc,
317 OPM_TIMEFUNCTION_LOCAL();
318 switch (params.approach()) {
319 case EclMultiplexerApproach::Stone1:
320 Stone1Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
321 params.template getRealParams<EclMultiplexerApproach::Stone1>());
324 case EclMultiplexerApproach::Stone2:
325 Stone2Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
326 params.template getRealParams<EclMultiplexerApproach::Stone2>());
329 case EclMultiplexerApproach::Default:
330 DefaultMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
331 params.template getRealParams<EclMultiplexerApproach::Default>());
334 case EclMultiplexerApproach::TwoPhase:
335 TwoPhaseMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
336 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
339 case EclMultiplexerApproach::OnePhase:
354 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
355 static Evaluation
pcgn(
const Params& ,
358 throw std::logic_error(
"Not implemented: pcgn()");
370 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
371 static Evaluation
pcnw(
const Params& ,
374 throw std::logic_error(
"Not implemented: pcnw()");
380 template <
class ContainerT,
class Flu
idState>
385 throw std::logic_error(
"Not implemented: saturations()");
391 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
392 static Evaluation
Sg(
const Params& ,
395 throw std::logic_error(
"Not implemented: Sg()");
401 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
402 static Evaluation
Sn(
const Params& ,
405 throw std::logic_error(
"Not implemented: Sn()");
411 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
412 static Evaluation
Sw(
const Params& ,
415 throw std::logic_error(
"Not implemented: Sw()");
433 template <
class ContainerT,
class Flu
idState>
435 const Params& params,
436 const FluidState& fluidState)
438 OPM_TIMEFUNCTION_LOCAL();
439 switch (params.approach()) {
440 case EclMultiplexerApproach::Stone1:
442 params.template getRealParams<EclMultiplexerApproach::Stone1>(),
446 case EclMultiplexerApproach::Stone2:
448 params.template getRealParams<EclMultiplexerApproach::Stone2>(),
452 case EclMultiplexerApproach::Default:
454 params.template getRealParams<EclMultiplexerApproach::Default>(),
458 case EclMultiplexerApproach::TwoPhase:
460 params.template getRealParams<EclMultiplexerApproach::TwoPhase>(),
464 case EclMultiplexerApproach::OnePhase:
469 throw std::logic_error(
"Not implemented: relativePermeabilities() option for unknown EclMultiplexerApproach (="
470 + std::to_string(
static_cast<int>(params.approach())) +
")");
477 template <
class Evaluation,
class Flu
idState>
479 const FluidState& fluidState)
481 OPM_TIMEFUNCTION_LOCAL();
482 switch (params.approach()) {
483 case EclMultiplexerApproach::Stone1:
484 return Stone1Material::template relpermOilInOilGasSystem<Evaluation>
485 (params.template getRealParams<EclMultiplexerApproach::Stone1>(),
488 case EclMultiplexerApproach::Stone2:
489 return Stone2Material::template relpermOilInOilGasSystem<Evaluation>
490 (params.template getRealParams<EclMultiplexerApproach::Stone2>(),
493 case EclMultiplexerApproach::Default:
494 return DefaultMaterial::template relpermOilInOilGasSystem<Evaluation>
495 (params.template getRealParams<EclMultiplexerApproach::Default>(),
499 throw std::logic_error {
500 "relpermOilInOilGasSystem() is specific to three phases"
508 template <
class Evaluation,
class Flu
idState>
510 const FluidState& fluidState)
512 OPM_TIMEFUNCTION_LOCAL();
513 switch (params.approach()) {
514 case EclMultiplexerApproach::Stone1:
515 return Stone1Material::template relpermOilInOilWaterSystem<Evaluation>
516 (params.template getRealParams<EclMultiplexerApproach::Stone1>(),
519 case EclMultiplexerApproach::Stone2:
520 return Stone2Material::template relpermOilInOilWaterSystem<Evaluation>
521 (params.template getRealParams<EclMultiplexerApproach::Stone2>(),
524 case EclMultiplexerApproach::Default:
525 return DefaultMaterial::template relpermOilInOilWaterSystem<Evaluation>
526 (params.template getRealParams<EclMultiplexerApproach::Default>(),
530 throw std::logic_error {
531 "relpermOilInOilWaterSystem() is specific to three phases"
539 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
540 static Evaluation
krg(
const Params& ,
543 throw std::logic_error(
"Not implemented: krg()");
549 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
550 static Evaluation
krw(
const Params& ,
553 throw std::logic_error(
"Not implemented: krw()");
559 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
560 static Evaluation
krn(
const Params& ,
563 throw std::logic_error(
"Not implemented: krn()");
574 template <
class Flu
idState>
577 OPM_TIMEFUNCTION_LOCAL();
578 switch (params.approach()) {
579 case EclMultiplexerApproach::Stone1:
584 case EclMultiplexerApproach::Stone2:
589 case EclMultiplexerApproach::Default:
594 case EclMultiplexerApproach::TwoPhase:
598 case EclMultiplexerApproach::OnePhase:
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 bool updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclDefaultMaterial.hpp:439
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition EclDefaultMaterial.hpp:319
Implements a multiplexer class that provides all three phase capillary pressure laws used by the ECLi...
Definition EclMultiplexerMaterial.hpp:58
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:135
static Evaluation relpermOilInOilWaterSystem(const Params ¶ms, const FluidState &fluidState)
The relative permeability of oil in oil/water system.
Definition EclMultiplexerMaterial.hpp:509
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition EclMultiplexerMaterial.hpp:550
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition EclMultiplexerMaterial.hpp:392
static Evaluation relpermOilInOilGasSystem(const Params ¶ms, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition EclMultiplexerMaterial.hpp:478
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition EclMultiplexerMaterial.hpp:355
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:371
static constexpr bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition EclMultiplexerMaterial.hpp:110
static constexpr bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition EclMultiplexerMaterial.hpp:102
static constexpr bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition EclMultiplexerMaterial.hpp:106
static constexpr bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition EclMultiplexerMaterial.hpp:114
static constexpr bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition EclMultiplexerMaterial.hpp:98
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition EclMultiplexerMaterial.hpp:434
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition EclMultiplexerMaterial.hpp:402
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition EclMultiplexerMaterial.hpp:540
static constexpr bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition EclMultiplexerMaterial.hpp:118
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition EclMultiplexerMaterial.hpp:381
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition EclMultiplexerMaterial.hpp:560
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition EclMultiplexerMaterial.hpp:412
static bool updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclMultiplexerMaterial.hpp:575
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 bool 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 bool 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:343
static bool updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclTwoPhaseMaterial.hpp:420
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:145
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30