45 using EffectiveLaw = EffectiveLawT;
46 using EffectiveLawParams =
typename EffectiveLaw::Params;
48 using Traits =
typename EffectiveLaw::Traits;
49 using Params = ParamsT;
50 using Scalar =
typename EffectiveLaw::Scalar;
52 enum { wettingPhaseIdx = Traits::wettingPhaseIdx };
53 enum { nonWettingPhaseIdx = Traits::nonWettingPhaseIdx };
56 static constexpr int numPhases = EffectiveLaw::numPhases;
58 "The endpoint scaling applies to the nested twophase laws, not to "
59 "the threephase one!");
65 static_assert(EffectiveLaw::implementsTwoPhaseApi,
66 "The material laws put into EclEpsTwoPhaseLaw must implement the "
67 "two-phase material law API!");
73 static_assert(EffectiveLaw::implementsTwoPhaseSatApi,
74 "The material laws put into EclEpsTwoPhaseLaw must implement the "
75 "two-phase material law saturation API!");
103 template <
class Container,
class Flu
idState>
108 throw std::invalid_argument(
"The capillaryPressures(fs) method is not yet implemented");
121 template <
class Container,
class Flu
idState>
126 throw std::invalid_argument(
"The pcnw(fs) method is not yet implemented");
140 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
141 static Evaluation
pcnw(
const Params& ,
144 throw std::invalid_argument(
"The pcnw(fs) method is not yet implemented");
147 template <
class Evaluation>
148 static Evaluation twoPhaseSatPcnw(
const Params& params,
const Evaluation&
Sw)
150 OPM_TIMEFUNCTION_LOCAL();
152 if (!params.config().enableHysteresis() || params.config().pcHysteresisModel() < 0)
153 return EffectiveLaw::twoPhaseSatPcnw(params.drainageParams(),
Sw);
156 if (params.initialImb()) {
157 if (
Sw >= params.pcSwMic()) {
158 return EffectiveLaw::twoPhaseSatPcnw(params.imbibitionParams(),
Sw);
161 const Evaluation& F = (1.0/(params.pcSwMic()-
Sw+params.curvatureCapPrs())-1.0/params.curvatureCapPrs())
162 / (1.0/(params.pcSwMic()-params.Swcrd()+params.curvatureCapPrs())-1.0/params.curvatureCapPrs());
164 const Evaluation& Pcd = EffectiveLaw::twoPhaseSatPcnw(params.drainageParams(),
Sw);
165 const Evaluation& Pci = EffectiveLaw::twoPhaseSatPcnw(params.imbibitionParams(),
Sw);
166 const Evaluation& pc_Killough = Pci+F*(Pcd-Pci);
173 if (
Sw <= params.pcSwMdc())
174 return EffectiveLaw::twoPhaseSatPcnw(params.drainageParams(),
Sw);
177 Scalar Swma = 1.0-params.Sncrt();
179 const Evaluation& Pci = EffectiveLaw::twoPhaseSatPcnw(params.imbibitionParams(),
Sw);
183 Scalar pciwght = params.pcWght();
185 const Evaluation& SwEff =
Sw;
186 const Evaluation& Pci = pciwght*EffectiveLaw::twoPhaseSatPcnw(params.imbibitionParams(), SwEff);
188 const Evaluation& Pcd = EffectiveLaw::twoPhaseSatPcnw(params.drainageParams(),
Sw);
193 const Evaluation& F = (1.0/(
Sw-params.pcSwMdc()+params.curvatureCapPrs())-1.0/params.curvatureCapPrs())
194 / (1.0/(Swma-params.pcSwMdc()+params.curvatureCapPrs())-1.0/params.curvatureCapPrs());
196 const Evaluation& pc_Killough = Pcd+F*(Pci-Pcd);
207 template <
class Container,
class Flu
idState>
212 throw std::invalid_argument(
"The saturations(fs) method is not yet implemented");
219 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
220 static Evaluation
Sw(
const Params& ,
223 throw std::invalid_argument(
"The Sw(fs) method is not yet implemented");
226 template <
class Evaluation>
227 static Evaluation twoPhaseSatSw(
const Params& ,
230 throw std::invalid_argument(
"The twoPhaseSatSw(pc) method is not yet implemented");
237 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
238 static Evaluation
Sn(
const Params& ,
241 throw std::invalid_argument(
"The Sn(pc) method is not yet implemented");
244 template <
class Evaluation>
245 static Evaluation twoPhaseSatSn(
const Params& ,
248 throw std::invalid_argument(
"The twoPhaseSatSn(pc) method is not yet implemented");
260 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
261 static Evaluation
krw(
const Params& ,
264 throw std::invalid_argument(
"The krw(fs) method is not yet implemented");
267 template <
class Evaluation>
268 static Evaluation twoPhaseSatKrw(
const Params& params,
const Evaluation&
Sw)
271 OPM_TIMEFUNCTION_LOCAL();
273 if (!params.config().enableHysteresis() || params.config().krHysteresisModel() < 0)
274 return EffectiveLaw::twoPhaseSatKrw(params.drainageParams(),
Sw);
276 if (params.config().krHysteresisModel() == 0 || params.config().krHysteresisModel() == 2)
278 return EffectiveLaw::twoPhaseSatKrw(params.drainageParams(),
Sw);
281 if (params.config().krHysteresisModel() == 1 || params.config().krHysteresisModel() == 3)
282 return EffectiveLaw::twoPhaseSatKrw(params.imbibitionParams(),
Sw);
284 if (
Sw >= params.krwSwMdc())
285 return EffectiveLaw::twoPhaseSatKrw(params.drainageParams(),
Sw);
288 assert(params.config().krHysteresisModel() == 4);
289 Evaluation Snorm = params.Swcri()+(
Sw-params.Swcrt())*(params.Swmaxd()-params.Swcri())/(params.Swhy()-params.Swcrt());
290 return params.krwWght()*EffectiveLaw::twoPhaseSatKrw(params.imbibitionParams(), Snorm);
296 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
297 static Evaluation
krn(
const Params& ,
300 throw std::invalid_argument(
"The krn(fs) method is not yet implemented");
303 template <
class Evaluation>
304 static Evaluation twoPhaseSatKrn(
const Params& params,
const Evaluation&
Sw)
306 OPM_TIMEFUNCTION_LOCAL();
309 if (params.gasOilHysteresisWAG()) {
312 if (
Sw <= params.krnSwMdc()+params.tolWAG() && params.nState()==1) {
313 Evaluation
krn = EffectiveLaw::twoPhaseSatKrn(params.drainageParams(),
Sw);
319 if (params.nState()==1) {
320 Evaluation Swf = params.computeSwf(
Sw);
321 Evaluation
krn = EffectiveLaw::twoPhaseSatKrn(params.drainageParams(), Swf);
326 if (
Sw <= params.krnSwDrainRevert()+params.tolWAG() ) {
327 Evaluation Krg = EffectiveLaw::twoPhaseSatKrn(params.drainageParams(),
Sw);
328 Evaluation KrgDrain2 = (Krg-params.krnDrainStart())*params.reductionDrain() + params.krnImbStart();
333 if (
Sw >= params.krnSwWAG()-params.tolWAG() ) {
334 Evaluation KrgImb2 = params.computeKrImbWAG(
Sw);
338 Evaluation Krg = EffectiveLaw::twoPhaseSatKrn(params.drainageParams(),
Sw);
339 Evaluation KrgDrainNxt = (Krg-params.krnDrainStartNxt())*params.reductionDrainNxt() + params.krnImbStartNxt();
345 if (!params.config().enableHysteresis() || params.config().krHysteresisModel() < 0)
346 return EffectiveLaw::twoPhaseSatKrn(params.drainageParams(),
Sw);
350 if (
Sw <= params.krnSwMdc()) {
351 return EffectiveLaw::twoPhaseSatKrn(params.drainageParams(),
Sw);
354 if (params.config().krHysteresisModel() <= 1) {
355 return EffectiveLaw::twoPhaseSatKrn(params.imbibitionParams(),
356 Sw + params.deltaSwImbKrn());
360 assert(params.config().krHysteresisModel() == 2 || params.config().krHysteresisModel() == 3 || params.config().krHysteresisModel() == 4);
361 Evaluation Snorm = params.Sncri()+(1.0-
Sw-params.Sncrt())*(params.Snmaxd()-params.Sncri())/(params.Snhy()-params.Sncrt());
362 return params.krnWght()*EffectiveLaw::twoPhaseSatKrn(params.imbibitionParams(),1.0-Snorm);