27 #ifndef OPM_ECL_DEFAULT_MATERIAL_HPP
28 #define OPM_ECL_DEFAULT_MATERIAL_HPP
52 template <
class TraitsT,
53 class GasOilMaterialLawT,
54 class OilWaterMaterialLawT,
55 class ParamsT = EclDefaultMaterialParams<TraitsT,
56 typename GasOilMaterialLawT::Params,
57 typename OilWaterMaterialLawT::Params> >
61 typedef GasOilMaterialLawT GasOilMaterialLaw;
62 typedef OilWaterMaterialLawT OilWaterMaterialLaw;
65 static_assert(TraitsT::numPhases == 3,
66 "The number of phases considered by this capillary pressure "
67 "law is always three!");
68 static_assert(GasOilMaterialLaw::numPhases == 2,
69 "The number of phases considered by the gas-oil capillary "
70 "pressure law must be two!");
71 static_assert(OilWaterMaterialLaw::numPhases == 2,
72 "The number of phases considered by the oil-water capillary "
73 "pressure law must be two!");
74 static_assert(std::is_same<
typename GasOilMaterialLaw::Scalar,
75 typename OilWaterMaterialLaw::Scalar>::value,
76 "The two two-phase capillary pressure laws must use the same "
77 "type of floating point values.");
79 static_assert(GasOilMaterialLaw::implementsTwoPhaseSatApi,
80 "The gas-oil material law must implement the two-phase saturation "
81 "only API to for the default Ecl capillary pressure law!");
82 static_assert(OilWaterMaterialLaw::implementsTwoPhaseSatApi,
83 "The oil-water material law must implement the two-phase saturation "
84 "only API to for the default Ecl capillary pressure law!");
86 typedef TraitsT Traits;
87 typedef ParamsT Params;
88 typedef typename Traits::Scalar Scalar;
90 static const int numPhases = 3;
91 static const int waterPhaseIdx = Traits::wettingPhaseIdx;
92 static const int oilPhaseIdx = Traits::nonWettingPhaseIdx;
93 static const int gasPhaseIdx = Traits::gasPhaseIdx;
133 template <
class ContainerT,
class Flu
idState>
135 const Params& params,
136 const FluidState& state)
138 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
139 values[gasPhaseIdx] = pcgn<FluidState, Evaluation>(params, state);
140 values[oilPhaseIdx] = 0;
141 values[waterPhaseIdx] = - pcnw<FluidState, Evaluation>(params, state);
143 Valgrind::CheckDefined(values[gasPhaseIdx]);
144 Valgrind::CheckDefined(values[oilPhaseIdx]);
145 Valgrind::CheckDefined(values[waterPhaseIdx]);
154 static void oilWaterHysteresisParams(Scalar& pcSwMdc,
156 const Params& params)
158 pcSwMdc = params.oilWaterParams().pcSwMdc();
159 krnSwMdc = params.oilWaterParams().krnSwMdc();
161 Valgrind::CheckDefined(pcSwMdc);
162 Valgrind::CheckDefined(krnSwMdc);
171 static void setOilWaterHysteresisParams(
const Scalar& pcSwMdc,
172 const Scalar& krnSwMdc,
175 const double krwSw = 2.0;
176 params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
185 static void gasOilHysteresisParams(Scalar& pcSwMdc,
187 const Params& params)
189 const auto Swco = params.Swl();
193 pcSwMdc = std::min(params.gasOilParams().pcSwMdc() + Swco, Scalar{2.0});
194 krnSwMdc = std::min(params.gasOilParams().krnSwMdc() + Swco, Scalar{2.0});
196 Valgrind::CheckDefined(pcSwMdc);
197 Valgrind::CheckDefined(krnSwMdc);
206 static void setGasOilHysteresisParams(
const Scalar& pcSwMdc,
207 const Scalar& krnSwMdc,
211 const auto Swco = params.Swl();
212 const double krwSw = 2.0;
213 params.gasOilParams().update(pcSwMdc - Swco, krwSw, krnSwMdc - Swco);
225 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
226 static Evaluation
pcgn(
const Params& params,
227 const FluidState& fs)
230 const auto Sw = 1.0 - params.Swl() - decay<Evaluation>(fs.saturation(gasPhaseIdx));
231 return GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(),
Sw);
243 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
244 static Evaluation
pcnw(
const Params& params,
245 const FluidState& fs)
247 const auto Sw = decay<Evaluation>(fs.saturation(waterPhaseIdx));
248 return OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(),
Sw);
254 template <
class ContainerT,
class Flu
idState>
259 throw std::logic_error(
"Not implemented: saturations()");
265 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
266 static Evaluation
Sg(
const Params& ,
269 throw std::logic_error(
"Not implemented: Sg()");
275 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
276 static Evaluation
Sn(
const Params& ,
279 throw std::logic_error(
"Not implemented: Sn()");
285 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
286 static Evaluation
Sw(
const Params& ,
289 throw std::logic_error(
"Not implemented: Sw()");
307 template <
class ContainerT,
class Flu
idState>
309 const Params& params,
310 const FluidState& fluidState)
312 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
314 values[waterPhaseIdx] = krw<FluidState, Evaluation>(params, fluidState);
315 values[oilPhaseIdx] = krn<FluidState, Evaluation>(params, fluidState);
316 values[gasPhaseIdx] = krg<FluidState, Evaluation>(params, fluidState);
322 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
323 static Evaluation
krg(
const Params& params,
324 const FluidState& fluidState)
327 const Evaluation
Sw = 1.0 - params.Swl() - decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
328 return GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(),
Sw);
334 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
335 static Evaluation
krw(
const Params& params,
336 const FluidState& fluidState)
338 const Evaluation
Sw = decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
339 return OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(),
Sw);
345 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
346 static Evaluation
krn(
const Params& params,
347 const FluidState& fluidState)
349 const Scalar Swco = params.Swl();
351 const Evaluation
Sw =
352 max(Evaluation(Swco),
353 decay<Evaluation>(fluidState.saturation(waterPhaseIdx)));
355 const Evaluation
Sg = decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
357 const Evaluation Sw_ow =
Sg +
Sw;
358 const Evaluation kro_ow = relpermOilInOilWaterSystem<Evaluation>(params, fluidState);
359 const Evaluation kro_go = relpermOilInOilGasSystem<Evaluation>(params, fluidState);
364 const Scalar epsilon = 1e-5;
365 if (scalarValue(Sw_ow) - Swco < epsilon) {
366 const Evaluation kro2 = (kro_ow + kro_go)/2;
367 if (scalarValue(Sw_ow) - Swco > epsilon/2) {
368 const Evaluation kro1 = (
Sg*kro_go + (
Sw - Swco)*kro_ow)/(Sw_ow - Swco);
369 const Evaluation alpha = (epsilon - (Sw_ow - Swco))/(epsilon/2);
371 return kro2*alpha + kro1*(1 - alpha);
377 return (
Sg*kro_go + (
Sw - Swco)*kro_ow) / (Sw_ow - Swco);
383 template <
class Evaluation,
class Flu
idState>
385 const FluidState& fluidState)
387 const Evaluation
Sw =
388 max(Evaluation{ params.Swl() },
389 decay<Evaluation>(fluidState.saturation(waterPhaseIdx)));
391 const Evaluation
Sg = decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
392 const Evaluation So_go = 1.0 - (
Sg +
Sw);
394 return GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), So_go);
400 template <
class Evaluation,
class Flu
idState>
402 const FluidState& fluidState)
404 const Evaluation
Sw =
405 max(Evaluation{ params.Swl() },
406 decay<Evaluation>(fluidState.saturation(waterPhaseIdx)));
408 const Evaluation
Sg = decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
409 const Evaluation Sw_ow =
Sg +
Sw;
411 return OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), Sw_ow);
421 template <
class Flu
idState>
424 const Scalar Swco = params.Swl();
426 const Scalar
Sw = clampSaturation(fluidState, waterPhaseIdx);
427 const Scalar So = clampSaturation(fluidState, oilPhaseIdx);
428 const Scalar
Sg = clampSaturation(fluidState, gasPhaseIdx);
430 if (params.inconsistentHysteresisUpdate()) {
440 params.oilWaterParams().update( 1.0 - So,
444 params.gasOilParams().update( 1.0 - Swco -
Sg,
449 const Scalar Sw_ow =
Sg + std::max(Swco,
Sw);
450 const Scalar So_go = 1.0 - Sw_ow;
452 params.oilWaterParams().update(
Sw,
456 params.gasOilParams().update( 1.0 - Swco -
Sg,
462 template <
class Flu
idState>
463 static Scalar clampSaturation(
const FluidState& fluidState,
const int phaseIndex)
465 const auto sat = scalarValue(fluidState.saturation(phaseIndex));
466 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:59
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclDefaultMaterial.hpp:101
static Evaluation relpermOilInOilWaterSystem(const Params ¶ms, const FluidState &fluidState)
The relative permeability of oil in oil/water system.
Definition: EclDefaultMaterial.hpp:401
static Evaluation krg(const Params ¶ms, const FluidState &fluidState)
The relative permeability of the gas phase.
Definition: EclDefaultMaterial.hpp:323
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:226
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: EclDefaultMaterial.hpp:109
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition: EclDefaultMaterial.hpp:266
static Evaluation relpermOilInOilGasSystem(const Params ¶ms, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition: EclDefaultMaterial.hpp:384
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 saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: EclDefaultMaterial.hpp:255
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: EclDefaultMaterial.hpp:117
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition: EclDefaultMaterial.hpp:276
static Evaluation krn(const Params ¶ms, const FluidState &fluidState)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition: EclDefaultMaterial.hpp:346
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:244
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclDefaultMaterial.hpp:97
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclDefaultMaterial.hpp:105
static void relativePermeabilities(ContainerT &values, const Params ¶ms, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclDefaultMaterial.hpp:308
static Evaluation krw(const Params ¶ms, const FluidState &fluidState)
The relative permeability of the wetting phase.
Definition: EclDefaultMaterial.hpp:335
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition: EclDefaultMaterial.hpp:286
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclDefaultMaterial.hpp:113
static void updateHysteresis(Params ¶ms, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclDefaultMaterial.hpp:422