My Project
EclTwoPhaseMaterial.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_TWO_PHASE_MATERIAL_HPP
28 #define OPM_ECL_TWO_PHASE_MATERIAL_HPP
29 
31 
34 
35 #include <algorithm>
36 
37 namespace Opm {
38 
48 template <class TraitsT,
49  class GasOilMaterialLawT,
50  class OilWaterMaterialLawT,
51  class GasWaterMaterialLawT,
52  class ParamsT = EclTwoPhaseMaterialParams<TraitsT,
53  typename GasOilMaterialLawT::Params,
54  typename OilWaterMaterialLawT::Params,
55  typename GasWaterMaterialLawT::Params> >
56 class EclTwoPhaseMaterial : public TraitsT
57 {
58 public:
59  using GasOilMaterialLaw = GasOilMaterialLawT;
60  using OilWaterMaterialLaw = OilWaterMaterialLawT;
61  using GasWaterMaterialLaw = GasWaterMaterialLawT;
62 
63  // some safety checks
64  static_assert(TraitsT::numPhases == 3,
65  "The number of phases considered by this capillary pressure "
66  "law is always three!");
67  static_assert(GasOilMaterialLaw::numPhases == 2,
68  "The number of phases considered by the gas-oil capillary "
69  "pressure law must be two!");
70  static_assert(OilWaterMaterialLaw::numPhases == 2,
71  "The number of phases considered by the oil-water capillary "
72  "pressure law must be two!");
73  static_assert(GasWaterMaterialLaw::numPhases == 2,
74  "The number of phases considered by the gas-water capillary "
75  "pressure law must be two!");
76  static_assert(std::is_same<typename GasOilMaterialLaw::Scalar,
77  typename OilWaterMaterialLaw::Scalar>::value,
78  "The two two-phase capillary pressure laws must use the same "
79  "type of floating point values.");
80 
81  using Traits = TraitsT;
82  using Params = ParamsT;
83  using Scalar = typename Traits::Scalar;
84 
85  static constexpr int numPhases = 3;
86  static constexpr int waterPhaseIdx = Traits::wettingPhaseIdx;
87  static constexpr int oilPhaseIdx = Traits::nonWettingPhaseIdx;
88  static constexpr int gasPhaseIdx = Traits::gasPhaseIdx;
89 
92  static constexpr bool implementsTwoPhaseApi = false;
93 
96  static constexpr bool implementsTwoPhaseSatApi = false;
97 
100  static constexpr bool isSaturationDependent = true;
101 
104  static constexpr bool isPressureDependent = false;
105 
108  static constexpr bool isTemperatureDependent = false;
109 
112  static constexpr bool isCompositionDependent = false;
113 
128  template <class ContainerT, class FluidState>
129  static void capillaryPressures(ContainerT& values,
130  const Params& params,
131  const FluidState& fluidState)
132  {
133  using Evaluation = typename std::remove_reference<decltype(values[0])>::type;
134 
135  switch (params.approach()) {
136  case EclTwoPhaseApproach::EclTwoPhaseGasOil: {
137  const Evaluation& So =
138  decay<Evaluation>(fluidState.saturation(oilPhaseIdx));
139 
140  values[oilPhaseIdx] = 0.0;
141  values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(), So);
142  break;
143  }
144 
145  case EclTwoPhaseApproach::EclTwoPhaseOilWater: {
146  const Evaluation& Sw =
147  decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
148 
149  values[waterPhaseIdx] = 0.0;
150  values[oilPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(), Sw);
151  break;
152  }
153 
154  case EclTwoPhaseApproach::EclTwoPhaseGasWater: {
155  const Evaluation& Sw =
156  decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
157 
158  values[waterPhaseIdx] = 0.0;
159  values[gasPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatPcnw(params.gasWaterParams(), Sw);
160  break;
161  }
162 
163  }
164  }
165 
166  /*
167  * Hysteresis parameters for oil-water
168  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
169  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
170  * \param params Parameters
171  */
172  static void oilWaterHysteresisParams(Scalar& pcSwMdc,
173  Scalar& krnSwMdc,
174  const Params& params)
175  {
176  pcSwMdc = params.oilWaterParams().pcSwMdc();
177  krnSwMdc = params.oilWaterParams().krnSwMdc();
178 
179  Valgrind::CheckDefined(pcSwMdc);
180  Valgrind::CheckDefined(krnSwMdc);
181  }
182 
183  /*
184  * Hysteresis parameters for oil-water
185  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
186  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
187  * \param params Parameters
188  */
189  static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
190  const Scalar& krnSwMdc,
191  Params& params)
192  {
193  constexpr const Scalar krwSw = 2.0; //Should not be used
194  params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
195  }
196 
197  /*
198  * Hysteresis parameters for gas-oil
199  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
200  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
201  * \param params Parameters
202  */
203  static void gasOilHysteresisParams(Scalar& pcSwMdc,
204  Scalar& krnSwMdc,
205  const Params& params)
206  {
207  pcSwMdc = params.gasOilParams().pcSwMdc();
208  krnSwMdc = params.gasOilParams().krnSwMdc();
209 
210  Valgrind::CheckDefined(pcSwMdc);
211  Valgrind::CheckDefined(krnSwMdc);
212  }
213 
214  /*
215  * Hysteresis parameters for gas-oil
216  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
217  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
218  * \param params Parameters
219  */
220  static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
221  const Scalar& krnSwMdc,
222  Params& params)
223  {
224  constexpr const Scalar krwSw = 2.0; //Should not be used
225  params.gasOilParams().update(pcSwMdc, krwSw, krnSwMdc);
226  }
227 
237  template <class FluidState, class Evaluation = typename FluidState::Scalar>
238  static Evaluation pcgn(const Params& /* params */,
239  const FluidState& /* fs */)
240  {
241  throw std::logic_error("Not implemented: pcgn()");
242  }
243 
253  template <class FluidState, class Evaluation = typename FluidState::Scalar>
254  static Evaluation pcnw(const Params& /* params */,
255  const FluidState& /* fs */)
256  {
257  throw std::logic_error("Not implemented: pcnw()");
258  }
259 
263  template <class ContainerT, class FluidState>
264  static void saturations(ContainerT& /* values */,
265  const Params& /* params */,
266  const FluidState& /* fs */)
267  {
268  throw std::logic_error("Not implemented: saturations()");
269  }
270 
274  template <class FluidState, class Evaluation = typename FluidState::Scalar>
275  static Evaluation Sg(const Params& /* params */,
276  const FluidState& /* fluidState */)
277  {
278  throw std::logic_error("Not implemented: Sg()");
279  }
280 
284  template <class FluidState, class Evaluation = typename FluidState::Scalar>
285  static Evaluation Sn(const Params& /* params */,
286  const FluidState& /* fluidState */)
287  {
288  throw std::logic_error("Not implemented: Sn()");
289  }
290 
294  template <class FluidState, class Evaluation = typename FluidState::Scalar>
295  static Evaluation Sw(const Params& /* params */,
296  const FluidState& /* fluidState */)
297  {
298  throw std::logic_error("Not implemented: Sw()");
299  }
300 
316  template <class ContainerT, class FluidState>
317  static void relativePermeabilities(ContainerT& values,
318  const Params& params,
319  const FluidState& fluidState)
320  {
321  using Evaluation = typename std::remove_reference<decltype(values[0])>::type;
322 
323  switch (params.approach()) {
324  case EclTwoPhaseApproach::EclTwoPhaseGasOil: {
325  const Evaluation& So =
326  decay<Evaluation>(fluidState.saturation(oilPhaseIdx));
327 
328  values[oilPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), So);
329  values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), So);
330  break;
331  }
332 
333  case EclTwoPhaseApproach::EclTwoPhaseOilWater: {
334  const Evaluation& Sw =
335  decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
336 
337  values[waterPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), Sw);
338  values[oilPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), Sw);
339  break;
340  }
341 
342  case EclTwoPhaseApproach::EclTwoPhaseGasWater: {
343  const Evaluation& Sw =
344  decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
345 
346  values[waterPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatKrw(params.gasWaterParams(), Sw);
347  values[gasPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatKrn(params.gasWaterParams(), Sw);
348 
349  break;
350  }
351  }
352  }
353 
357  template <class FluidState, class Evaluation = typename FluidState::Scalar>
358  static Evaluation krg(const Params& /* params */,
359  const FluidState& /* fluidState */)
360  {
361  throw std::logic_error("Not implemented: krg()");
362  }
363 
367  template <class FluidState, class Evaluation = typename FluidState::Scalar>
368  static Evaluation krw(const Params& /* params */,
369  const FluidState& /* fluidState */)
370  {
371  throw std::logic_error("Not implemented: krw()");
372  }
373 
377  template <class FluidState, class Evaluation = typename FluidState::Scalar>
378  static Evaluation krn(const Params& /* params */,
379  const FluidState& /* fluidState */)
380  {
381  throw std::logic_error("Not implemented: krn()");
382  }
383 
384 
392  template <class FluidState>
393  static void updateHysteresis(Params& params, const FluidState& fluidState)
394  {
395  switch (params.approach()) {
396  case EclTwoPhaseApproach::EclTwoPhaseGasOil: {
397  Scalar So = scalarValue(fluidState.saturation(oilPhaseIdx));
398 
399  params.gasOilParams().update(/*pcSw=*/So, /*krwSw=*/So, /*krnSw=*/So);
400  break;
401  }
402 
403  case EclTwoPhaseApproach::EclTwoPhaseOilWater: {
404  Scalar Sw = scalarValue(fluidState.saturation(waterPhaseIdx));
405 
406  params.oilWaterParams().update(/*pcSw=*/Sw, /*krwSw=*/Sw, /*krnSw=*/Sw);
407  break;
408  }
409 
410  case EclTwoPhaseApproach::EclTwoPhaseGasWater: {
411  Scalar Sw = scalarValue(fluidState.saturation(waterPhaseIdx));
412 
413  params.gasWaterParams().update(/*pcSw=*/1.0, /*krwSw=*/0.0, /*krnSw=*/Sw);
414  break;
415  }
416  }
417  }
418 };
419 
420 } // namespace Opm
421 
422 #endif
Implementation for the parameters required by the material law for two-phase simulations.
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
Some templates to wrap the valgrind client request macros.
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Definition: EclTwoPhaseMaterial.hpp:57
static constexpr bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclTwoPhaseMaterial.hpp:108
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclTwoPhaseMaterial.hpp:317
static constexpr bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclTwoPhaseMaterial.hpp:92
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition: EclTwoPhaseMaterial.hpp:295
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition: EclTwoPhaseMaterial.hpp:238
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: EclTwoPhaseMaterial.hpp:264
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition: EclTwoPhaseMaterial.hpp:378
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclTwoPhaseMaterial.hpp:393
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition: EclTwoPhaseMaterial.hpp:285
static constexpr bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclTwoPhaseMaterial.hpp:100
static constexpr bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: EclTwoPhaseMaterial.hpp:104
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition: EclTwoPhaseMaterial.hpp:358
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition: EclTwoPhaseMaterial.hpp:275
static void capillaryPressures(ContainerT &values, const Params &params, const FluidState &fluidState)
Implements the multiplexer three phase capillary pressure law used by the ECLipse simulator.
Definition: EclTwoPhaseMaterial.hpp:129
static constexpr bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: EclTwoPhaseMaterial.hpp:112
static Evaluation pcnw(const Params &, const FluidState &)
Capillary pressure between the non-wetting liquid (i.e., oil) and the wetting liquid (i....
Definition: EclTwoPhaseMaterial.hpp:254
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition: EclTwoPhaseMaterial.hpp:368
static constexpr bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclTwoPhaseMaterial.hpp:96