My Project
EclEpsTwoPhaseLaw.hpp
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
27 #ifndef OPM_ECL_EPS_TWO_PHASE_LAW_HPP
28 #define OPM_ECL_EPS_TWO_PHASE_LAW_HPP
29 
31 
33 
34 #include <algorithm>
35 #include <cstddef>
36 #include <type_traits>
37 
38 namespace Opm {
50 template <class EffLawT,
51  class ParamsT = EclEpsTwoPhaseLawParams<EffLawT> >
52 class EclEpsTwoPhaseLaw : public EffLawT::Traits
53 {
54  typedef EffLawT EffLaw;
55 
56 public:
57  typedef typename EffLaw::Traits Traits;
58  typedef ParamsT Params;
59  typedef typename EffLaw::Scalar Scalar;
60 
61  enum { wettingPhaseIdx = Traits::wettingPhaseIdx };
62  enum { nonWettingPhaseIdx = Traits::nonWettingPhaseIdx };
63 
65  static const int numPhases = EffLaw::numPhases;
66  static_assert(numPhases == 2,
67  "The endpoint scaling applies to the nested twophase laws, not to "
68  "the threephase one!");
69 
72  static const bool implementsTwoPhaseApi = true;
73 
74  static_assert(EffLaw::implementsTwoPhaseApi,
75  "The material laws put into EclEpsTwoPhaseLaw must implement the "
76  "two-phase material law API!");
77 
80  static const bool implementsTwoPhaseSatApi = true;
81 
82  static_assert(EffLaw::implementsTwoPhaseSatApi,
83  "The material laws put into EclEpsTwoPhaseLaw must implement the "
84  "two-phase material law saturation API!");
85 
88  static const bool isSaturationDependent = true;
89 
92  static const bool isPressureDependent = false;
93 
96  static const bool isTemperatureDependent = false;
97 
100  static const bool isCompositionDependent = false;
101 
112  template <class Container, class FluidState>
113  static void capillaryPressures(Container& /*values*/, const Params& /*params*/, const FluidState& /*fluidState*/)
114  {
115  throw std::invalid_argument("The capillaryPressures(fs) method is not yet implemented");
116  }
117 
128  template <class Container, class FluidState>
129  static void relativePermeabilities(Container& /*values*/, const Params& /*params*/, const FluidState& /*fluidState*/)
130  {
131  throw std::invalid_argument("The pcnw(fs) method is not yet implemented");
132  }
133 
145  template <class FluidState, class Evaluation = typename FluidState::Scalar>
146  static Evaluation pcnw(const Params& /*params*/, const FluidState& /*fluidState*/)
147  {
148  throw std::invalid_argument("The pcnw(fs) method is not yet implemented");
149  }
150 
151  template <class Evaluation>
152  static Evaluation twoPhaseSatPcnw(const Params& params, const Evaluation& SwScaled)
153  {
154  const Evaluation SwUnscaled = scaledToUnscaledSatPc(params, SwScaled);
155  const Evaluation pcUnscaled = EffLaw::twoPhaseSatPcnw(params.effectiveLawParams(), SwUnscaled);
156  return unscaledToScaledPcnw_(params, pcUnscaled);
157  }
158 
159  template <class Evaluation>
160  static Evaluation twoPhaseSatPcnwInv(const Params& params, const Evaluation& pcnwScaled)
161  {
162  const Evaluation pcnwUnscaled = scaledToUnscaledPcnw_(params, pcnwScaled);
163  const Evaluation SwUnscaled = EffLaw::twoPhaseSatPcnwInv(params.effectiveLawParams(), pcnwUnscaled);
164  return unscaledToScaledSatPc(params, SwUnscaled);
165  }
166 
170  template <class Container, class FluidState>
171  static void saturations(Container& /*values*/, const Params& /*params*/, const FluidState& /*fluidState*/)
172  {
173  throw std::invalid_argument("The saturations(fs) method is not yet implemented");
174  }
175 
180  template <class FluidState, class Evaluation = typename FluidState::Scalar>
181  static Evaluation Sw(const Params& /*params*/, const FluidState& /*fluidState*/)
182  {
183  throw std::invalid_argument("The Sw(fs) method is not yet implemented");
184  }
185 
186  template <class Evaluation>
187  static Evaluation twoPhaseSatSw(const Params& /*params*/, const Evaluation& /*pc*/)
188  {
189  throw std::invalid_argument("The twoPhaseSatSw(pc) method is not yet implemented");
190  }
191 
196  template <class FluidState, class Evaluation = typename FluidState::Scalar>
197  static Evaluation Sn(const Params& /*params*/, const FluidState& /*fluidState*/)
198  {
199  throw std::invalid_argument("The Sn(pc) method is not yet implemented");
200  }
201 
202  template <class Evaluation>
203  static Evaluation twoPhaseSatSn(const Params& /*params*/, const Evaluation& /*pc*/)
204  {
205  throw std::invalid_argument("The twoPhaseSatSn(pc) method is not yet implemented");
206  }
207 
217  template <class FluidState, class Evaluation = typename FluidState::Scalar>
218  static Evaluation krw(const Params& /*params*/, const FluidState& /*fluidState*/)
219  {
220  throw std::invalid_argument("The krw(fs) method is not yet implemented");
221  }
222 
223  template <class Evaluation>
224  static Evaluation twoPhaseSatKrw(const Params& params, const Evaluation& SwScaled)
225  {
226  const Evaluation SwUnscaled = scaledToUnscaledSatKrw(params, SwScaled);
227  const Evaluation krwUnscaled = EffLaw::twoPhaseSatKrw(params.effectiveLawParams(), SwUnscaled);
228  return unscaledToScaledKrw_(SwScaled, params, krwUnscaled);
229  }
230 
231  template <class Evaluation>
232  static Evaluation twoPhaseSatKrwInv(const Params& params, const Evaluation& krwScaled)
233  {
234  const Evaluation krwUnscaled = scaledToUnscaledKrw_(params, krwScaled);
235  const Evaluation SwUnscaled = EffLaw::twoPhaseSatKrwInv(params.effectiveLawParams(), krwUnscaled);
236  return unscaledToScaledSatKrw(params, SwUnscaled);
237  }
238 
242  template <class FluidState, class Evaluation = typename FluidState::Scalar>
243  static Evaluation krn(const Params& /*params*/, const FluidState& /*fluidState*/)
244  {
245  throw std::invalid_argument("The krn(fs) method is not yet implemented");
246  }
247 
248  template <class Evaluation>
249  static Evaluation twoPhaseSatKrn(const Params& params, const Evaluation& SwScaled)
250  {
251  const Evaluation SwUnscaled = scaledToUnscaledSatKrn(params, SwScaled);
252  const Evaluation krnUnscaled = EffLaw::twoPhaseSatKrn(params.effectiveLawParams(), SwUnscaled);
253  return unscaledToScaledKrn_(SwScaled, params, krnUnscaled);
254  }
255 
256  template <class Evaluation>
257  static Evaluation twoPhaseSatKrnInv(const Params& params, const Evaluation& krnScaled)
258  {
259  const Evaluation krnUnscaled = scaledToUnscaledKrn_(params, krnScaled);
260  const Evaluation SwUnscaled = EffLaw::twoPhaseSatKrnInv(params.effectiveLawParams(), krnUnscaled);
261  return unscaledToScaledSatKrn(params, SwUnscaled);
262  }
263 
269  template <class Evaluation>
270  static Evaluation scaledToUnscaledSatPc(const Params& params, const Evaluation& SwScaled)
271  {
272  if (!params.config().enableSatScaling())
273  return SwScaled;
274 
275  // the saturations of capillary pressure are always scaled using two-point
276  // scaling
277  return scaledToUnscaledSatTwoPoint_(SwScaled,
278  params.unscaledPoints().saturationPcPoints(),
279  params.scaledPoints().saturationPcPoints());
280  }
281 
282  template <class Evaluation>
283  static Evaluation unscaledToScaledSatPc(const Params& params, const Evaluation& SwUnscaled)
284  {
285  if (!params.config().enableSatScaling())
286  return SwUnscaled;
287 
288  // the saturations of capillary pressure are always scaled using two-point
289  // scaling
290  return unscaledToScaledSatTwoPoint_(SwUnscaled,
291  params.unscaledPoints().saturationPcPoints(),
292  params.scaledPoints().saturationPcPoints());
293  }
294 
299  template <class Evaluation>
300  static Evaluation scaledToUnscaledSatKrw(const Params& params, const Evaluation& SwScaled)
301  {
302  if (!params.config().enableSatScaling())
303  return SwScaled;
304 
305  if (params.config().enableThreePointKrSatScaling()) {
306  return scaledToUnscaledSatThreePoint_(SwScaled,
307  params.unscaledPoints().saturationKrwPoints(),
308  params.scaledPoints().saturationKrwPoints());
309  }
310  else { // two-point relperm saturation scaling
311  return scaledToUnscaledSatTwoPoint_(SwScaled,
312  params.unscaledPoints().saturationKrwPoints(),
313  params.scaledPoints().saturationKrwPoints());
314  }
315  }
316 
317  template <class Evaluation>
318  static Evaluation unscaledToScaledSatKrw(const Params& params, const Evaluation& SwUnscaled)
319  {
320  if (!params.config().enableSatScaling())
321  return SwUnscaled;
322 
323  if (params.config().enableThreePointKrSatScaling()) {
324  return unscaledToScaledSatThreePoint_(SwUnscaled,
325  params.unscaledPoints().saturationKrwPoints(),
326  params.scaledPoints().saturationKrwPoints());
327  }
328  else { // two-point relperm saturation scaling
329  return unscaledToScaledSatTwoPoint_(SwUnscaled,
330  params.unscaledPoints().saturationKrwPoints(),
331  params.scaledPoints().saturationKrwPoints());
332  }
333  }
334 
339  template <class Evaluation>
340  static Evaluation scaledToUnscaledSatKrn(const Params& params, const Evaluation& SwScaled)
341  {
342  if (!params.config().enableSatScaling())
343  return SwScaled;
344 
345  if (params.config().enableThreePointKrSatScaling())
346  return scaledToUnscaledSatThreePoint_(SwScaled,
347  params.unscaledPoints().saturationKrnPoints(),
348  params.scaledPoints().saturationKrnPoints());
349  else // two-point relperm saturation scaling
350  return scaledToUnscaledSatTwoPoint_(SwScaled,
351  params.unscaledPoints().saturationKrnPoints(),
352  params.scaledPoints().saturationKrnPoints());
353  }
354 
355 
356  template <class Evaluation>
357  static Evaluation unscaledToScaledSatKrn(const Params& params, const Evaluation& SwUnscaled)
358  {
359  if (!params.config().enableSatScaling())
360  return SwUnscaled;
361 
362  if (params.config().enableThreePointKrSatScaling()) {
363  return unscaledToScaledSatThreePoint_(SwUnscaled,
364  params.unscaledPoints().saturationKrnPoints(),
365  params.scaledPoints().saturationKrnPoints());
366  }
367  else { // two-point relperm saturation scaling
368  return unscaledToScaledSatTwoPoint_(SwUnscaled,
369  params.unscaledPoints().saturationKrnPoints(),
370  params.scaledPoints().saturationKrnPoints());
371  }
372  }
373 
374 private:
375  template <class Evaluation, class PointsContainer>
376  static Evaluation scaledToUnscaledSatTwoPoint_(const Evaluation& scaledSat,
377  const PointsContainer& unscaledSats,
378  const PointsContainer& scaledSats)
379  {
380  return
381  unscaledSats[0]
382  +
383  (scaledSat - scaledSats[0])*((unscaledSats[2] - unscaledSats[0])/(scaledSats[2] - scaledSats[0]));
384  }
385 
386  template <class Evaluation, class PointsContainer>
387  static Evaluation unscaledToScaledSatTwoPoint_(const Evaluation& unscaledSat,
388  const PointsContainer& unscaledSats,
389  const PointsContainer& scaledSats)
390  {
391  return
392  scaledSats[0]
393  +
394  (unscaledSat - unscaledSats[0])*((scaledSats[2] - scaledSats[0])/(unscaledSats[2] - unscaledSats[0]));
395  }
396 
397  template <class Evaluation, class PointsContainer>
398  static Evaluation scaledToUnscaledSatThreePoint_(const Evaluation& scaledSat,
399  const PointsContainer& unscaledSats,
400  const PointsContainer& scaledSats)
401  {
402  using UnscaledSat = std::remove_cv_t<std::remove_reference_t<decltype(unscaledSats[0])>>;
403 
404  auto map = [&scaledSat, &unscaledSats, &scaledSats](const std::size_t i)
405  {
406  const auto distance = (scaledSat - scaledSats[i])
407  / (scaledSats[i + 1] - scaledSats[i]);
408 
409  const auto displacement =
410  std::max(unscaledSats[i + 1] - unscaledSats[i], UnscaledSat{ 0 });
411 
412  return std::min(unscaledSats[i] + distance*displacement,
413  Evaluation { unscaledSats[i + 1] });
414  };
415 
416  if (! (scaledSat > scaledSats[0])) {
417  // s <= sL
418  return unscaledSats[0];
419  }
420  else if (scaledSat < std::min(scaledSats[1], scaledSats[2])) {
421  // Scaled saturation in interval [sL, sR).
422  // Map to tabulated saturation in [unscaledSats[0], unscaledSats[1]).
423  return map(0);
424  }
425  else if (scaledSat < scaledSats[2]) {
426  // Scaled saturation in interval [sR, sU); sR guaranteed to be less
427  // than sU from previous condition. Map to tabulated saturation in
428  // [unscaledSats[1], unscaledSats[2]).
429  return map(1);
430  }
431  else {
432  // s >= sU
433  return unscaledSats[2];
434  }
435  }
436 
437  template <class Evaluation, class PointsContainer>
438  static Evaluation unscaledToScaledSatThreePoint_(const Evaluation& unscaledSat,
439  const PointsContainer& unscaledSats,
440  const PointsContainer& scaledSats)
441  {
442  using ScaledSat = std::remove_cv_t<std::remove_reference_t<decltype(scaledSats[0])>>;
443 
444  auto map = [&unscaledSat, &unscaledSats, &scaledSats](const std::size_t i)
445  {
446  const auto distance = (unscaledSat - unscaledSats[i])
447  / (unscaledSats[i + 1] - unscaledSats[i]);
448 
449  const auto displacement =
450  std::max(scaledSats[i + 1] - scaledSats[i], ScaledSat{ 0 });
451 
452  return std::min(scaledSats[i] + distance*displacement,
453  Evaluation { scaledSats[i + 1] });
454  };
455 
456  if (! (unscaledSat > unscaledSats[0])) {
457  return scaledSats[0];
458  }
459  else if (unscaledSat < unscaledSats[1]) {
460  // Tabulated saturation in interval [unscaledSats[0], unscaledSats[1]).
461  // Map to scaled saturation in [sL, sR).
462  return map(0);
463  }
464  else if (unscaledSat < unscaledSats[2]) {
465  // Tabulated saturation in interval [unscaledSats[1], unscaledSats[2]).
466  // Map to scaled saturation in [sR, sU).
467  return map(1);
468  }
469  else {
470  return scaledSats[2];
471  }
472  }
473 
477  template <class Evaluation>
478  static Evaluation unscaledToScaledPcnw_(const Params& params, const Evaluation& unscaledPcnw)
479  {
480  if (params.config().enableLeverettScaling()) {
481  Scalar alpha = params.scaledPoints().leverettFactor();
482  return unscaledPcnw*alpha;
483  }
484  else if (params.config().enablePcScaling()) {
485  Scalar alpha = params.scaledPoints().maxPcnw()/params.unscaledPoints().maxPcnw();
486  return unscaledPcnw*alpha;
487  }
488 
489  return unscaledPcnw;
490  }
491 
492  template <class Evaluation>
493  static Evaluation scaledToUnscaledPcnw_(const Params& params, const Evaluation& scaledPcnw)
494  {
495  if (params.config().enableLeverettScaling()) {
496  Scalar alpha = params.scaledPoints().leverettFactor();
497  return scaledPcnw/alpha;
498  }
499  else if (params.config().enablePcScaling()) {
500  Scalar alpha = params.unscaledPoints().maxPcnw()/params.scaledPoints().maxPcnw();
501  return scaledPcnw/alpha;
502  }
503 
504  return scaledPcnw;
505  }
506 
510  template <class Evaluation>
511  static Evaluation unscaledToScaledKrw_(const Evaluation& SwScaled,
512  const Params& params,
513  const Evaluation& unscaledKrw)
514  {
515  const auto& cfg = params.config();
516 
517  if (! cfg.enableKrwScaling())
518  return unscaledKrw;
519 
520  const auto& scaled = params.scaledPoints();
521  const auto& unscaled = params.unscaledPoints();
522 
523  if (! cfg.enableThreePointKrwScaling()) {
524  // Simple case: Run uses pure vertical scaling of water relperm (keyword KRW)
525  const Scalar alpha = scaled.maxKrw() / unscaled.maxKrw();
526  return unscaledKrw * alpha;
527  }
528 
529  // Otherwise, run uses three-point vertical scaling (keywords KRWR and KRW)
530  const auto fdisp = unscaled.krwr();
531  const auto fmax = unscaled.maxKrw();
532 
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();
537 
538  if (! (SwScaled > sr)) {
539  // Pure vertical scaling in left interval ([SWL, SR])
540  return unscaledKrw * (fr / fdisp);
541  }
542  else if (fmax > fdisp) {
543  // s \in [sr, sm), sm > sr; normal case: Kr(Smax) > Kr(Sr).
544  //
545  // Linear function between (sr,fr) and (sm,fm) in terms of
546  // function value 'unscaledKrw'. This usually alters the shape
547  // of the relative permeability function in this interval (e.g.,
548  // roughly quadratic to linear).
549  const auto t = (unscaledKrw - fdisp) / (fmax - fdisp);
550 
551  return fr + t*(fm - fr);
552  }
553  else if (sr < sm) {
554  // s \in [sr, sm), sm > sr; special case: Kr(Smax) == Kr(Sr).
555  //
556  // Linear function between (sr,fr) and (sm,fm) in terms of
557  // saturation value 'SwScaled'. This usually alters the shape
558  // of the relative permeability function in this interval (e.g.,
559  // roughly quadratic to linear).
560  const auto t = (SwScaled - sr) / (sm - sr);
561 
562  return fr + t*(fm - fr);
563  }
564  else {
565  // sm == sr (pure scaling). Almost arbitrarily pick 'fm'.
566  return fm;
567  }
568  }
569 
570  template <class Evaluation>
571  static Evaluation scaledToUnscaledKrw_(const Params& params, const Evaluation& scaledKrw)
572  {
573  if (!params.config().enableKrwScaling())
574  return scaledKrw;
575 
576  Scalar alpha = params.unscaledPoints().maxKrw()/params.scaledPoints().maxKrw();
577  return scaledKrw*alpha;
578  }
579 
583  template <class Evaluation>
584  static Evaluation unscaledToScaledKrn_(const Evaluation& SwScaled,
585  const Params& params,
586  const Evaluation& unscaledKrn)
587  {
588  const auto& cfg = params.config();
589 
590  if (! cfg.enableKrnScaling())
591  return unscaledKrn;
592 
593  const auto& scaled = params.scaledPoints();
594  const auto& unscaled = params.unscaledPoints();
595 
596  if (! cfg.enableThreePointKrnScaling()) {
597  // Simple case: Run uses pure vertical scaling of non-wetting
598  // phase's relative permeability (e.g., KRG)
599  const Scalar alpha = scaled.maxKrn() / unscaled.maxKrn();
600  return unscaledKrn * alpha;
601  }
602 
603  // Otherwise, run uses three-point vertical scaling (e.g., keywords KRGR and KRG)
604  const auto fdisp = unscaled.krnr();
605  const auto fmax = unscaled.maxKrn();
606 
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();
611 
612  // Note logic here. Krn is a decreasing function of Sw (dKrn/dSw <=
613  // 0) so the roles of left and right intervals are reversed viz
614  // unscaledToScaledKrw_().
615 
616  if (! (SwScaled < sr)) {
617  // Pure vertical scaling in right-hand interval ([SR, SWU])
618  return unscaledKrn * (fr / fdisp);
619  }
620  else if (fmax > fdisp) {
621  // s \in [SWL, SR), SWL < SR; normal case: Kr(Swl) > Kr(Sr).
622  //
623  // Linear function between (sr,fr) and (sl,fm) in terms of
624  // function value 'unscaledKrn'. This usually alters the shape
625  // of the relative permeability function in this interval (e.g.,
626  // roughly quadratic to linear).
627  const auto t = (unscaledKrn - fdisp) / (fmax - fdisp);
628 
629  return fr + t*(fm - fr);
630  }
631  else if (sr > sl) {
632  // s \in [SWL, SR), SWL < SR; special case: Kr(Swl) == Kr(Sr).
633  //
634  // Linear function between (sr,fr) and (sl,fm) in terms of
635  // saturation value 'SwScaled'. This usually alters the shape
636  // of the relative permeability function in this interval (e.g.,
637  // roughly quadratic to linear).
638  const auto t = (sr - SwScaled) / (sr - sl);
639 
640  return fr + t*(fm - fr);
641  }
642  else {
643  // sl == sr (pure scaling). Almost arbitrarily pick 'fm'.
644  return fm;
645  }
646  }
647 
648  template <class Evaluation>
649  static Evaluation scaledToUnscaledKrn_(const Params& params, const Evaluation& scaledKrn)
650  {
651  if (!params.config().enableKrnScaling())
652  return scaledKrn;
653 
654  Scalar alpha = params.unscaledPoints().maxKrn()/params.scaledPoints().maxKrn();
655  return scaledKrn*alpha;
656  }
657 };
658 } // namespace Opm
659 
660 #endif
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 &params, 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 &params, 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 &params, 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