27 #ifndef OPM_ECL_EPS_TWO_PHASE_LAW_HPP
28 #define OPM_ECL_EPS_TWO_PHASE_LAW_HPP
33 #include <opm/common/ErrorMacros.hpp>
34 #include <opm/common/Exceptions.hpp>
48 template <
class EffLawT,
49 class ParamsT = EclEpsTwoPhaseLawParams<EffLawT> >
52 typedef EffLawT EffLaw;
55 typedef typename EffLaw::Traits Traits;
56 typedef ParamsT Params;
57 typedef typename EffLaw::Scalar Scalar;
59 enum { wettingPhaseIdx = Traits::wettingPhaseIdx };
60 enum { nonWettingPhaseIdx = Traits::nonWettingPhaseIdx };
65 "The endpoint scaling applies to the nested twophase laws, not to "
66 "the threephase one!");
72 static_assert(EffLaw::implementsTwoPhaseApi,
73 "The material laws put into EclEpsTwoPhaseLaw must implement the "
74 "two-phase material law API!");
80 static_assert(EffLaw::implementsTwoPhaseSatApi,
81 "The material laws put into EclEpsTwoPhaseLaw must implement the "
82 "two-phase material law saturation API!");
110 template <
class Container,
class Flu
idState>
113 OPM_THROW(NotImplemented,
114 "The capillaryPressures(fs) method is not yet implemented");
127 template <
class Container,
class Flu
idState>
130 OPM_THROW(NotImplemented,
131 "The pcnw(fs) method is not yet implemented");
145 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
146 static Evaluation
pcnw(
const Params& ,
const FluidState& )
148 OPM_THROW(NotImplemented,
149 "The pcnw(fs) method is not yet implemented");
152 template <
class Evaluation>
153 static Evaluation twoPhaseSatPcnw(
const Params& params,
const Evaluation& SwScaled)
156 const Evaluation& pcUnscaled = EffLaw::twoPhaseSatPcnw(params.effectiveLawParams(), SwUnscaled);
157 return unscaledToScaledPcnw_(params, pcUnscaled);
160 template <
class Evaluation>
161 static Evaluation twoPhaseSatPcnwInv(
const Params& params,
const Evaluation& pcnwScaled)
163 Evaluation pcnwUnscaled = scaledToUnscaledPcnw_(params, pcnwScaled);
164 Evaluation SwUnscaled = EffLaw::twoPhaseSatPcnwInv(params.effectiveLawParams(), pcnwUnscaled);
165 return unscaledToScaledSatPc(params, SwUnscaled);
171 template <
class Container,
class Flu
idState>
172 static void saturations(Container& ,
const Params& ,
const FluidState& )
174 OPM_THROW(NotImplemented,
175 "The saturations(fs) method is not yet implemented");
182 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
183 static Evaluation
Sw(
const Params& ,
const FluidState& )
185 OPM_THROW(NotImplemented,
186 "The Sw(fs) method is not yet implemented");
189 template <
class Evaluation>
190 static Evaluation twoPhaseSatSw(
const Params& ,
const Evaluation& )
192 OPM_THROW(NotImplemented,
193 "The twoPhaseSatSw(pc) method is not yet implemented");
200 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
201 static Evaluation
Sn(
const Params& ,
const FluidState& )
203 OPM_THROW(NotImplemented,
204 "The Sn(pc) method is not yet implemented");
207 template <
class Evaluation>
208 static Evaluation twoPhaseSatSn(
const Params& ,
const Evaluation& )
210 OPM_THROW(NotImplemented,
211 "The twoPhaseSatSn(pc) method is not yet implemented");
223 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
224 static Evaluation
krw(
const Params& ,
const FluidState& )
226 OPM_THROW(NotImplemented,
227 "The krw(fs) method is not yet implemented");
230 template <
class Evaluation>
231 static Evaluation twoPhaseSatKrw(
const Params& params,
const Evaluation& SwScaled)
234 const Evaluation& krwUnscaled = EffLaw::twoPhaseSatKrw(params.effectiveLawParams(), SwUnscaled);
235 return unscaledToScaledKrw_(params, krwUnscaled);
238 template <
class Evaluation>
239 static Evaluation twoPhaseSatKrwInv(
const Params& params,
const Evaluation& krwScaled)
241 Evaluation krwUnscaled = scaledToUnscaledKrw_(params, krwScaled);
242 Evaluation SwUnscaled = EffLaw::twoPhaseSatKrwInv(params.effectiveLawParams(), krwUnscaled);
243 return unscaledToScaledSatKrw(params, SwUnscaled);
249 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
250 static Evaluation
krn(
const Params& ,
const FluidState& )
252 OPM_THROW(NotImplemented,
253 "The krn(fs) method is not yet implemented");
256 template <
class Evaluation>
257 static Evaluation twoPhaseSatKrn(
const Params& params,
const Evaluation& SwScaled)
260 const Evaluation& krnUnscaled = EffLaw::twoPhaseSatKrn(params.effectiveLawParams(), SwUnscaled);
261 return unscaledToScaledKrn_(params, krnUnscaled);
264 template <
class Evaluation>
265 static Evaluation twoPhaseSatKrnInv(
const Params& params,
const Evaluation& krnScaled)
267 Evaluation krnUnscaled = scaledToUnscaledKrn_(params, krnScaled);
268 Evaluation SwUnscaled = EffLaw::twoPhaseSatKrnInv(params.effectiveLawParams(), krnUnscaled);
269 return unscaledToScaledSatKrn(params, SwUnscaled);
277 template <
class Evaluation>
280 if (!params.config().enableSatScaling())
285 return scaledToUnscaledSatTwoPoint_(SwScaled,
286 params.unscaledPoints().saturationPcPoints(),
287 params.scaledPoints().saturationPcPoints());
290 template <
class Evaluation>
291 static Evaluation unscaledToScaledSatPc(
const Params& params,
const Evaluation& SwUnscaled)
293 if (!params.config().enableSatScaling())
298 return unscaledToScaledSatTwoPoint_(SwUnscaled,
299 params.unscaledPoints().saturationPcPoints(),
300 params.scaledPoints().saturationPcPoints());
307 template <
class Evaluation>
310 if (!params.config().enableSatScaling())
313 if (params.config().enableThreePointKrSatScaling()) {
314 return scaledToUnscaledSatThreePoint_(SwScaled,
315 params.unscaledPoints().saturationKrwPoints(),
316 params.scaledPoints().saturationKrwPoints());
319 return scaledToUnscaledSatTwoPoint_(SwScaled,
320 params.unscaledPoints().saturationKrwPoints(),
321 params.scaledPoints().saturationKrwPoints());
325 template <
class Evaluation>
326 static Evaluation unscaledToScaledSatKrw(
const Params& params,
const Evaluation& SwUnscaled)
328 if (!params.config().enableSatScaling())
331 if (params.config().enableThreePointKrSatScaling()) {
332 return unscaledToScaledSatThreePoint_(SwUnscaled,
333 params.unscaledPoints().saturationKrwPoints(),
334 params.scaledPoints().saturationKrwPoints());
337 return unscaledToScaledSatTwoPoint_(SwUnscaled,
338 params.unscaledPoints().saturationKrwPoints(),
339 params.scaledPoints().saturationKrwPoints());
347 template <
class Evaluation>
350 if (!params.config().enableSatScaling())
353 if (params.config().enableThreePointKrSatScaling())
354 return scaledToUnscaledSatThreePoint_(SwScaled,
355 params.unscaledPoints().saturationKrnPoints(),
356 params.scaledPoints().saturationKrnPoints());
358 return scaledToUnscaledSatTwoPoint_(SwScaled,
359 params.unscaledPoints().saturationKrnPoints(),
360 params.scaledPoints().saturationKrnPoints());
364 template <
class Evaluation>
365 static Evaluation unscaledToScaledSatKrn(
const Params& params,
const Evaluation& SwUnscaled)
367 if (!params.config().enableSatScaling())
370 if (params.config().enableThreePointKrSatScaling()) {
371 return unscaledToScaledSatThreePoint_(SwUnscaled,
372 params.unscaledPoints().saturationKrnPoints(),
373 params.scaledPoints().saturationKrnPoints());
376 return unscaledToScaledSatTwoPoint_(SwUnscaled,
377 params.unscaledPoints().saturationKrnPoints(),
378 params.scaledPoints().saturationKrnPoints());
383 template <
class Evaluation,
class Po
intsContainer>
384 static Evaluation scaledToUnscaledSatTwoPoint_(
const Evaluation& scaledSat,
385 const PointsContainer& unscaledSats,
386 const PointsContainer& scaledSats)
391 (scaledSat - scaledSats[0])*((unscaledSats[1] - unscaledSats[0])/(scaledSats[1] - scaledSats[0]));
394 template <
class Evaluation,
class Po
intsContainer>
395 static Evaluation unscaledToScaledSatTwoPoint_(
const Evaluation& unscaledSat,
396 const PointsContainer& unscaledSats,
397 const PointsContainer& scaledSats)
402 (unscaledSat - unscaledSats[0])*((scaledSats[1] - scaledSats[0])/(unscaledSats[1] - unscaledSats[0]));
405 template <
class Evaluation,
class Po
intsContainer>
406 static Evaluation scaledToUnscaledSatThreePoint_(
const Evaluation& scaledSat,
407 const PointsContainer& unscaledSats,
408 const PointsContainer& scaledSats)
410 if (unscaledSats[1] >= unscaledSats[2])
411 return scaledToUnscaledSatTwoPoint_(scaledSat, unscaledSats, scaledSats);
413 if (scaledSat < scaledSats[1]) {
414 Scalar delta = scaledSats[1] - scaledSats[0];
421 (scaledSat - scaledSats[0])*((unscaledSats[1] - unscaledSats[0])/delta);
424 Scalar delta = scaledSats[2] - scaledSats[1];
431 (scaledSat - scaledSats[1])*((unscaledSats[2] - unscaledSats[1])/delta);
435 template <
class Evaluation,
class Po
intsContainer>
436 static Evaluation unscaledToScaledSatThreePoint_(
const Evaluation& unscaledSat,
437 const PointsContainer& unscaledSats,
438 const PointsContainer& scaledSats)
440 if (unscaledSats[1] >= unscaledSats[2])
441 return unscaledToScaledSatTwoPoint_(unscaledSat, unscaledSats, scaledSats);
443 if (unscaledSat < unscaledSats[1]) {
444 Scalar delta = unscaledSats[1] - unscaledSats[0];
451 (unscaledSat - unscaledSats[0])*((scaledSats[1] - scaledSats[0])/delta);
454 Scalar delta = unscaledSats[2] - unscaledSats[1];
461 (unscaledSat - unscaledSats[1])*((scaledSats[2] - scaledSats[1])/delta);
468 template <
class Evaluation>
469 static Evaluation unscaledToScaledPcnw_(
const Params& params,
const Evaluation& unscaledPcnw)
471 if (params.config().enableLeverettScaling()) {
472 Scalar alpha = params.scaledPoints().leverettFactor();
473 return unscaledPcnw*alpha;
475 else if (params.config().enablePcScaling()) {
476 Scalar alpha = params.scaledPoints().maxPcnw()/params.unscaledPoints().maxPcnw();
477 return unscaledPcnw*alpha;
483 template <
class Evaluation>
484 static Evaluation scaledToUnscaledPcnw_(
const Params& params,
const Evaluation& scaledPcnw)
486 if (params.config().enableLeverettScaling()) {
487 Scalar alpha = params.scaledPoints().leverettFactor();
488 return scaledPcnw/alpha;
490 else if (params.config().enablePcScaling()) {
491 Scalar alpha = params.unscaledPoints().maxPcnw()/params.scaledPoints().maxPcnw();
492 return scaledPcnw/alpha;
501 template <
class Evaluation>
502 static Evaluation unscaledToScaledKrw_(
const Params& params,
const Evaluation& unscaledKrw)
504 if (!params.config().enableKrwScaling())
508 Scalar alpha = params.scaledPoints().maxKrw()/params.unscaledPoints().maxKrw();
509 return unscaledKrw*alpha;
512 template <
class Evaluation>
513 static Evaluation scaledToUnscaledKrw_(
const Params& params,
const Evaluation& scaledKrw)
515 if (!params.config().enableKrwScaling())
518 Scalar alpha = params.unscaledPoints().maxKrw()/params.scaledPoints().maxKrw();
519 return scaledKrw*alpha;
525 template <
class Evaluation>
526 static Evaluation unscaledToScaledKrn_(
const Params& params,
const Evaluation& unscaledKrn)
528 if (!params.config().enableKrnScaling())
532 Scalar alpha = params.scaledPoints().maxKrn()/params.unscaledPoints().maxKrn();
533 return unscaledKrn*alpha;
536 template <
class Evaluation>
537 static Evaluation scaledToUnscaledKrn_(
const Params& params,
const Evaluation& scaledKrn)
539 if (!params.config().enableKrnScaling())
542 Scalar alpha = params.unscaledPoints().maxKrn()/params.scaledPoints().maxKrn();
543 return scaledKrn*alpha;
static void relativePermeabilities(Container &, const Params &, const FluidState &)
The relative permeability-saturation curves depending on absolute saturations.
Definition: EclEpsTwoPhaseLaw.hpp:128
static Evaluation pcnw(const Params &, const FluidState &)
The capillary pressure-saturation curve.
Definition: EclEpsTwoPhaseLaw.hpp:146
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure...
Definition: EclEpsTwoPhaseLaw.hpp:90
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclEpsTwoPhaseLaw.hpp:70
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition...
Definition: EclEpsTwoPhaseLaw.hpp:98
static Evaluation Sw(const Params &, const FluidState &)
Calculate wetting liquid phase saturation given that the rest of the fluid state has been initialized...
Definition: EclEpsTwoPhaseLaw.hpp:183
static void capillaryPressures(Container &, const Params &, const FluidState &)
The capillary pressure-saturation curves depending on absolute saturations.
Definition: EclEpsTwoPhaseLaw.hpp:111
static Evaluation scaledToUnscaledSatKrn(const Params ¶ms, const Evaluation &SwScaled)
Convert an absolute saturation to an effective one for the scaling of the relperm of the non-wetting ...
Definition: EclEpsTwoPhaseLaw.hpp:348
A default implementation of the parameters for the material law adapter class which implements ECL en...
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclEpsTwoPhaseLaw.hpp:86
static Evaluation scaledToUnscaledSatKrw(const Params ¶ms, const Evaluation &SwScaled)
Convert an absolute saturation to an effective one for the scaling of the relperm of the wetting phas...
Definition: EclEpsTwoPhaseLaw.hpp:308
static Evaluation krw(const Params &, const FluidState &)
The relative permeability for the wetting phase.
Definition: EclEpsTwoPhaseLaw.hpp:224
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclEpsTwoPhaseLaw.hpp:94
This is a fluid state which allows to set the fluid saturations and takes all other quantities from a...
static void saturations(Container &, const Params &, const FluidState &)
The saturation-capillary pressure curves.
Definition: EclEpsTwoPhaseLaw.hpp:172
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclEpsTwoPhaseLaw.hpp:78
static const int numPhases
The number of fluid phases.
Definition: EclEpsTwoPhaseLaw.hpp:63
This material law takes a material law defined for unscaled saturation and converts it to a material ...
Definition: EclEpsTwoPhaseLaw.hpp:50
static Evaluation Sn(const Params &, const FluidState &)
Calculate non-wetting liquid phase saturation given that the rest of the fluid state has been initial...
Definition: EclEpsTwoPhaseLaw.hpp:201
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting phase.
Definition: EclEpsTwoPhaseLaw.hpp:250
static Evaluation scaledToUnscaledSatPc(const Params ¶ms, const Evaluation &SwScaled)
Convert an absolute saturation to an effective one for capillary pressure.
Definition: EclEpsTwoPhaseLaw.hpp:278