27 #ifndef OPM_ECL_MULTIPLEXER_MATERIAL_HPP
28 #define OPM_ECL_MULTIPLEXER_MATERIAL_HPP
49 template <
class TraitsT,
50 class GasOilMaterialLawT,
51 class OilWaterMaterialLawT,
52 class GasWaterMaterialLawT,
53 class ParamsT = EclMultiplexerMaterialParams<TraitsT,
56 GasWaterMaterialLawT> >
60 typedef GasOilMaterialLawT GasOilMaterialLaw;
61 typedef OilWaterMaterialLawT OilWaterMaterialLaw;
62 typedef GasWaterMaterialLawT GasWaterMaterialLaw;
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 typedef TraitsT Traits;
88 typedef ParamsT Params;
89 typedef typename Traits::Scalar Scalar;
91 static const int numPhases = 3;
92 static const int waterPhaseIdx = Traits::wettingPhaseIdx;
93 static const int oilPhaseIdx = Traits::nonWettingPhaseIdx;
94 static const int gasPhaseIdx = Traits::gasPhaseIdx;
134 template <
class ContainerT,
class Flu
idState>
136 const Params& params,
137 const FluidState& fluidState)
139 switch (params.approach()) {
140 case EclMultiplexerApproach::EclStone1Approach:
142 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
146 case EclMultiplexerApproach::EclStone2Approach:
148 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
152 case EclMultiplexerApproach::EclDefaultApproach:
154 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
158 case EclMultiplexerApproach::EclTwoPhaseApproach:
160 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
164 case EclMultiplexerApproach::EclOnePhaseApproach:
176 static void oilWaterHysteresisParams(Scalar& pcSwMdc,
178 const Params& params)
180 switch (params.approach()) {
181 case EclMultiplexerApproach::EclStone1Approach:
182 Stone1Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
183 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
186 case EclMultiplexerApproach::EclStone2Approach:
187 Stone2Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
188 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
191 case EclMultiplexerApproach::EclDefaultApproach:
192 DefaultMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
193 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
196 case EclMultiplexerApproach::EclTwoPhaseApproach:
197 TwoPhaseMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
198 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
201 case EclMultiplexerApproach::EclOnePhaseApproach:
213 static void setOilWaterHysteresisParams(
const Scalar& pcSwMdc,
214 const Scalar& krnSwMdc,
217 switch (params.approach()) {
218 case EclMultiplexerApproach::EclStone1Approach:
219 Stone1Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
220 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
223 case EclMultiplexerApproach::EclStone2Approach:
224 Stone2Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
225 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
228 case EclMultiplexerApproach::EclDefaultApproach:
229 DefaultMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
230 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
233 case EclMultiplexerApproach::EclTwoPhaseApproach:
234 TwoPhaseMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
235 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
238 case EclMultiplexerApproach::EclOnePhaseApproach:
250 static void gasOilHysteresisParams(Scalar& pcSwMdc,
252 const Params& params)
254 switch (params.approach()) {
255 case EclMultiplexerApproach::EclStone1Approach:
256 Stone1Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
257 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
260 case EclMultiplexerApproach::EclStone2Approach:
261 Stone2Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
262 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
265 case EclMultiplexerApproach::EclDefaultApproach:
266 DefaultMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
267 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
270 case EclMultiplexerApproach::EclTwoPhaseApproach:
271 TwoPhaseMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
272 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
275 case EclMultiplexerApproach::EclOnePhaseApproach:
287 static void setGasOilHysteresisParams(
const Scalar& pcSwMdc,
288 const Scalar& krnSwMdc,
291 switch (params.approach()) {
292 case EclMultiplexerApproach::EclStone1Approach:
293 Stone1Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
294 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
297 case EclMultiplexerApproach::EclStone2Approach:
298 Stone2Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
299 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
302 case EclMultiplexerApproach::EclDefaultApproach:
303 DefaultMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
304 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
307 case EclMultiplexerApproach::EclTwoPhaseApproach:
308 TwoPhaseMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
309 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
312 case EclMultiplexerApproach::EclOnePhaseApproach:
327 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
328 static Evaluation
pcgn(
const Params& ,
331 throw std::logic_error(
"Not implemented: pcgn()");
343 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
344 static Evaluation
pcnw(
const Params& ,
347 throw std::logic_error(
"Not implemented: pcnw()");
353 template <
class ContainerT,
class Flu
idState>
358 throw std::logic_error(
"Not implemented: saturations()");
364 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
365 static Evaluation
Sg(
const Params& ,
368 throw std::logic_error(
"Not implemented: Sg()");
374 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
375 static Evaluation
Sn(
const Params& ,
378 throw std::logic_error(
"Not implemented: Sn()");
384 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
385 static Evaluation
Sw(
const Params& ,
388 throw std::logic_error(
"Not implemented: Sw()");
406 template <
class ContainerT,
class Flu
idState>
408 const Params& params,
409 const FluidState& fluidState)
411 switch (params.approach()) {
412 case EclMultiplexerApproach::EclStone1Approach:
414 params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
418 case EclMultiplexerApproach::EclStone2Approach:
420 params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
424 case EclMultiplexerApproach::EclDefaultApproach:
426 params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
430 case EclMultiplexerApproach::EclTwoPhaseApproach:
432 params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
436 case EclMultiplexerApproach::EclOnePhaseApproach:
441 throw std::logic_error(
"Not implemented: relativePermeabilities() option for unknown EclMultiplexerApproach (="
442 + std::to_string(
static_cast<int>(params.approach())) +
")");
449 template <
class Evaluation,
class Flu
idState>
451 const FluidState& fluidState)
453 switch (params.approach()) {
454 case EclMultiplexerApproach::EclStone1Approach:
455 return Stone1Material::template relpermOilInOilGasSystem<Evaluation>
456 (params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
459 case EclMultiplexerApproach::EclStone2Approach:
460 return Stone2Material::template relpermOilInOilGasSystem<Evaluation>
461 (params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
464 case EclMultiplexerApproach::EclDefaultApproach:
465 return DefaultMaterial::template relpermOilInOilGasSystem<Evaluation>
466 (params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
470 throw std::logic_error {
471 "relpermOilInOilGasSystem() is specific to three phases"
479 template <
class Evaluation,
class Flu
idState>
481 const FluidState& fluidState)
483 switch (params.approach()) {
484 case EclMultiplexerApproach::EclStone1Approach:
485 return Stone1Material::template relpermOilInOilWaterSystem<Evaluation>
486 (params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
489 case EclMultiplexerApproach::EclStone2Approach:
490 return Stone2Material::template relpermOilInOilWaterSystem<Evaluation>
491 (params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
494 case EclMultiplexerApproach::EclDefaultApproach:
495 return DefaultMaterial::template relpermOilInOilWaterSystem<Evaluation>
496 (params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
500 throw std::logic_error {
501 "relpermOilInOilWaterSystem() is specific to three phases"
509 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
510 static Evaluation
krg(
const Params& ,
513 throw std::logic_error(
"Not implemented: krg()");
519 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
520 static Evaluation
krw(
const Params& ,
523 throw std::logic_error(
"Not implemented: krw()");
529 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
530 static Evaluation
krn(
const Params& ,
533 throw std::logic_error(
"Not implemented: krn()");
544 template <
class Flu
idState>
547 switch (params.approach()) {
548 case EclMultiplexerApproach::EclStone1Approach:
553 case EclMultiplexerApproach::EclStone2Approach:
558 case EclMultiplexerApproach::EclDefaultApproach:
563 case EclMultiplexerApproach::EclTwoPhaseApproach:
567 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.
Some templates to wrap the valgrind client request macros.
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclDefaultMaterial.hpp:59
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:134
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclDefaultMaterial.hpp:308
static void updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclDefaultMaterial.hpp:422
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:480
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclMultiplexerMaterial.hpp:98
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclMultiplexerMaterial.hpp:106
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition: EclMultiplexerMaterial.hpp:520
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition: EclMultiplexerMaterial.hpp:365
static Evaluation relpermOilInOilGasSystem(const Params ¶ms, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition: EclMultiplexerMaterial.hpp:450
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: EclMultiplexerMaterial.hpp:110
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:328
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:344
static void updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclMultiplexerMaterial.hpp:545
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: EclMultiplexerMaterial.hpp:118
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclMultiplexerMaterial.hpp:407
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:375
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition: EclMultiplexerMaterial.hpp:510
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: EclMultiplexerMaterial.hpp:354
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclMultiplexerMaterial.hpp:102
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:530
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition: EclMultiplexerMaterial.hpp:385
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclMultiplexerMaterial.hpp:114
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