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 #include <opm/common/ErrorMacros.hpp>
34 #include <opm/common/Exceptions.hpp>
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  OPM_THROW(NotImplemented,
114  "The capillaryPressures(fs) method is not yet implemented");
115  }
116 
127  template <class Container, class FluidState>
128  static void relativePermeabilities(Container& /*values*/, const Params& /*params*/, const FluidState& /*fluidState*/)
129  {
130  OPM_THROW(NotImplemented,
131  "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  OPM_THROW(NotImplemented,
149  "The pcnw(fs) method is not yet implemented");
150  }
151 
152  template <class Evaluation>
153  static Evaluation twoPhaseSatPcnw(const Params& params, const Evaluation& SwScaled)
154  {
155  const Evaluation& SwUnscaled = scaledToUnscaledSatPc(params, SwScaled);
156  const Evaluation& pcUnscaled = EffLaw::twoPhaseSatPcnw(params.effectiveLawParams(), SwUnscaled);
157  return unscaledToScaledPcnw_(params, pcUnscaled);
158  }
159 
160  template <class Evaluation>
161  static Evaluation twoPhaseSatPcnwInv(const Params& params, const Evaluation& pcnwScaled)
162  {
163  Evaluation pcnwUnscaled = scaledToUnscaledPcnw_(params, pcnwScaled);
164  Evaluation SwUnscaled = EffLaw::twoPhaseSatPcnwInv(params.effectiveLawParams(), pcnwUnscaled);
165  return unscaledToScaledSatPc(params, SwUnscaled);
166  }
167 
171  template <class Container, class FluidState>
172  static void saturations(Container& /*values*/, const Params& /*params*/, const FluidState& /*fluidState*/)
173  {
174  OPM_THROW(NotImplemented,
175  "The saturations(fs) method is not yet implemented");
176  }
177 
182  template <class FluidState, class Evaluation = typename FluidState::Scalar>
183  static Evaluation Sw(const Params& /*params*/, const FluidState& /*fluidState*/)
184  {
185  OPM_THROW(NotImplemented,
186  "The Sw(fs) method is not yet implemented");
187  }
188 
189  template <class Evaluation>
190  static Evaluation twoPhaseSatSw(const Params& /*params*/, const Evaluation& /*pc*/)
191  {
192  OPM_THROW(NotImplemented,
193  "The twoPhaseSatSw(pc) method is not yet implemented");
194  }
195 
200  template <class FluidState, class Evaluation = typename FluidState::Scalar>
201  static Evaluation Sn(const Params& /*params*/, const FluidState& /*fluidState*/)
202  {
203  OPM_THROW(NotImplemented,
204  "The Sn(pc) method is not yet implemented");
205  }
206 
207  template <class Evaluation>
208  static Evaluation twoPhaseSatSn(const Params& /*params*/, const Evaluation& /*pc*/)
209  {
210  OPM_THROW(NotImplemented,
211  "The twoPhaseSatSn(pc) method is not yet implemented");
212  }
213 
223  template <class FluidState, class Evaluation = typename FluidState::Scalar>
224  static Evaluation krw(const Params& /*params*/, const FluidState& /*fluidState*/)
225  {
226  OPM_THROW(NotImplemented,
227  "The krw(fs) method is not yet implemented");
228  }
229 
230  template <class Evaluation>
231  static Evaluation twoPhaseSatKrw(const Params& params, const Evaluation& SwScaled)
232  {
233  const Evaluation& SwUnscaled = scaledToUnscaledSatKrw(params, SwScaled);
234  const Evaluation& krwUnscaled = EffLaw::twoPhaseSatKrw(params.effectiveLawParams(), SwUnscaled);
235  return unscaledToScaledKrw_(params, krwUnscaled);
236  }
237 
238  template <class Evaluation>
239  static Evaluation twoPhaseSatKrwInv(const Params& params, const Evaluation& krwScaled)
240  {
241  Evaluation krwUnscaled = scaledToUnscaledKrw_(params, krwScaled);
242  Evaluation SwUnscaled = EffLaw::twoPhaseSatKrwInv(params.effectiveLawParams(), krwUnscaled);
243  return unscaledToScaledSatKrw(params, SwUnscaled);
244  }
245 
249  template <class FluidState, class Evaluation = typename FluidState::Scalar>
250  static Evaluation krn(const Params& /*params*/, const FluidState& /*fluidState*/)
251  {
252  OPM_THROW(NotImplemented,
253  "The krn(fs) method is not yet implemented");
254  }
255 
256  template <class Evaluation>
257  static Evaluation twoPhaseSatKrn(const Params& params, const Evaluation& SwScaled)
258  {
259  const Evaluation& SwUnscaled = scaledToUnscaledSatKrn(params, SwScaled);
260  const Evaluation& krnUnscaled = EffLaw::twoPhaseSatKrn(params.effectiveLawParams(), SwUnscaled);
261  return unscaledToScaledKrn_(params, krnUnscaled);
262  }
263 
264  template <class Evaluation>
265  static Evaluation twoPhaseSatKrnInv(const Params& params, const Evaluation& krnScaled)
266  {
267  Evaluation krnUnscaled = scaledToUnscaledKrn_(params, krnScaled);
268  Evaluation SwUnscaled = EffLaw::twoPhaseSatKrnInv(params.effectiveLawParams(), krnUnscaled);
269  return unscaledToScaledSatKrn(params, SwUnscaled);
270  }
271 
277  template <class Evaluation>
278  static Evaluation scaledToUnscaledSatPc(const Params& params, const Evaluation& SwScaled)
279  {
280  if (!params.config().enableSatScaling())
281  return SwScaled;
282 
283  // the saturations of capillary pressure are always scaled using two-point
284  // scaling
285  return scaledToUnscaledSatTwoPoint_(SwScaled,
286  params.unscaledPoints().saturationPcPoints(),
287  params.scaledPoints().saturationPcPoints());
288  }
289 
290  template <class Evaluation>
291  static Evaluation unscaledToScaledSatPc(const Params& params, const Evaluation& SwUnscaled)
292  {
293  if (!params.config().enableSatScaling())
294  return SwUnscaled;
295 
296  // the saturations of capillary pressure are always scaled using two-point
297  // scaling
298  return unscaledToScaledSatTwoPoint_(SwUnscaled,
299  params.unscaledPoints().saturationPcPoints(),
300  params.scaledPoints().saturationPcPoints());
301  }
302 
307  template <class Evaluation>
308  static Evaluation scaledToUnscaledSatKrw(const Params& params, const Evaluation& SwScaled)
309  {
310  if (!params.config().enableSatScaling())
311  return SwScaled;
312 
313  if (params.config().enableThreePointKrSatScaling()) {
314  return scaledToUnscaledSatThreePoint_(SwScaled,
315  params.unscaledPoints().saturationKrwPoints(),
316  params.scaledPoints().saturationKrwPoints());
317  }
318  else { // two-point relperm saturation scaling
319  return scaledToUnscaledSatTwoPoint_(SwScaled,
320  params.unscaledPoints().saturationKrwPoints(),
321  params.scaledPoints().saturationKrwPoints());
322  }
323  }
324 
325  template <class Evaluation>
326  static Evaluation unscaledToScaledSatKrw(const Params& params, const Evaluation& SwUnscaled)
327  {
328  if (!params.config().enableSatScaling())
329  return SwUnscaled;
330 
331  if (params.config().enableThreePointKrSatScaling()) {
332  return unscaledToScaledSatThreePoint_(SwUnscaled,
333  params.unscaledPoints().saturationKrwPoints(),
334  params.scaledPoints().saturationKrwPoints());
335  }
336  else { // two-point relperm saturation scaling
337  return unscaledToScaledSatTwoPoint_(SwUnscaled,
338  params.unscaledPoints().saturationKrwPoints(),
339  params.scaledPoints().saturationKrwPoints());
340  }
341  }
342 
347  template <class Evaluation>
348  static Evaluation scaledToUnscaledSatKrn(const Params& params, const Evaluation& SwScaled)
349  {
350  if (!params.config().enableSatScaling())
351  return SwScaled;
352 
353  if (params.config().enableThreePointKrSatScaling())
354  return scaledToUnscaledSatThreePoint_(SwScaled,
355  params.unscaledPoints().saturationKrnPoints(),
356  params.scaledPoints().saturationKrnPoints());
357  else // two-point relperm saturation scaling
358  return scaledToUnscaledSatTwoPoint_(SwScaled,
359  params.unscaledPoints().saturationKrnPoints(),
360  params.scaledPoints().saturationKrnPoints());
361  }
362 
363 
364  template <class Evaluation>
365  static Evaluation unscaledToScaledSatKrn(const Params& params, const Evaluation& SwUnscaled)
366  {
367  if (!params.config().enableSatScaling())
368  return SwUnscaled;
369 
370  if (params.config().enableThreePointKrSatScaling()) {
371  return unscaledToScaledSatThreePoint_(SwUnscaled,
372  params.unscaledPoints().saturationKrnPoints(),
373  params.scaledPoints().saturationKrnPoints());
374  }
375  else { // two-point relperm saturation scaling
376  return unscaledToScaledSatTwoPoint_(SwUnscaled,
377  params.unscaledPoints().saturationKrnPoints(),
378  params.scaledPoints().saturationKrnPoints());
379  }
380  }
381 
382 private:
383  template <class Evaluation, class PointsContainer>
384  static Evaluation scaledToUnscaledSatTwoPoint_(const Evaluation& scaledSat,
385  const PointsContainer& unscaledSats,
386  const PointsContainer& scaledSats)
387  {
388  return
389  unscaledSats[0]
390  +
391  (scaledSat - scaledSats[0])*((unscaledSats[1] - unscaledSats[0])/(scaledSats[1] - scaledSats[0]));
392  }
393 
394  template <class Evaluation, class PointsContainer>
395  static Evaluation unscaledToScaledSatTwoPoint_(const Evaluation& unscaledSat,
396  const PointsContainer& unscaledSats,
397  const PointsContainer& scaledSats)
398  {
399  return
400  scaledSats[0]
401  +
402  (unscaledSat - unscaledSats[0])*((scaledSats[1] - scaledSats[0])/(unscaledSats[1] - unscaledSats[0]));
403  }
404 
405  template <class Evaluation, class PointsContainer>
406  static Evaluation scaledToUnscaledSatThreePoint_(const Evaluation& scaledSat,
407  const PointsContainer& unscaledSats,
408  const PointsContainer& scaledSats)
409  {
410  if (unscaledSats[1] >= unscaledSats[2])
411  return scaledToUnscaledSatTwoPoint_(scaledSat, unscaledSats, scaledSats);
412 
413  if (scaledSat < scaledSats[1]) {
414  Scalar delta = scaledSats[1] - scaledSats[0];
415  if (delta <= 1e-20)
416  delta = 1.0; // prevent division by zero for (possibly) incorrect input data
417 
418  return
419  unscaledSats[0]
420  +
421  (scaledSat - scaledSats[0])*((unscaledSats[1] - unscaledSats[0])/delta);
422  }
423  else {
424  Scalar delta = scaledSats[2] - scaledSats[1];
425  if (delta <= 1e-20)
426  delta = 1.0; // prevent division by zero for (possibly) incorrect input data
427 
428  return
429  unscaledSats[1]
430  +
431  (scaledSat - scaledSats[1])*((unscaledSats[2] - unscaledSats[1])/delta);
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  if (unscaledSats[1] >= unscaledSats[2])
441  return unscaledToScaledSatTwoPoint_(unscaledSat, unscaledSats, scaledSats);
442 
443  if (unscaledSat < unscaledSats[1]) {
444  Scalar delta = unscaledSats[1] - unscaledSats[0];
445  if (delta <= 1e-20)
446  delta = 1.0; // prevent division by zero for (possibly) incorrect input data
447 
448  return
449  scaledSats[0]
450  +
451  (unscaledSat - unscaledSats[0])*((scaledSats[1] - scaledSats[0])/delta);
452  }
453  else {
454  Scalar delta = unscaledSats[2] - unscaledSats[1];
455  if (delta <= 1e-20)
456  delta = 1.0; // prevent division by zero for (possibly) incorrect input data
457 
458  return
459  scaledSats[1]
460  +
461  (unscaledSat - unscaledSats[1])*((scaledSats[2] - scaledSats[1])/delta);
462  }
463  }
464 
468  template <class Evaluation>
469  static Evaluation unscaledToScaledPcnw_(const Params& params, const Evaluation& unscaledPcnw)
470  {
471  if (params.config().enableLeverettScaling()) {
472  Scalar alpha = params.scaledPoints().leverettFactor();
473  return unscaledPcnw*alpha;
474  }
475  else if (params.config().enablePcScaling()) {
476  Scalar alpha = params.scaledPoints().maxPcnw()/params.unscaledPoints().maxPcnw();
477  return unscaledPcnw*alpha;
478  }
479 
480  return unscaledPcnw;
481  }
482 
483  template <class Evaluation>
484  static Evaluation scaledToUnscaledPcnw_(const Params& params, const Evaluation& scaledPcnw)
485  {
486  if (params.config().enableLeverettScaling()) {
487  Scalar alpha = params.scaledPoints().leverettFactor();
488  return scaledPcnw/alpha;
489  }
490  else if (params.config().enablePcScaling()) {
491  Scalar alpha = params.unscaledPoints().maxPcnw()/params.scaledPoints().maxPcnw();
492  return scaledPcnw/alpha;
493  }
494 
495  return scaledPcnw;
496  }
497 
501  template <class Evaluation>
502  static Evaluation unscaledToScaledKrw_(const Params& params, const Evaluation& unscaledKrw)
503  {
504  if (!params.config().enableKrwScaling())
505  return unscaledKrw;
506 
507  // TODO: three point krw y-scaling
508  Scalar alpha = params.scaledPoints().maxKrw()/params.unscaledPoints().maxKrw();
509  return unscaledKrw*alpha;
510  }
511 
512  template <class Evaluation>
513  static Evaluation scaledToUnscaledKrw_(const Params& params, const Evaluation& scaledKrw)
514  {
515  if (!params.config().enableKrwScaling())
516  return scaledKrw;
517 
518  Scalar alpha = params.unscaledPoints().maxKrw()/params.scaledPoints().maxKrw();
519  return scaledKrw*alpha;
520  }
521 
525  template <class Evaluation>
526  static Evaluation unscaledToScaledKrn_(const Params& params, const Evaluation& unscaledKrn)
527  {
528  if (!params.config().enableKrnScaling())
529  return unscaledKrn;
530 
531  //TODO: three point krn y-scaling
532  Scalar alpha = params.scaledPoints().maxKrn()/params.unscaledPoints().maxKrn();
533  return unscaledKrn*alpha;
534  }
535 
536  template <class Evaluation>
537  static Evaluation scaledToUnscaledKrn_(const Params& params, const Evaluation& scaledKrn)
538  {
539  if (!params.config().enableKrnScaling())
540  return scaledKrn;
541 
542  Scalar alpha = params.unscaledPoints().maxKrn()/params.scaledPoints().maxKrn();
543  return scaledKrn*alpha;
544  }
545 };
546 } // namespace Opm
547 
548 #endif
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 &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:348
A default implementation of the parameters for the material law adapter class which implements ECL en...
Definition: Air_Mesitylene.hpp:33
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 &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: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 &params, const Evaluation &SwScaled)
Convert an absolute saturation to an effective one for capillary pressure.
Definition: EclEpsTwoPhaseLaw.hpp:278