25 #ifndef OPM_ECL_HYSTERESIS_TWO_PHASE_LAW_HPP
26 #define OPM_ECL_HYSTERESIS_TWO_PHASE_LAW_HPP
36 template <
class EffectiveLawT,
37 class ParamsT = EclHysteresisTwoPhaseLawParams<EffectiveLawT> >
44 typedef typename EffectiveLaw::Traits
Traits;
46 typedef typename EffectiveLaw::Scalar
Scalar;
52 static const int numPhases = EffectiveLaw::numPhases;
53 static_assert(numPhases == 2,
54 "The endpoint scaling applies to the nested twophase laws, not to "
55 "the threephase one!");
59 static const bool implementsTwoPhaseApi =
true;
61 static_assert(EffectiveLaw::implementsTwoPhaseApi,
62 "The material laws put into EclEpsTwoPhaseLaw must implement the "
63 "two-phase material law API!");
67 static const bool implementsTwoPhaseSatApi =
true;
69 static_assert(EffectiveLaw::implementsTwoPhaseSatApi,
70 "The material laws put into EclEpsTwoPhaseLaw must implement the "
71 "two-phase material law saturation API!");
75 static const bool isSaturationDependent =
true;
79 static const bool isPressureDependent =
false;
83 static const bool isTemperatureDependent =
false;
87 static const bool isCompositionDependent =
false;
99 template <
class Container,
class Flu
idState>
100 static void capillaryPressures(Container& ,
104 OPM_THROW(NotImplemented,
105 "The capillaryPressures(fs) method is not yet implemented");
118 template <
class Container,
class Flu
idState>
119 static void relativePermeabilities(Container& ,
123 OPM_THROW(NotImplemented,
124 "The pcnw(fs) method is not yet implemented");
138 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
139 static Evaluation pcnw(
const Params& ,
142 OPM_THROW(NotImplemented,
143 "The pcnw(fs) method is not yet implemented");
146 template <
class Evaluation>
147 static Evaluation twoPhaseSatPcnw(
const Params ¶ms,
const Evaluation& Sw)
150 return EffectiveLaw::twoPhaseSatPcnw(params.drainageParams(), Sw);
168 template <
class Container,
class Flu
idState>
169 static void saturations(Container& ,
173 OPM_THROW(NotImplemented,
174 "The saturations(fs) method is not yet implemented");
181 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
182 static Evaluation Sw(
const Params& ,
185 OPM_THROW(NotImplemented,
186 "The Sw(fs) method is not yet implemented");
189 template <
class Evaluation>
190 static Evaluation twoPhaseSatSw(
const Params& ,
193 OPM_THROW(NotImplemented,
194 "The twoPhaseSatSw(pc) method is not yet implemented");
201 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
202 static Evaluation Sn(
const Params& ,
205 OPM_THROW(NotImplemented,
206 "The Sn(pc) method is not yet implemented");
209 template <
class Evaluation>
210 static Evaluation twoPhaseSatSn(
const Params& ,
213 OPM_THROW(NotImplemented,
214 "The twoPhaseSatSn(pc) method is not yet implemented");
226 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
227 static Evaluation krw(
const Params& ,
230 OPM_THROW(NotImplemented,
231 "The krw(fs) method is not yet implemented");
234 template <
class Evaluation>
235 static Evaluation twoPhaseSatKrw(
const Params ¶ms,
const Evaluation& Sw)
239 if (!params.config().enableHysteresis() || params.config().krHysteresisModel() < 0)
240 return EffectiveLaw::twoPhaseSatKrw(params.drainageParams(), Sw);
244 if (Sw <= params.krwSwMdc())
245 return EffectiveLaw::twoPhaseSatKrw(params.drainageParams(), Sw);
247 return EffectiveLaw::twoPhaseSatKrw(params.imbibitionParams(),
248 Sw + params.deltaSwImbKrw());
254 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
255 static Evaluation krn(
const Params& ,
258 OPM_THROW(NotImplemented,
259 "The krn(fs) method is not yet implemented");
262 template <
class Evaluation>
263 static Evaluation twoPhaseSatKrn(
const Params ¶ms,
const Evaluation& Sw)
266 if (!params.config().enableHysteresis() || params.config().krHysteresisModel() < 0)
267 return EffectiveLaw::twoPhaseSatKrn(params.drainageParams(), Sw);
271 if (Sw <= params.krnSwMdc())
272 return EffectiveLaw::twoPhaseSatKrn(params.drainageParams(), Sw);
274 return EffectiveLaw::twoPhaseSatKrn(params.imbibitionParams(),
275 Sw + params.deltaSwImbKrn());
EffectiveLaw::Traits Traits
Definition: EclHysteresisTwoPhaseLaw.hpp:44
Definition: Air_Mesitylene.hpp:31
EffectiveLawT EffectiveLaw
Definition: EclHysteresisTwoPhaseLaw.hpp:41
Definition: EclHysteresisTwoPhaseLaw.hpp:48
EffectiveLaw::Params EffectiveLawParams
Definition: EclHysteresisTwoPhaseLaw.hpp:42
ParamsT Params
Definition: EclHysteresisTwoPhaseLaw.hpp:45
A default implementation of the parameters for the material law which implements the ECL relative per...
EffectiveLaw::Scalar Scalar
Definition: EclHysteresisTwoPhaseLaw.hpp:46
Definition: EclHysteresisTwoPhaseLaw.hpp:49
This material law implements the hysteresis model of the ECL file format.
Definition: EclHysteresisTwoPhaseLaw.hpp:38
static const int numPhases
The number of fluid phases.
Definition: EclHysteresisTwoPhaseLaw.hpp:52