27#ifndef OPM_ECL_DEFAULT_MATERIAL_HPP
28#define OPM_ECL_DEFAULT_MATERIAL_HPP
30#include <opm/common/TimingMacros.hpp>
55template <
class TraitsT,
56 class GasOilMaterialLawT,
57 class OilWaterMaterialLawT,
58 class ParamsT = EclDefaultMaterialParams<TraitsT,
59 typename GasOilMaterialLawT::Params,
60 typename OilWaterMaterialLawT::Params> >
64 using GasOilMaterialLaw = GasOilMaterialLawT;
65 using OilWaterMaterialLaw = OilWaterMaterialLawT;
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(std::is_same<
typename GasOilMaterialLaw::Scalar,
78 typename OilWaterMaterialLaw::Scalar>::value,
79 "The two two-phase capillary pressure laws must use the same "
80 "type of floating point values.");
82 static_assert(GasOilMaterialLaw::implementsTwoPhaseSatApi,
83 "The gas-oil material law must implement the two-phase saturation "
84 "only API to for the default Ecl capillary pressure law!");
85 static_assert(OilWaterMaterialLaw::implementsTwoPhaseSatApi,
86 "The oil-water material law must implement the two-phase saturation "
87 "only API to for the default Ecl capillary pressure law!");
89 using Traits = TraitsT;
90 using Params = ParamsT;
91 using Scalar =
typename Traits::Scalar;
93 static constexpr int numPhases = 3;
94 static constexpr int waterPhaseIdx = Traits::wettingPhaseIdx;
95 static constexpr int oilPhaseIdx = Traits::nonWettingPhaseIdx;
96 static constexpr int gasPhaseIdx = Traits::gasPhaseIdx;
136 template <
class ContainerT,
class Flu
idState>
138 const Params& params,
139 const FluidState& state)
141 OPM_TIMEFUNCTION_LOCAL();
142 using Evaluation =
typename std::remove_reference<
decltype(values[0])>::type;
143 values[gasPhaseIdx] = pcgn<FluidState, Evaluation>(params, state);
144 values[oilPhaseIdx] = 0;
145 values[waterPhaseIdx] = - pcnw<FluidState, Evaluation>(params, state);
147 Valgrind::CheckDefined(values[gasPhaseIdx]);
148 Valgrind::CheckDefined(values[oilPhaseIdx]);
149 Valgrind::CheckDefined(values[waterPhaseIdx]);
159 static void oilWaterHysteresisParams(Scalar& soMax,
162 const Params& params)
164 soMax = 1.0 - params.oilWaterParams().krnSwMdc();
165 swMax = params.oilWaterParams().krwSwMdc();
166 swMin = params.oilWaterParams().pcSwMdc();
167 Valgrind::CheckDefined(soMax);
168 Valgrind::CheckDefined(swMax);
169 Valgrind::CheckDefined(swMin);
179 static void setOilWaterHysteresisParams(
const Scalar& soMax,
184 params.oilWaterParams().update(swMin, swMax, 1.0 - soMax);
194 static void gasOilHysteresisParams(Scalar& sgMax,
197 const Params& params)
199 const auto Swco = params.Swl();
200 sgMax = 1.0 - params.gasOilParams().krnSwMdc() - Swco;
201 shMax = params.gasOilParams().krwSwMdc();
202 soMin = params.gasOilParams().pcSwMdc();
204 Valgrind::CheckDefined(sgMax);
205 Valgrind::CheckDefined(shMax);
206 Valgrind::CheckDefined(soMin);
216 static void setGasOilHysteresisParams(
const Scalar& sgMax,
221 const auto Swco = params.Swl();
222 params.gasOilParams().update(soMin, shMax, 1.0 - sgMax - Swco);
225 static Scalar trappedGasSaturation(
const Params& params,
bool maximumTrapping)
227 const auto Swco = params.Swl();
228 return params.gasOilParams().SnTrapped(maximumTrapping) - Swco;
231 static Scalar trappedOilSaturation(
const Params& params,
bool maximumTrapping)
233 return params.oilWaterParams().SnTrapped(maximumTrapping) + params.gasOilParams().SwTrapped();
236 static Scalar trappedWaterSaturation(
const Params& params)
238 return params.oilWaterParams().SwTrapped();
241 static Scalar strandedGasSaturation(
const Params& params, Scalar
Sg, Scalar Kg)
243 const auto Swco = params.Swl();
244 return params.gasOilParams().SnStranded(
Sg, Kg) - Swco;
256 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
257 static Evaluation
pcgn(
const Params& params,
258 const FluidState& fs)
260 OPM_TIMEFUNCTION_LOCAL();
262 const auto Sw = 1.0 - params.Swl() - decay<Evaluation>(fs.saturation(gasPhaseIdx));
263 return GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(),
Sw);
275 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
276 static Evaluation
pcnw(
const Params& params,
277 const FluidState& fs)
279 OPM_TIMEFUNCTION_LOCAL();
280 const auto Sw = decay<Evaluation>(fs.saturation(waterPhaseIdx));
281 return OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(),
Sw);
287 template <
class ContainerT,
class Flu
idState>
292 throw std::logic_error(
"Not implemented: saturations()");
298 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
299 static Evaluation
Sg(
const Params& ,
302 throw std::logic_error(
"Not implemented: Sg()");
308 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
309 static Evaluation
Sn(
const Params& ,
312 throw std::logic_error(
"Not implemented: Sn()");
318 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
319 static Evaluation
Sw(
const Params& ,
322 throw std::logic_error(
"Not implemented: Sw()");
340 template <
class ContainerT,
class Flu
idState>
342 const Params& params,
343 const FluidState& fluidState)
345 OPM_TIMEFUNCTION_LOCAL();
346 using Evaluation =
typename std::remove_reference<
decltype(values[0])>::type;
348 values[waterPhaseIdx] = krw<FluidState, Evaluation>(params, fluidState);
349 values[oilPhaseIdx] = krn<FluidState, Evaluation>(params, fluidState);
350 values[gasPhaseIdx] = krg<FluidState, Evaluation>(params, fluidState);
356 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
357 static Evaluation
krg(
const Params& params,
358 const FluidState& fluidState)
360 OPM_TIMEFUNCTION_LOCAL();
362 const Evaluation sw = 1.0 - params.Swl() - decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
363 return GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), sw);
369 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
370 static Evaluation
krw(
const Params& params,
371 const FluidState& fluidState)
373 OPM_TIMEFUNCTION_LOCAL();
374 const Evaluation sw = decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
375 return OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), sw);
381 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
382 static Evaluation
krn(
const Params& params,
383 const FluidState& fluidState)
385 OPM_TIMEFUNCTION_LOCAL();
386 const Scalar Swco = params.Swl();
388 const Evaluation sw =
389 max(Evaluation(Swco),
390 decay<Evaluation>(fluidState.saturation(waterPhaseIdx)));
392 const Evaluation sg = decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
394 const Evaluation Sw_ow = sg + sw;
395 const Evaluation kro_ow = relpermOilInOilWaterSystem<Evaluation>(params, fluidState);
396 const Evaluation kro_go = relpermOilInOilGasSystem<Evaluation>(params, fluidState);
401 constexpr const Scalar epsilon = 1e-5;
402 if (scalarValue(Sw_ow) - Swco < epsilon) {
403 const Evaluation kro2 = (kro_ow + kro_go)/2;
404 if (scalarValue(Sw_ow) - Swco > epsilon/2) {
405 const Evaluation kro1 = (sg * kro_go + (sw - Swco) * kro_ow) / (Sw_ow - Swco);
406 const Evaluation alpha = (epsilon - (Sw_ow - Swco)) / (epsilon / 2);
408 return kro2 * alpha + kro1 * (1 - alpha);
414 return (sg * kro_go + (sw - Swco) * kro_ow) / (Sw_ow - Swco);
420 template <
class Evaluation,
class Flu
idState>
422 const FluidState& fluidState)
424 OPM_TIMEFUNCTION_LOCAL();
425 const Evaluation sw =
426 max(Evaluation{ params.Swl() },
427 decay<Evaluation>(fluidState.saturation(waterPhaseIdx)));
429 const Evaluation sg = decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
430 const Evaluation So_go = 1.0 - (sg + sw);
432 return GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), So_go);
438 template <
class Evaluation,
class Flu
idState>
440 const FluidState& fluidState)
442 OPM_TIMEFUNCTION_LOCAL();
443 const Evaluation sw =
444 max(Evaluation{ params.Swl() },
445 decay<Evaluation>(fluidState.saturation(waterPhaseIdx)));
447 const Evaluation sg = decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
448 const Evaluation Sw_ow = sg + sw;
450 return OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), Sw_ow);
460 template <
class Flu
idState>
463 const Scalar Swco = params.Swl();
464 const Scalar sw = clampSaturation(fluidState, waterPhaseIdx);
465 const Scalar So = clampSaturation(fluidState, oilPhaseIdx);
466 const Scalar sg = clampSaturation(fluidState, gasPhaseIdx);
467 bool owChanged = params.oilWaterParams().update(sw, sw, 1 - So);
468 bool gochanged = params.gasOilParams().update( So,
471 return owChanged || gochanged;
474 template <
class Flu
idState>
475 static Scalar clampSaturation(
const FluidState& fluidState,
const int phaseIndex)
477 OPM_TIMEFUNCTION_LOCAL();
478 const auto sat = scalarValue(fluidState.saturation(phaseIndex));
479 return std::clamp(sat, Scalar{0.0}, Scalar{1.0});
Default implementation for the parameters required by the default three-phase capillary pressure mode...
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:62
static Evaluation relpermOilInOilWaterSystem(const Params ¶ms, const FluidState &fluidState)
The relative permeability of oil in oil/water system.
Definition EclDefaultMaterial.hpp:439
static Evaluation krg(const Params ¶ms, const FluidState &fluidState)
The relative permeability of the gas phase.
Definition EclDefaultMaterial.hpp:357
static Evaluation pcgn(const Params ¶ms, const FluidState &fs)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition EclDefaultMaterial.hpp:257
static constexpr bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition EclDefaultMaterial.hpp:120
static constexpr bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition EclDefaultMaterial.hpp:116
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition EclDefaultMaterial.hpp:299
static Evaluation relpermOilInOilGasSystem(const Params ¶ms, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition EclDefaultMaterial.hpp:421
static constexpr bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition EclDefaultMaterial.hpp:100
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:137
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition EclDefaultMaterial.hpp:288
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition EclDefaultMaterial.hpp:309
static Evaluation krn(const Params ¶ms, const FluidState &fluidState)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition EclDefaultMaterial.hpp:382
static Evaluation pcnw(const Params ¶ms, const FluidState &fs)
Capillary pressure between the non-wetting liquid (i.e., oil) and the wetting liquid (i....
Definition EclDefaultMaterial.hpp:276
static constexpr bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition EclDefaultMaterial.hpp:108
static bool updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclDefaultMaterial.hpp:461
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition EclDefaultMaterial.hpp:341
static Evaluation krw(const Params ¶ms, const FluidState &fluidState)
The relative permeability of the wetting phase.
Definition EclDefaultMaterial.hpp:370
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition EclDefaultMaterial.hpp:319
static constexpr bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition EclDefaultMaterial.hpp:112
static constexpr bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition EclDefaultMaterial.hpp:104
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30