25 #ifndef OPM_ECL_STONE1_MATERIAL_HPP
26 #define OPM_ECL_STONE1_MATERIAL_HPP
33 #include <opm/common/Exceptions.hpp>
34 #include <opm/common/ErrorMacros.hpp>
54 template <
class TraitsT,
55 class GasOilMaterialLawT,
56 class OilWaterMaterialLawT,
57 class ParamsT = EclStone1MaterialParams<TraitsT,
58 typename GasOilMaterialLawT::Params,
59 typename OilWaterMaterialLawT::Params> >
67 static_assert(TraitsT::numPhases == 3,
68 "The number of phases considered by this capillary pressure "
69 "law is always three!");
70 static_assert(GasOilMaterialLaw::numPhases == 2,
71 "The number of phases considered by the gas-oil capillary "
72 "pressure law must be two!");
73 static_assert(OilWaterMaterialLaw::numPhases == 2,
74 "The number of phases considered by the oil-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 static_assert(GasOilMaterialLaw::implementsTwoPhaseSatApi,
82 "The gas-oil material law must implement the two-phase saturation "
83 "only API to for the default Ecl capillary pressure law!");
84 static_assert(OilWaterMaterialLaw::implementsTwoPhaseSatApi,
85 "The oil-water material law must implement the two-phase saturation "
86 "only API to for the default Ecl capillary pressure law!");
88 typedef TraitsT Traits;
89 typedef ParamsT Params;
90 typedef typename Traits::Scalar Scalar;
92 static const int numPhases = 3;
93 static const int waterPhaseIdx = Traits::wettingPhaseIdx;
94 static const int oilPhaseIdx = Traits::nonWettingPhaseIdx;
95 static const int gasPhaseIdx = Traits::gasPhaseIdx;
99 static const bool implementsTwoPhaseApi =
false;
103 static const bool implementsTwoPhaseSatApi =
false;
107 static const bool isSaturationDependent =
true;
111 static const bool isPressureDependent =
false;
115 static const bool isTemperatureDependent =
false;
119 static const bool isCompositionDependent =
false;
135 template <
class ContainerT,
class Flu
idState>
136 static void capillaryPressures(ContainerT &values,
137 const Params ¶ms,
138 const FluidState &state)
140 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
141 values[gasPhaseIdx] = pcgn<FluidState, Evaluation>(params, state);
142 values[oilPhaseIdx] = 0;
143 values[waterPhaseIdx] = - pcnw<FluidState, Evaluation>(params, state);
158 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
159 static Evaluation pcgn(
const Params ¶ms,
160 const FluidState &fs)
164 const auto& Sw = 1.0 - FsToolbox::template toLhs<Evaluation>(fs.saturation(gasPhaseIdx));
165 return GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(), Sw);
177 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
178 static Evaluation pcnw(
const Params ¶ms,
179 const FluidState &fs)
181 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
183 const auto& Sw = FsToolbox::template toLhs<Evaluation>(fs.saturation(waterPhaseIdx));
185 const auto& result = OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(), Sw);
193 template <
class ContainerT,
class Flu
idState>
194 static void saturations(ContainerT& ,
198 OPM_THROW(std::logic_error,
"Not implemented: saturations()");
204 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
205 static Evaluation Sg(
const Params& ,
208 OPM_THROW(std::logic_error,
"Not implemented: Sg()");
214 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
215 static Evaluation Sn(
const Params& ,
218 OPM_THROW(std::logic_error,
"Not implemented: Sn()");
224 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
225 static Evaluation Sw(
const Params& ,
228 OPM_THROW(std::logic_error,
"Not implemented: Sw()");
246 template <
class ContainerT,
class Flu
idState>
247 static void relativePermeabilities(ContainerT &values,
248 const Params ¶ms,
249 const FluidState &fluidState)
251 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
253 values[waterPhaseIdx] = krw<FluidState, Evaluation>(params, fluidState);
254 values[oilPhaseIdx] = krn<FluidState, Evaluation>(params, fluidState);
255 values[gasPhaseIdx] = krg<FluidState, Evaluation>(params, fluidState);
261 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
262 static Evaluation krg(
const Params ¶ms,
263 const FluidState &fluidState)
265 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
267 const Evaluation& Sw = 1 - FsToolbox::template toLhs<Evaluation>(fluidState.saturation(gasPhaseIdx));
268 return GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), Sw);
274 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
275 static Evaluation krw(
const Params ¶ms,
276 const FluidState &fluidState)
278 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
280 const Evaluation& Sw = FsToolbox::template toLhs<Evaluation>(fluidState.saturation(waterPhaseIdx));
281 return OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), Sw);
287 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
288 static Evaluation krn(
const Params ¶ms,
289 const FluidState &fluidState)
291 typedef MathToolbox<Evaluation> Toolbox;
292 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
296 Scalar Swco = params.Swl();
298 Scalar Sowcr = params.Sowcr();
299 Scalar Sogcr = params.Sogcr();
300 Scalar Som =
std::min(Sowcr, Sogcr);
302 Scalar eta = params.eta();
304 const Evaluation& Sw = FsToolbox::template toLhs<Evaluation>(fluidState.saturation(waterPhaseIdx));
305 const Evaluation& So = FsToolbox::template toLhs<Evaluation>(fluidState.saturation(oilPhaseIdx));
306 const Evaluation& Sg = FsToolbox::template toLhs<Evaluation>(fluidState.saturation(gasPhaseIdx));
310 SSw = (Sw - Swco)/(1 - Swco - Som);
316 SSo = (So - Som)/(1 - Swco - Som);
320 Evaluation SSg = Sg/(1 - Swco - Som);
322 Scalar krocw = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), Swco);
323 Evaluation krow = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), Sw);
324 Evaluation krog = GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), 1 - Sg);
326 Evaluation beta =
Toolbox::pow(SSo/((1 - SSw)*(1 - SSg)), eta);
327 return beta*krow*krog/krocw;
337 template <
class Flu
idState>
338 static void updateHysteresis(Params ¶ms,
const FluidState &fluidState)
340 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
342 Scalar Sw = FsToolbox::value(fluidState.saturation(waterPhaseIdx));
343 Scalar Sg = FsToolbox::value(fluidState.saturation(gasPhaseIdx));
345 params.oilWaterParams().update(Sw, Sw, Sw);
346 params.gasOilParams().update(1 - Sg, 1 - Sg, 1 - Sg);
bool CheckDefined(const T &value OPM_UNUSED)
Make valgrind complain if any of the memory occupied by an object is undefined.
Definition: Valgrind.hpp:74
Definition: Air_Mesitylene.hpp:31
GasOilMaterialLawT GasOilMaterialLaw
Definition: EclStone1Material.hpp:63
OilWaterMaterialLawT OilWaterMaterialLaw
Definition: EclStone1Material.hpp:64
Some templates to wrap the valgrind client request macros.
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition: EclStone1Material.hpp:60
Evaluation< Scalar, VarSetTag, numVars > min(const Evaluation< Scalar, VarSetTag, numVars > &x1, const Evaluation< Scalar, VarSetTag, numVars > &x2)
Definition: Math.hpp:61
Evaluation< Scalar, VarSetTag, numVars > pow(const Evaluation< Scalar, VarSetTag, numVars > &base, Scalar exp)
Definition: Math.hpp:312
Default implementation for the parameters required by the three-phase capillary pressure/relperm Ston...