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 
32 #include <algorithm>
33 #include <cstddef>
34 #include <type_traits>
35 
36 namespace Opm {
48 template <class EffLawT,
49  class ParamsT = EclEpsTwoPhaseLawParams<EffLawT> >
50 class EclEpsTwoPhaseLaw : public EffLawT::Traits
51 {
52  typedef EffLawT EffLaw;
53 
54 public:
55  typedef typename EffLaw::Traits Traits;
56  typedef ParamsT Params;
57  typedef typename EffLaw::Scalar Scalar;
58 
59  enum { wettingPhaseIdx = Traits::wettingPhaseIdx };
60  enum { nonWettingPhaseIdx = Traits::nonWettingPhaseIdx };
61 
63  static const int numPhases = EffLaw::numPhases;
64  static_assert(numPhases == 2,
65  "The endpoint scaling applies to the nested twophase laws, not to "
66  "the threephase one!");
67 
70  static const bool implementsTwoPhaseApi = true;
71 
72  static_assert(EffLaw::implementsTwoPhaseApi,
73  "The material laws put into EclEpsTwoPhaseLaw must implement the "
74  "two-phase material law API!");
75 
78  static const bool implementsTwoPhaseSatApi = true;
79 
80  static_assert(EffLaw::implementsTwoPhaseSatApi,
81  "The material laws put into EclEpsTwoPhaseLaw must implement the "
82  "two-phase material law saturation API!");
83 
86  static const bool isSaturationDependent = true;
87 
90  static const bool isPressureDependent = false;
91 
94  static const bool isTemperatureDependent = false;
95 
98  static const bool isCompositionDependent = false;
99 
110  template <class Container, class FluidState>
111  static void capillaryPressures(Container& /*values*/, const Params& /*params*/, const FluidState& /*fluidState*/)
112  {
113  throw std::invalid_argument("The capillaryPressures(fs) method is not yet implemented");
114  }
115 
126  template <class Container, class FluidState>
127  static void relativePermeabilities(Container& /*values*/, const Params& /*params*/, const FluidState& /*fluidState*/)
128  {
129  throw std::invalid_argument("The pcnw(fs) method is not yet implemented");
130  }
131 
143  template <class FluidState, class Evaluation = typename FluidState::Scalar>
144  static Evaluation pcnw(const Params& /*params*/, const FluidState& /*fluidState*/)
145  {
146  throw std::invalid_argument("The pcnw(fs) method is not yet implemented");
147  }
148 
149  template <class Evaluation>
150  static Evaluation twoPhaseSatPcnw(const Params& params, const Evaluation& SwScaled)
151  {
152  const Evaluation SwUnscaled = scaledToUnscaledSatPc(params, SwScaled);
153  const Evaluation pcUnscaled = EffLaw::twoPhaseSatPcnw(params.effectiveLawParams(), SwUnscaled);
154  return unscaledToScaledPcnw_(params, pcUnscaled);
155  }
156 
157  template <class Evaluation>
158  static Evaluation twoPhaseSatPcnwInv(const Params& params, const Evaluation& pcnwScaled)
159  {
160  const Evaluation pcnwUnscaled = scaledToUnscaledPcnw_(params, pcnwScaled);
161  const Evaluation SwUnscaled = EffLaw::twoPhaseSatPcnwInv(params.effectiveLawParams(), pcnwUnscaled);
162  return unscaledToScaledSatPc(params, SwUnscaled);
163  }
164 
168  template <class Container, class FluidState>
169  static void saturations(Container& /*values*/, const Params& /*params*/, const FluidState& /*fluidState*/)
170  {
171  throw std::invalid_argument("The saturations(fs) method is not yet implemented");
172  }
173 
178  template <class FluidState, class Evaluation = typename FluidState::Scalar>
179  static Evaluation Sw(const Params& /*params*/, const FluidState& /*fluidState*/)
180  {
181  throw std::invalid_argument("The Sw(fs) method is not yet implemented");
182  }
183 
184  template <class Evaluation>
185  static Evaluation twoPhaseSatSw(const Params& /*params*/, const Evaluation& /*pc*/)
186  {
187  throw std::invalid_argument("The twoPhaseSatSw(pc) method is not yet implemented");
188  }
189 
194  template <class FluidState, class Evaluation = typename FluidState::Scalar>
195  static Evaluation Sn(const Params& /*params*/, const FluidState& /*fluidState*/)
196  {
197  throw std::invalid_argument("The Sn(pc) method is not yet implemented");
198  }
199 
200  template <class Evaluation>
201  static Evaluation twoPhaseSatSn(const Params& /*params*/, const Evaluation& /*pc*/)
202  {
203  throw std::invalid_argument("The twoPhaseSatSn(pc) method is not yet implemented");
204  }
205 
215  template <class FluidState, class Evaluation = typename FluidState::Scalar>
216  static Evaluation krw(const Params& /*params*/, const FluidState& /*fluidState*/)
217  {
218  throw std::invalid_argument("The krw(fs) method is not yet implemented");
219  }
220 
221  template <class Evaluation>
222  static Evaluation twoPhaseSatKrw(const Params& params, const Evaluation& SwScaled)
223  {
224  const Evaluation SwUnscaled = scaledToUnscaledSatKrw(params, SwScaled);
225  const Evaluation krwUnscaled = EffLaw::twoPhaseSatKrw(params.effectiveLawParams(), SwUnscaled);
226  return unscaledToScaledKrw_(SwScaled, params, krwUnscaled);
227  }
228 
229  template <class Evaluation>
230  static Evaluation twoPhaseSatKrwInv(const Params& params, const Evaluation& krwScaled)
231  {
232  const Evaluation krwUnscaled = scaledToUnscaledKrw_(params, krwScaled);
233  const Evaluation SwUnscaled = EffLaw::twoPhaseSatKrwInv(params.effectiveLawParams(), krwUnscaled);
234  return unscaledToScaledSatKrw(params, SwUnscaled);
235  }
236 
240  template <class FluidState, class Evaluation = typename FluidState::Scalar>
241  static Evaluation krn(const Params& /*params*/, const FluidState& /*fluidState*/)
242  {
243  throw std::invalid_argument("The krn(fs) method is not yet implemented");
244  }
245 
246  template <class Evaluation>
247  static Evaluation twoPhaseSatKrn(const Params& params, const Evaluation& SwScaled)
248  {
249  const Evaluation SwUnscaled = scaledToUnscaledSatKrn(params, SwScaled);
250  const Evaluation krnUnscaled = EffLaw::twoPhaseSatKrn(params.effectiveLawParams(), SwUnscaled);
251  return unscaledToScaledKrn_(SwScaled, params, krnUnscaled);
252  }
253 
254  template <class Evaluation>
255  static Evaluation twoPhaseSatKrnInv(const Params& params, const Evaluation& krnScaled)
256  {
257  const Evaluation krnUnscaled = scaledToUnscaledKrn_(params, krnScaled);
258  const Evaluation SwUnscaled = EffLaw::twoPhaseSatKrnInv(params.effectiveLawParams(), krnUnscaled);
259  return unscaledToScaledSatKrn(params, SwUnscaled);
260  }
261 
267  template <class Evaluation>
268  static Evaluation scaledToUnscaledSatPc(const Params& params, const Evaluation& SwScaled)
269  {
270  if (!params.config().enableSatScaling())
271  return SwScaled;
272 
273  // the saturations of capillary pressure are always scaled using two-point
274  // scaling
275  return scaledToUnscaledSatTwoPoint_(SwScaled,
276  params.unscaledPoints().saturationPcPoints(),
277  params.scaledPoints().saturationPcPoints());
278  }
279 
280  template <class Evaluation>
281  static Evaluation unscaledToScaledSatPc(const Params& params, const Evaluation& SwUnscaled)
282  {
283  if (!params.config().enableSatScaling())
284  return SwUnscaled;
285 
286  // the saturations of capillary pressure are always scaled using two-point
287  // scaling
288  return unscaledToScaledSatTwoPoint_(SwUnscaled,
289  params.unscaledPoints().saturationPcPoints(),
290  params.scaledPoints().saturationPcPoints());
291  }
292 
297  template <class Evaluation>
298  static Evaluation scaledToUnscaledSatKrw(const Params& params, const Evaluation& SwScaled)
299  {
300  if (!params.config().enableSatScaling())
301  return SwScaled;
302 
303  if (params.config().enableThreePointKrSatScaling()) {
304  return scaledToUnscaledSatThreePoint_(SwScaled,
305  params.unscaledPoints().saturationKrwPoints(),
306  params.scaledPoints().saturationKrwPoints());
307  }
308  else { // two-point relperm saturation scaling
309  return scaledToUnscaledSatTwoPoint_(SwScaled,
310  params.unscaledPoints().saturationKrwPoints(),
311  params.scaledPoints().saturationKrwPoints());
312  }
313  }
314 
315  template <class Evaluation>
316  static Evaluation unscaledToScaledSatKrw(const Params& params, const Evaluation& SwUnscaled)
317  {
318  if (!params.config().enableSatScaling())
319  return SwUnscaled;
320 
321  if (params.config().enableThreePointKrSatScaling()) {
322  return unscaledToScaledSatThreePoint_(SwUnscaled,
323  params.unscaledPoints().saturationKrwPoints(),
324  params.scaledPoints().saturationKrwPoints());
325  }
326  else { // two-point relperm saturation scaling
327  return unscaledToScaledSatTwoPoint_(SwUnscaled,
328  params.unscaledPoints().saturationKrwPoints(),
329  params.scaledPoints().saturationKrwPoints());
330  }
331  }
332 
337  template <class Evaluation>
338  static Evaluation scaledToUnscaledSatKrn(const Params& params, const Evaluation& SwScaled)
339  {
340  if (!params.config().enableSatScaling())
341  return SwScaled;
342 
343  if (params.config().enableThreePointKrSatScaling())
344  return scaledToUnscaledSatThreePoint_(SwScaled,
345  params.unscaledPoints().saturationKrnPoints(),
346  params.scaledPoints().saturationKrnPoints());
347  else // two-point relperm saturation scaling
348  return scaledToUnscaledSatTwoPoint_(SwScaled,
349  params.unscaledPoints().saturationKrnPoints(),
350  params.scaledPoints().saturationKrnPoints());
351  }
352 
353 
354  template <class Evaluation>
355  static Evaluation unscaledToScaledSatKrn(const Params& params, const Evaluation& SwUnscaled)
356  {
357  if (!params.config().enableSatScaling())
358  return SwUnscaled;
359 
360  if (params.config().enableThreePointKrSatScaling()) {
361  return unscaledToScaledSatThreePoint_(SwUnscaled,
362  params.unscaledPoints().saturationKrnPoints(),
363  params.scaledPoints().saturationKrnPoints());
364  }
365  else { // two-point relperm saturation scaling
366  return unscaledToScaledSatTwoPoint_(SwUnscaled,
367  params.unscaledPoints().saturationKrnPoints(),
368  params.scaledPoints().saturationKrnPoints());
369  }
370  }
371 
372 private:
373  template <class Evaluation, class PointsContainer>
374  static Evaluation scaledToUnscaledSatTwoPoint_(const Evaluation& scaledSat,
375  const PointsContainer& unscaledSats,
376  const PointsContainer& scaledSats)
377  {
378  return
379  unscaledSats[0]
380  +
381  (scaledSat - scaledSats[0])*((unscaledSats[2] - unscaledSats[0])/(scaledSats[2] - scaledSats[0]));
382  }
383 
384  template <class Evaluation, class PointsContainer>
385  static Evaluation unscaledToScaledSatTwoPoint_(const Evaluation& unscaledSat,
386  const PointsContainer& unscaledSats,
387  const PointsContainer& scaledSats)
388  {
389  return
390  scaledSats[0]
391  +
392  (unscaledSat - unscaledSats[0])*((scaledSats[2] - scaledSats[0])/(unscaledSats[2] - unscaledSats[0]));
393  }
394 
395  template <class Evaluation, class PointsContainer>
396  static Evaluation scaledToUnscaledSatThreePoint_(const Evaluation& scaledSat,
397  const PointsContainer& unscaledSats,
398  const PointsContainer& scaledSats)
399  {
400  using UnscaledSat = std::remove_cv_t<std::remove_reference_t<decltype(unscaledSats[0])>>;
401 
402  auto map = [&scaledSat, &unscaledSats, &scaledSats](const std::size_t i)
403  {
404  const auto distance = (scaledSat - scaledSats[i])
405  / (scaledSats[i + 1] - scaledSats[i]);
406 
407  const auto displacement =
408  std::max(unscaledSats[i + 1] - unscaledSats[i], UnscaledSat{ 0 });
409 
410  return std::min(unscaledSats[i] + distance*displacement,
411  Evaluation { unscaledSats[i + 1] });
412  };
413 
414  if (! (scaledSat > scaledSats[0])) {
415  // s <= sL
416  return unscaledSats[0];
417  }
418  else if (scaledSat < std::min(scaledSats[1], scaledSats[2])) {
419  // Scaled saturation in interval [sL, sR).
420  // Map to tabulated saturation in [unscaledSats[0], unscaledSats[1]).
421  return map(0);
422  }
423  else if (scaledSat < scaledSats[2]) {
424  // Scaled saturation in interval [sR, sU); sR guaranteed to be less
425  // than sU from previous condition. Map to tabulated saturation in
426  // [unscaledSats[1], unscaledSats[2]).
427  return map(1);
428  }
429  else {
430  // s >= sU
431  return unscaledSats[2];
432  }
433  }
434 
435  template <class Evaluation, class PointsContainer>
436  static Evaluation unscaledToScaledSatThreePoint_(const Evaluation& unscaledSat,
437  const PointsContainer& unscaledSats,
438  const PointsContainer& scaledSats)
439  {
440  using ScaledSat = std::remove_cv_t<std::remove_reference_t<decltype(scaledSats[0])>>;
441 
442  auto map = [&unscaledSat, &unscaledSats, &scaledSats](const std::size_t i)
443  {
444  const auto distance = (unscaledSat - unscaledSats[i])
445  / (unscaledSats[i + 1] - unscaledSats[i]);
446 
447  const auto displacement =
448  std::max(scaledSats[i + 1] - scaledSats[i], ScaledSat{ 0 });
449 
450  return std::min(scaledSats[i] + distance*displacement,
451  Evaluation { scaledSats[i + 1] });
452  };
453 
454  if (! (unscaledSat > unscaledSats[0])) {
455  return scaledSats[0];
456  }
457  else if (unscaledSat < unscaledSats[1]) {
458  // Tabulated saturation in interval [unscaledSats[0], unscaledSats[1]).
459  // Map to scaled saturation in [sL, sR).
460  return map(0);
461  }
462  else if (unscaledSat < unscaledSats[2]) {
463  // Tabulated saturation in interval [unscaledSats[1], unscaledSats[2]).
464  // Map to scaled saturation in [sR, sU).
465  return map(1);
466  }
467  else {
468  return scaledSats[2];
469  }
470  }
471 
475  template <class Evaluation>
476  static Evaluation unscaledToScaledPcnw_(const Params& params, const Evaluation& unscaledPcnw)
477  {
478  if (params.config().enableLeverettScaling()) {
479  Scalar alpha = params.scaledPoints().leverettFactor();
480  return unscaledPcnw*alpha;
481  }
482  else if (params.config().enablePcScaling()) {
483  const auto& scaled_maxPcnw = params.scaledPoints().maxPcnw();
484  const auto& unscaled_maxPcnw = params.unscaledPoints().maxPcnw();
485 
486  Scalar alpha;
487  if (scaled_maxPcnw == unscaled_maxPcnw)
488  alpha = 1.0;
489  else
490  alpha = params.scaledPoints().maxPcnw()/params.unscaledPoints().maxPcnw();
491 
492  return unscaledPcnw*alpha;
493  }
494 
495  return unscaledPcnw;
496  }
497 
498  template <class Evaluation>
499  static Evaluation scaledToUnscaledPcnw_(const Params& params, const Evaluation& scaledPcnw)
500  {
501  if (params.config().enableLeverettScaling()) {
502  Scalar alpha = params.scaledPoints().leverettFactor();
503  return scaledPcnw/alpha;
504  }
505  else if (params.config().enablePcScaling()) {
506  const auto& scaled_maxPcnw = params.scaledPoints().maxPcnw();
507  const auto& unscaled_maxPcnw = params.unscaledPoints().maxPcnw();
508 
509  Scalar alpha;
510  if (scaled_maxPcnw == unscaled_maxPcnw)
511  alpha = 1.0;
512  else
513  alpha = params.scaledPoints().maxPcnw()/params.unscaledPoints().maxPcnw();
514 
515  return scaledPcnw/alpha;
516  }
517 
518  return scaledPcnw;
519  }
520 
524  template <class Evaluation>
525  static Evaluation unscaledToScaledKrw_(const Evaluation& SwScaled,
526  const Params& params,
527  const Evaluation& unscaledKrw)
528  {
529  const auto& cfg = params.config();
530 
531  if (! cfg.enableKrwScaling())
532  return unscaledKrw;
533 
534  const auto& scaled = params.scaledPoints();
535  const auto& unscaled = params.unscaledPoints();
536 
537  if (! cfg.enableThreePointKrwScaling()) {
538  // Simple case: Run uses pure vertical scaling of water relperm (keyword KRW)
539  const Scalar alpha = scaled.maxKrw() / unscaled.maxKrw();
540  return unscaledKrw * alpha;
541  }
542 
543  // Otherwise, run uses three-point vertical scaling (keywords KRWR and KRW)
544  const auto fdisp = unscaled.krwr();
545  const auto fmax = unscaled.maxKrw();
546 
547  const auto sm = scaled.saturationKrwPoints()[2];
548  const auto sr = std::min(scaled.saturationKrwPoints()[1], sm);
549  const auto fr = scaled.krwr();
550  const auto fm = scaled.maxKrw();
551 
552  if (! (SwScaled > sr)) {
553  // Pure vertical scaling in left interval ([SWL, SR])
554  return unscaledKrw * (fr / fdisp);
555  }
556  else if (fmax > fdisp) {
557  // s \in [sr, sm), sm > sr; normal case: Kr(Smax) > Kr(Sr).
558  //
559  // Linear function between (sr,fr) and (sm,fm) in terms of
560  // function value 'unscaledKrw'. This usually alters the shape
561  // of the relative permeability function in this interval (e.g.,
562  // roughly quadratic to linear).
563  const auto t = (unscaledKrw - fdisp) / (fmax - fdisp);
564 
565  return fr + t*(fm - fr);
566  }
567  else if (sr < sm) {
568  // s \in [sr, sm), sm > sr; special case: Kr(Smax) == Kr(Sr).
569  //
570  // Linear function between (sr,fr) and (sm,fm) in terms of
571  // saturation value 'SwScaled'. This usually alters the shape
572  // of the relative permeability function in this interval (e.g.,
573  // roughly quadratic to linear).
574  const auto t = (SwScaled - sr) / (sm - sr);
575 
576  return fr + t*(fm - fr);
577  }
578  else {
579  // sm == sr (pure scaling). Almost arbitrarily pick 'fm'.
580  return fm;
581  }
582  }
583 
584  template <class Evaluation>
585  static Evaluation scaledToUnscaledKrw_(const Params& params, const Evaluation& scaledKrw)
586  {
587  if (!params.config().enableKrwScaling())
588  return scaledKrw;
589 
590  Scalar alpha = params.unscaledPoints().maxKrw()/params.scaledPoints().maxKrw();
591  return scaledKrw*alpha;
592  }
593 
597  template <class Evaluation>
598  static Evaluation unscaledToScaledKrn_(const Evaluation& SwScaled,
599  const Params& params,
600  const Evaluation& unscaledKrn)
601  {
602  const auto& cfg = params.config();
603 
604  if (! cfg.enableKrnScaling())
605  return unscaledKrn;
606 
607  const auto& scaled = params.scaledPoints();
608  const auto& unscaled = params.unscaledPoints();
609 
610  if (! cfg.enableThreePointKrnScaling()) {
611  // Simple case: Run uses pure vertical scaling of non-wetting
612  // phase's relative permeability (e.g., KRG)
613  const Scalar alpha = scaled.maxKrn() / unscaled.maxKrn();
614  return unscaledKrn * alpha;
615  }
616 
617  // Otherwise, run uses three-point vertical scaling (e.g., keywords KRGR and KRG)
618  const auto fdisp = unscaled.krnr();
619  const auto fmax = unscaled.maxKrn();
620 
621  const auto sl = scaled.saturationKrnPoints()[0];
622  const auto sr = std::max(scaled.saturationKrnPoints()[1], sl);
623  const auto fr = scaled.krnr();
624  const auto fm = scaled.maxKrn();
625 
626  // Note logic here. Krn is a decreasing function of Sw (dKrn/dSw <=
627  // 0) so the roles of left and right intervals are reversed viz
628  // unscaledToScaledKrw_().
629 
630  if (! (SwScaled < sr)) {
631  // Pure vertical scaling in right-hand interval ([SR, SWU])
632  return unscaledKrn * (fr / fdisp);
633  }
634  else if (fmax > fdisp) {
635  // s \in [SWL, SR), SWL < SR; normal case: Kr(Swl) > Kr(Sr).
636  //
637  // Linear function between (sr,fr) and (sl,fm) in terms of
638  // function value 'unscaledKrn'. This usually alters the shape
639  // of the relative permeability function in this interval (e.g.,
640  // roughly quadratic to linear).
641  const auto t = (unscaledKrn - fdisp) / (fmax - fdisp);
642 
643  return fr + t*(fm - fr);
644  }
645  else if (sr > sl) {
646  // s \in [SWL, SR), SWL < SR; special case: Kr(Swl) == Kr(Sr).
647  //
648  // Linear function between (sr,fr) and (sl,fm) in terms of
649  // saturation value 'SwScaled'. This usually alters the shape
650  // of the relative permeability function in this interval (e.g.,
651  // roughly quadratic to linear).
652  const auto t = (sr - SwScaled) / (sr - sl);
653 
654  return fr + t*(fm - fr);
655  }
656  else {
657  // sl == sr (pure scaling). Almost arbitrarily pick 'fm'.
658  return fm;
659  }
660  }
661 
662  template <class Evaluation>
663  static Evaluation scaledToUnscaledKrn_(const Params& params, const Evaluation& scaledKrn)
664  {
665  if (!params.config().enableKrnScaling())
666  return scaledKrn;
667 
668  Scalar alpha = params.unscaledPoints().maxKrn()/params.scaledPoints().maxKrn();
669  return scaledKrn*alpha;
670  }
671 };
672 } // namespace Opm
673 
674 #endif
A default implementation of the parameters for the material law adapter class which implements ECL en...
This material law takes a material law defined for unscaled saturation and converts it to a material ...
Definition: EclEpsTwoPhaseLaw.hpp:51
static Evaluation krw(const Params &, const FluidState &)
The relative permeability for the wetting phase.
Definition: EclEpsTwoPhaseLaw.hpp:216
static const int numPhases
The number of fluid phases.
Definition: EclEpsTwoPhaseLaw.hpp:63
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:338
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:179
static void saturations(Container &, const Params &, const FluidState &)
The saturation-capillary pressure curves.
Definition: EclEpsTwoPhaseLaw.hpp:169
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting phase.
Definition: EclEpsTwoPhaseLaw.hpp:241
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclEpsTwoPhaseLaw.hpp:70
static Evaluation scaledToUnscaledSatPc(const Params &params, const Evaluation &SwScaled)
Convert an absolute saturation to an effective one for capillary pressure.
Definition: EclEpsTwoPhaseLaw.hpp:268
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclEpsTwoPhaseLaw.hpp:86
static void relativePermeabilities(Container &, const Params &, const FluidState &)
The relative permeability-saturation curves depending on absolute saturations.
Definition: EclEpsTwoPhaseLaw.hpp:127
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 pcnw(const Params &, const FluidState &)
The capillary pressure-saturation curve.
Definition: EclEpsTwoPhaseLaw.hpp:144
static void capillaryPressures(Container &, const Params &, const FluidState &)
The capillary pressure-saturation curves depending on absolute saturations.
Definition: EclEpsTwoPhaseLaw.hpp:111
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:298
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:195
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclEpsTwoPhaseLaw.hpp:94
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 bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: EclEpsTwoPhaseLaw.hpp:90