27#ifndef OPM_ECL_TWO_PHASE_MATERIAL_HPP
28#define OPM_ECL_TWO_PHASE_MATERIAL_HPP
48template <
class TraitsT,
49 class GasOilMaterialLawT,
50 class OilWaterMaterialLawT,
51 class GasWaterMaterialLawT,
52 class ParamsT = EclTwoPhaseMaterialParams<TraitsT,
53 typename GasOilMaterialLawT::Params,
54 typename OilWaterMaterialLawT::Params,
55 typename GasWaterMaterialLawT::Params> >
59 using GasOilMaterialLaw = GasOilMaterialLawT;
60 using OilWaterMaterialLaw = OilWaterMaterialLawT;
61 using GasWaterMaterialLaw = GasWaterMaterialLawT;
64 static_assert(TraitsT::numPhases == 3,
65 "The number of phases considered by this capillary pressure "
66 "law is always three!");
67 static_assert(GasOilMaterialLaw::numPhases == 2,
68 "The number of phases considered by the gas-oil capillary "
69 "pressure law must be two!");
70 static_assert(OilWaterMaterialLaw::numPhases == 2,
71 "The number of phases considered by the oil-water capillary "
72 "pressure law must be two!");
73 static_assert(GasWaterMaterialLaw::numPhases == 2,
74 "The number of phases considered by the gas-water capillary "
75 "pressure law must be two!");
76 static_assert(std::is_same<
typename GasOilMaterialLaw::Scalar,
77 typename OilWaterMaterialLaw::Scalar>::value,
78 "The two two-phase capillary pressure laws must use the same "
79 "type of floating point values.");
81 using Traits = TraitsT;
82 using Params = ParamsT;
83 using Scalar =
typename Traits::Scalar;
85 static constexpr int numPhases = 3;
86 static constexpr int waterPhaseIdx = Traits::wettingPhaseIdx;
87 static constexpr int oilPhaseIdx = Traits::nonWettingPhaseIdx;
88 static constexpr int gasPhaseIdx = Traits::gasPhaseIdx;
114 template <
class ContainerT,
class Flu
idState>
115 static Scalar relpermOilInOilGasSystem(
const Params& ,
116 const FluidState& ) {
117 throw std::logic_error {
118 "relpermOilInOilGasSystem() is specific to three phases"
121 template <
class ContainerT,
class Flu
idState>
122 static Scalar relpermOilInOilWaterSystem(
const Params& ,
123 const FluidState& ) {
124 throw std::logic_error {
125 "relpermOilInOilWaterSystem() is specific to three phases"
144 template <
class ContainerT,
class Flu
idState>
146 const Params& params,
147 const FluidState& fluidState)
149 OPM_TIMEFUNCTION_LOCAL();
150 using Evaluation =
typename std::remove_reference<
decltype(values[0])>::type;
152 switch (params.approach()) {
153 case EclTwoPhaseApproach::GasOil: {
154 const Evaluation& So =
155 decay<Evaluation>(fluidState.saturation(oilPhaseIdx));
157 values[oilPhaseIdx] = 0.0;
158 values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(), So);
162 case EclTwoPhaseApproach::OilWater: {
163 const Evaluation&
Sw =
164 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
166 values[waterPhaseIdx] = 0.0;
167 values[oilPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(),
Sw);
171 case EclTwoPhaseApproach::GasWater: {
172 const Evaluation&
Sw =
173 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
175 values[waterPhaseIdx] = 0.0;
176 values[gasPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatPcnw(params.gasWaterParams(),
Sw);
189 static void oilWaterHysteresisParams(Scalar& pcSwMdc,
191 const Params& params)
193 OPM_TIMEFUNCTION_LOCAL();
194 pcSwMdc = params.oilWaterParams().pcSwMdc();
195 krnSwMdc = params.oilWaterParams().krnSwMdc();
197 Valgrind::CheckDefined(pcSwMdc);
198 Valgrind::CheckDefined(krnSwMdc);
207 static void setOilWaterHysteresisParams(
const Scalar& pcSwMdc,
208 const Scalar& krnSwMdc,
211 OPM_TIMEFUNCTION_LOCAL();
212 constexpr const Scalar krwSw = 2.0;
213 params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
222 static void gasOilHysteresisParams(Scalar& pcSwMdc,
224 const Params& params)
226 pcSwMdc = params.gasOilParams().pcSwMdc();
227 krnSwMdc = params.gasOilParams().krnSwMdc();
229 Valgrind::CheckDefined(pcSwMdc);
230 Valgrind::CheckDefined(krnSwMdc);
239 static void setGasOilHysteresisParams(
const Scalar& pcSwMdc,
240 const Scalar& krnSwMdc,
243 constexpr const Scalar krwSw = 2.0;
244 params.gasOilParams().update(pcSwMdc, krwSw, krnSwMdc);
247 static Scalar trappedGasSaturation(
const Params& params){
248 if(params.approach() == EclTwoPhaseApproach::GasOil)
249 return params.gasOilParams().SnTrapped();
250 if(params.approach() == EclTwoPhaseApproach::GasWater)
251 return params.gasWaterParams().SnTrapped();
255 static Scalar trappedOilSaturation(
const Params& params){
256 if(params.approach() == EclTwoPhaseApproach::GasOil)
257 return params.gasOilParams().SwTrapped();
258 if(params.approach() == EclTwoPhaseApproach::OilWater)
259 return params.oilWaterParams().SnTrapped();
271 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
272 static Evaluation
pcgn(
const Params& ,
275 throw std::logic_error(
"Not implemented: pcgn()");
287 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
288 static Evaluation
pcnw(
const Params& ,
291 throw std::logic_error(
"Not implemented: pcnw()");
297 template <
class ContainerT,
class Flu
idState>
302 throw std::logic_error(
"Not implemented: saturations()");
308 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
309 static Evaluation
Sg(
const Params& ,
312 throw std::logic_error(
"Not implemented: Sg()");
318 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
319 static Evaluation
Sn(
const Params& ,
322 throw std::logic_error(
"Not implemented: Sn()");
328 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
329 static Evaluation
Sw(
const Params& ,
332 throw std::logic_error(
"Not implemented: Sw()");
350 template <
class ContainerT,
class Flu
idState>
352 const Params& params,
353 const FluidState& fluidState)
355 OPM_TIMEFUNCTION_LOCAL();
356 using Evaluation =
typename std::remove_reference<
decltype(values[0])>::type;
358 switch (params.approach()) {
359 case EclTwoPhaseApproach::GasOil: {
360 const Evaluation& So =
361 decay<Evaluation>(fluidState.saturation(oilPhaseIdx));
363 values[oilPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), So);
364 values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), So);
368 case EclTwoPhaseApproach::OilWater: {
369 const Evaluation&
Sw =
370 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
372 values[waterPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(),
Sw);
373 values[oilPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(),
Sw);
377 case EclTwoPhaseApproach::GasWater: {
378 const Evaluation&
Sw =
379 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
381 values[waterPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatKrw(params.gasWaterParams(),
Sw);
382 values[gasPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatKrn(params.gasWaterParams(),
Sw);
392 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
393 static Evaluation
krg(
const Params& ,
396 throw std::logic_error(
"Not implemented: krg()");
402 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
403 static Evaluation
krw(
const Params& ,
406 throw std::logic_error(
"Not implemented: krw()");
412 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
413 static Evaluation
krn(
const Params& ,
416 throw std::logic_error(
"Not implemented: krn()");
427 template <
class Flu
idState>
430 OPM_TIMEFUNCTION_LOCAL();
431 switch (params.approach()) {
432 case EclTwoPhaseApproach::GasOil: {
433 Scalar So = scalarValue(fluidState.saturation(oilPhaseIdx));
435 return params.gasOilParams().update(So, So, So);
439 case EclTwoPhaseApproach::OilWater: {
440 Scalar
Sw = scalarValue(fluidState.saturation(waterPhaseIdx));
442 return params.oilWaterParams().update(
Sw,
Sw,
Sw);
446 case EclTwoPhaseApproach::GasWater: {
447 Scalar
Sw = scalarValue(fluidState.saturation(waterPhaseIdx));
449 return params.gasWaterParams().update(
Sw,
Sw,
Sw);
Implementation for the parameters required by the material law for two-phase simulations.
Some templates to wrap the valgrind client request macros.
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Definition EclTwoPhaseMaterial.hpp:57
static constexpr bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition EclTwoPhaseMaterial.hpp:108
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition EclTwoPhaseMaterial.hpp:351
static constexpr bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition EclTwoPhaseMaterial.hpp:92
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition EclTwoPhaseMaterial.hpp:329
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition EclTwoPhaseMaterial.hpp:272
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition EclTwoPhaseMaterial.hpp:298
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition EclTwoPhaseMaterial.hpp:413
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition EclTwoPhaseMaterial.hpp:319
static constexpr bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition EclTwoPhaseMaterial.hpp:100
static constexpr bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition EclTwoPhaseMaterial.hpp:104
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition EclTwoPhaseMaterial.hpp:393
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition EclTwoPhaseMaterial.hpp:309
static bool updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclTwoPhaseMaterial.hpp:428
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
static constexpr bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition EclTwoPhaseMaterial.hpp:112
static Evaluation pcnw(const Params &, const FluidState &)
Capillary pressure between the non-wetting liquid (i.e., oil) and the wetting liquid (i....
Definition EclTwoPhaseMaterial.hpp:288
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition EclTwoPhaseMaterial.hpp:403
static constexpr bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition EclTwoPhaseMaterial.hpp:96
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30