27 #ifndef OPM_ECL_EPS_TWO_PHASE_LAW_HPP
28 #define OPM_ECL_EPS_TWO_PHASE_LAW_HPP
36 #include <type_traits>
50 template <
class EffLawT,
51 class ParamsT = EclEpsTwoPhaseLawParams<EffLawT> >
54 typedef EffLawT EffLaw;
57 typedef typename EffLaw::Traits Traits;
58 typedef ParamsT Params;
59 typedef typename EffLaw::Scalar Scalar;
61 enum { wettingPhaseIdx = Traits::wettingPhaseIdx };
62 enum { nonWettingPhaseIdx = Traits::nonWettingPhaseIdx };
67 "The endpoint scaling applies to the nested twophase laws, not to "
68 "the threephase one!");
74 static_assert(EffLaw::implementsTwoPhaseApi,
75 "The material laws put into EclEpsTwoPhaseLaw must implement the "
76 "two-phase material law API!");
82 static_assert(EffLaw::implementsTwoPhaseSatApi,
83 "The material laws put into EclEpsTwoPhaseLaw must implement the "
84 "two-phase material law saturation API!");
112 template <
class Container,
class Flu
idState>
115 throw std::invalid_argument(
"The capillaryPressures(fs) method is not yet implemented");
128 template <
class Container,
class Flu
idState>
131 throw std::invalid_argument(
"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 throw std::invalid_argument(
"The pcnw(fs) method is not yet implemented");
151 template <
class Evaluation>
152 static Evaluation twoPhaseSatPcnw(
const Params& params,
const Evaluation& SwScaled)
155 const Evaluation pcUnscaled = EffLaw::twoPhaseSatPcnw(params.effectiveLawParams(), SwUnscaled);
156 return unscaledToScaledPcnw_(params, pcUnscaled);
159 template <
class Evaluation>
160 static Evaluation twoPhaseSatPcnwInv(
const Params& params,
const Evaluation& pcnwScaled)
162 const Evaluation pcnwUnscaled = scaledToUnscaledPcnw_(params, pcnwScaled);
163 const Evaluation SwUnscaled = EffLaw::twoPhaseSatPcnwInv(params.effectiveLawParams(), pcnwUnscaled);
164 return unscaledToScaledSatPc(params, SwUnscaled);
170 template <
class Container,
class Flu
idState>
171 static void saturations(Container& ,
const Params& ,
const FluidState& )
173 throw std::invalid_argument(
"The saturations(fs) method is not yet implemented");
180 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
181 static Evaluation
Sw(
const Params& ,
const FluidState& )
183 throw std::invalid_argument(
"The Sw(fs) method is not yet implemented");
186 template <
class Evaluation>
187 static Evaluation twoPhaseSatSw(
const Params& ,
const Evaluation& )
189 throw std::invalid_argument(
"The twoPhaseSatSw(pc) method is not yet implemented");
196 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
197 static Evaluation
Sn(
const Params& ,
const FluidState& )
199 throw std::invalid_argument(
"The Sn(pc) method is not yet implemented");
202 template <
class Evaluation>
203 static Evaluation twoPhaseSatSn(
const Params& ,
const Evaluation& )
205 throw std::invalid_argument(
"The twoPhaseSatSn(pc) method is not yet implemented");
217 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
218 static Evaluation
krw(
const Params& ,
const FluidState& )
220 throw std::invalid_argument(
"The krw(fs) method is not yet implemented");
223 template <
class Evaluation>
224 static Evaluation twoPhaseSatKrw(
const Params& params,
const Evaluation& SwScaled)
227 const Evaluation krwUnscaled = EffLaw::twoPhaseSatKrw(params.effectiveLawParams(), SwUnscaled);
228 return unscaledToScaledKrw_(SwScaled, params, krwUnscaled);
231 template <
class Evaluation>
232 static Evaluation twoPhaseSatKrwInv(
const Params& params,
const Evaluation& krwScaled)
234 const Evaluation krwUnscaled = scaledToUnscaledKrw_(params, krwScaled);
235 const Evaluation SwUnscaled = EffLaw::twoPhaseSatKrwInv(params.effectiveLawParams(), krwUnscaled);
236 return unscaledToScaledSatKrw(params, SwUnscaled);
242 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
243 static Evaluation
krn(
const Params& ,
const FluidState& )
245 throw std::invalid_argument(
"The krn(fs) method is not yet implemented");
248 template <
class Evaluation>
249 static Evaluation twoPhaseSatKrn(
const Params& params,
const Evaluation& SwScaled)
252 const Evaluation krnUnscaled = EffLaw::twoPhaseSatKrn(params.effectiveLawParams(), SwUnscaled);
253 return unscaledToScaledKrn_(SwScaled, params, krnUnscaled);
256 template <
class Evaluation>
257 static Evaluation twoPhaseSatKrnInv(
const Params& params,
const Evaluation& krnScaled)
259 const Evaluation krnUnscaled = scaledToUnscaledKrn_(params, krnScaled);
260 const Evaluation SwUnscaled = EffLaw::twoPhaseSatKrnInv(params.effectiveLawParams(), krnUnscaled);
261 return unscaledToScaledSatKrn(params, SwUnscaled);
269 template <
class Evaluation>
272 if (!params.config().enableSatScaling())
277 return scaledToUnscaledSatTwoPoint_(SwScaled,
278 params.unscaledPoints().saturationPcPoints(),
279 params.scaledPoints().saturationPcPoints());
282 template <
class Evaluation>
283 static Evaluation unscaledToScaledSatPc(
const Params& params,
const Evaluation& SwUnscaled)
285 if (!params.config().enableSatScaling())
290 return unscaledToScaledSatTwoPoint_(SwUnscaled,
291 params.unscaledPoints().saturationPcPoints(),
292 params.scaledPoints().saturationPcPoints());
299 template <
class Evaluation>
302 if (!params.config().enableSatScaling())
305 if (params.config().enableThreePointKrSatScaling()) {
306 return scaledToUnscaledSatThreePoint_(SwScaled,
307 params.unscaledPoints().saturationKrwPoints(),
308 params.scaledPoints().saturationKrwPoints());
311 return scaledToUnscaledSatTwoPoint_(SwScaled,
312 params.unscaledPoints().saturationKrwPoints(),
313 params.scaledPoints().saturationKrwPoints());
317 template <
class Evaluation>
318 static Evaluation unscaledToScaledSatKrw(
const Params& params,
const Evaluation& SwUnscaled)
320 if (!params.config().enableSatScaling())
323 if (params.config().enableThreePointKrSatScaling()) {
324 return unscaledToScaledSatThreePoint_(SwUnscaled,
325 params.unscaledPoints().saturationKrwPoints(),
326 params.scaledPoints().saturationKrwPoints());
329 return unscaledToScaledSatTwoPoint_(SwUnscaled,
330 params.unscaledPoints().saturationKrwPoints(),
331 params.scaledPoints().saturationKrwPoints());
339 template <
class Evaluation>
342 if (!params.config().enableSatScaling())
345 if (params.config().enableThreePointKrSatScaling())
346 return scaledToUnscaledSatThreePoint_(SwScaled,
347 params.unscaledPoints().saturationKrnPoints(),
348 params.scaledPoints().saturationKrnPoints());
350 return scaledToUnscaledSatTwoPoint_(SwScaled,
351 params.unscaledPoints().saturationKrnPoints(),
352 params.scaledPoints().saturationKrnPoints());
356 template <
class Evaluation>
357 static Evaluation unscaledToScaledSatKrn(
const Params& params,
const Evaluation& SwUnscaled)
359 if (!params.config().enableSatScaling())
362 if (params.config().enableThreePointKrSatScaling()) {
363 return unscaledToScaledSatThreePoint_(SwUnscaled,
364 params.unscaledPoints().saturationKrnPoints(),
365 params.scaledPoints().saturationKrnPoints());
368 return unscaledToScaledSatTwoPoint_(SwUnscaled,
369 params.unscaledPoints().saturationKrnPoints(),
370 params.scaledPoints().saturationKrnPoints());
375 template <
class Evaluation,
class Po
intsContainer>
376 static Evaluation scaledToUnscaledSatTwoPoint_(
const Evaluation& scaledSat,
377 const PointsContainer& unscaledSats,
378 const PointsContainer& scaledSats)
383 (scaledSat - scaledSats[0])*((unscaledSats[2] - unscaledSats[0])/(scaledSats[2] - scaledSats[0]));
386 template <
class Evaluation,
class Po
intsContainer>
387 static Evaluation unscaledToScaledSatTwoPoint_(
const Evaluation& unscaledSat,
388 const PointsContainer& unscaledSats,
389 const PointsContainer& scaledSats)
394 (unscaledSat - unscaledSats[0])*((scaledSats[2] - scaledSats[0])/(unscaledSats[2] - unscaledSats[0]));
397 template <
class Evaluation,
class Po
intsContainer>
398 static Evaluation scaledToUnscaledSatThreePoint_(
const Evaluation& scaledSat,
399 const PointsContainer& unscaledSats,
400 const PointsContainer& scaledSats)
402 using UnscaledSat = std::remove_cv_t<std::remove_reference_t<decltype(unscaledSats[0])>>;
404 auto map = [&scaledSat, &unscaledSats, &scaledSats](
const std::size_t i)
406 const auto distance = (scaledSat - scaledSats[i])
407 / (scaledSats[i + 1] - scaledSats[i]);
409 const auto displacement =
410 std::max(unscaledSats[i + 1] - unscaledSats[i], UnscaledSat{ 0 });
412 return std::min(unscaledSats[i] + distance*displacement,
413 Evaluation { unscaledSats[i + 1] });
416 if (! (scaledSat > scaledSats[0])) {
418 return unscaledSats[0];
420 else if (scaledSat < std::min(scaledSats[1], scaledSats[2])) {
425 else if (scaledSat < scaledSats[2]) {
433 return unscaledSats[2];
437 template <
class Evaluation,
class Po
intsContainer>
438 static Evaluation unscaledToScaledSatThreePoint_(
const Evaluation& unscaledSat,
439 const PointsContainer& unscaledSats,
440 const PointsContainer& scaledSats)
442 using ScaledSat = std::remove_cv_t<std::remove_reference_t<decltype(scaledSats[0])>>;
444 auto map = [&unscaledSat, &unscaledSats, &scaledSats](
const std::size_t i)
446 const auto distance = (unscaledSat - unscaledSats[i])
447 / (unscaledSats[i + 1] - unscaledSats[i]);
449 const auto displacement =
450 std::max(scaledSats[i + 1] - scaledSats[i], ScaledSat{ 0 });
452 return std::min(scaledSats[i] + distance*displacement,
453 Evaluation { scaledSats[i + 1] });
456 if (! (unscaledSat > unscaledSats[0])) {
457 return scaledSats[0];
459 else if (unscaledSat < unscaledSats[1]) {
464 else if (unscaledSat < unscaledSats[2]) {
470 return scaledSats[2];
477 template <
class Evaluation>
478 static Evaluation unscaledToScaledPcnw_(
const Params& params,
const Evaluation& unscaledPcnw)
480 if (params.config().enableLeverettScaling()) {
481 Scalar alpha = params.scaledPoints().leverettFactor();
482 return unscaledPcnw*alpha;
484 else if (params.config().enablePcScaling()) {
485 Scalar alpha = params.scaledPoints().maxPcnw()/params.unscaledPoints().maxPcnw();
486 return unscaledPcnw*alpha;
492 template <
class Evaluation>
493 static Evaluation scaledToUnscaledPcnw_(
const Params& params,
const Evaluation& scaledPcnw)
495 if (params.config().enableLeverettScaling()) {
496 Scalar alpha = params.scaledPoints().leverettFactor();
497 return scaledPcnw/alpha;
499 else if (params.config().enablePcScaling()) {
500 Scalar alpha = params.unscaledPoints().maxPcnw()/params.scaledPoints().maxPcnw();
501 return scaledPcnw/alpha;
510 template <
class Evaluation>
511 static Evaluation unscaledToScaledKrw_(
const Evaluation& SwScaled,
512 const Params& params,
513 const Evaluation& unscaledKrw)
515 const auto& cfg = params.config();
517 if (! cfg.enableKrwScaling())
520 const auto& scaled = params.scaledPoints();
521 const auto& unscaled = params.unscaledPoints();
523 if (! cfg.enableThreePointKrwScaling()) {
525 const Scalar alpha = scaled.maxKrw() / unscaled.maxKrw();
526 return unscaledKrw * alpha;
530 const auto fdisp = unscaled.krwr();
531 const auto fmax = unscaled.maxKrw();
533 const auto sm = scaled.saturationKrwPoints()[2];
534 const auto sr = std::min(scaled.saturationKrwPoints()[1], sm);
535 const auto fr = scaled.krwr();
536 const auto fm = scaled.maxKrw();
538 if (! (SwScaled > sr)) {
540 return unscaledKrw * (fr / fdisp);
542 else if (fmax > fdisp) {
549 const auto t = (unscaledKrw - fdisp) / (fmax - fdisp);
551 return fr + t*(fm - fr);
560 const auto t = (SwScaled - sr) / (sm - sr);
562 return fr + t*(fm - fr);
570 template <
class Evaluation>
571 static Evaluation scaledToUnscaledKrw_(
const Params& params,
const Evaluation& scaledKrw)
573 if (!params.config().enableKrwScaling())
576 Scalar alpha = params.unscaledPoints().maxKrw()/params.scaledPoints().maxKrw();
577 return scaledKrw*alpha;
583 template <
class Evaluation>
584 static Evaluation unscaledToScaledKrn_(
const Evaluation& SwScaled,
585 const Params& params,
586 const Evaluation& unscaledKrn)
588 const auto& cfg = params.config();
590 if (! cfg.enableKrnScaling())
593 const auto& scaled = params.scaledPoints();
594 const auto& unscaled = params.unscaledPoints();
596 if (! cfg.enableThreePointKrnScaling()) {
599 const Scalar alpha = scaled.maxKrn() / unscaled.maxKrn();
600 return unscaledKrn * alpha;
604 const auto fdisp = unscaled.krnr();
605 const auto fmax = unscaled.maxKrn();
607 const auto sl = scaled.saturationKrnPoints()[0];
608 const auto sr = std::max(scaled.saturationKrnPoints()[1], sl);
609 const auto fr = scaled.krnr();
610 const auto fm = scaled.maxKrn();
616 if (! (SwScaled < sr)) {
618 return unscaledKrn * (fr / fdisp);
620 else if (fmax > fdisp) {
627 const auto t = (unscaledKrn - fdisp) / (fmax - fdisp);
629 return fr + t*(fm - fr);
638 const auto t = (sr - SwScaled) / (sr - sl);
640 return fr + t*(fm - fr);
648 template <
class Evaluation>
649 static Evaluation scaledToUnscaledKrn_(
const Params& params,
const Evaluation& scaledKrn)
651 if (!params.config().enableKrnScaling())
654 Scalar alpha = params.unscaledPoints().maxKrn()/params.scaledPoints().maxKrn();
655 return scaledKrn*alpha;
A default implementation of the parameters for the material law adapter class which implements ECL en...
This is a fluid state which allows to set the fluid saturations and takes all other quantities from a...
This material law takes a material law defined for unscaled saturation and converts it to a material ...
Definition: EclEpsTwoPhaseLaw.hpp:53
static Evaluation krw(const Params &, const FluidState &)
The relative permeability for the wetting phase.
Definition: EclEpsTwoPhaseLaw.hpp:218
static const int numPhases
The number of fluid phases.
Definition: EclEpsTwoPhaseLaw.hpp:65
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:340
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:181
static void saturations(Container &, const Params &, const FluidState &)
The saturation-capillary pressure curves.
Definition: EclEpsTwoPhaseLaw.hpp:171
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting phase.
Definition: EclEpsTwoPhaseLaw.hpp:243
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclEpsTwoPhaseLaw.hpp:72
static Evaluation scaledToUnscaledSatPc(const Params ¶ms, const Evaluation &SwScaled)
Convert an absolute saturation to an effective one for capillary pressure.
Definition: EclEpsTwoPhaseLaw.hpp:270
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclEpsTwoPhaseLaw.hpp:88
static void relativePermeabilities(Container &, const Params &, const FluidState &)
The relative permeability-saturation curves depending on absolute saturations.
Definition: EclEpsTwoPhaseLaw.hpp:129
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: EclEpsTwoPhaseLaw.hpp:100
static Evaluation pcnw(const Params &, const FluidState &)
The capillary pressure-saturation curve.
Definition: EclEpsTwoPhaseLaw.hpp:146
static void capillaryPressures(Container &, const Params &, const FluidState &)
The capillary pressure-saturation curves depending on absolute saturations.
Definition: EclEpsTwoPhaseLaw.hpp:113
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:300
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:197
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclEpsTwoPhaseLaw.hpp:96
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclEpsTwoPhaseLaw.hpp:80
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: EclEpsTwoPhaseLaw.hpp:92