EclStone2Material.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_STONE2_MATERIAL_HPP
28 #define OPM_ECL_STONE2_MATERIAL_HPP
29 
31 
32 #include <opm/common/Valgrind.hpp>
34 
35 #include <opm/common/Exceptions.hpp>
36 #include <opm/common/ErrorMacros.hpp>
37 
38 #include <algorithm>
39 
40 namespace Opm {
41 
55 template <class TraitsT,
56  class GasOilMaterialLawT,
57  class OilWaterMaterialLawT,
58  class ParamsT = EclStone2MaterialParams<TraitsT,
59  typename GasOilMaterialLawT::Params,
60  typename OilWaterMaterialLawT::Params> >
61 class EclStone2Material : public TraitsT
62 {
63 public:
64  typedef GasOilMaterialLawT GasOilMaterialLaw;
65  typedef OilWaterMaterialLawT OilWaterMaterialLaw;
66 
67  // some safety checks
68  static_assert(TraitsT::numPhases == 3,
69  "The number of phases considered by this capillary pressure "
70  "law is always three!");
71  static_assert(GasOilMaterialLaw::numPhases == 2,
72  "The number of phases considered by the gas-oil capillary "
73  "pressure law must be two!");
74  static_assert(OilWaterMaterialLaw::numPhases == 2,
75  "The number of phases considered by the oil-water capillary "
76  "pressure law must be two!");
77  static_assert(std::is_same<typename GasOilMaterialLaw::Scalar,
78  typename OilWaterMaterialLaw::Scalar>::value,
79  "The two two-phase capillary pressure laws must use the same "
80  "type of floating point values.");
81 
82  static_assert(GasOilMaterialLaw::implementsTwoPhaseSatApi,
83  "The gas-oil material law must implement the two-phase saturation "
84  "only API to for the default Ecl capillary pressure law!");
85  static_assert(OilWaterMaterialLaw::implementsTwoPhaseSatApi,
86  "The oil-water material law must implement the two-phase saturation "
87  "only API to for the default Ecl capillary pressure law!");
88 
89  typedef TraitsT Traits;
90  typedef ParamsT Params;
91  typedef typename Traits::Scalar Scalar;
92 
93  static const int numPhases = 3;
94  static const int waterPhaseIdx = Traits::wettingPhaseIdx;
95  static const int oilPhaseIdx = Traits::nonWettingPhaseIdx;
96  static const int gasPhaseIdx = Traits::gasPhaseIdx;
97 
100  static const bool implementsTwoPhaseApi = false;
101 
104  static const bool implementsTwoPhaseSatApi = false;
105 
108  static const bool isSaturationDependent = true;
109 
112  static const bool isPressureDependent = false;
113 
116  static const bool isTemperatureDependent = false;
117 
120  static const bool isCompositionDependent = false;
121 
136  template <class ContainerT, class FluidState>
137  static void capillaryPressures(ContainerT& values,
138  const Params& params,
139  const FluidState& state)
140  {
141  typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
142  values[gasPhaseIdx] = pcgn<FluidState, Evaluation>(params, state);
143  values[oilPhaseIdx] = 0;
144  values[waterPhaseIdx] = - pcnw<FluidState, Evaluation>(params, state);
145  Valgrind::CheckDefined(values[gasPhaseIdx]);
146  Valgrind::CheckDefined(values[oilPhaseIdx]);
147  Valgrind::CheckDefined(values[waterPhaseIdx]);
148  }
149 
150  /*
151  * Hysteresis parameters for oil-water
152  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
153  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
154  * \param params Parameters
155  */
156  static void oilWaterHysteresisParams(Scalar& pcSwMdc,
157  Scalar& krnSwMdc,
158  const Params& params)
159  {
160  pcSwMdc = params.oilWaterParams().pcSwMdc();
161  krnSwMdc = params.oilWaterParams().krnSwMdc();
162 
163  Valgrind::CheckDefined(pcSwMdc);
164  Valgrind::CheckDefined(krnSwMdc);
165  }
166 
167  /*
168  * Hysteresis parameters for oil-water
169  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
170  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
171  * \param params Parameters
172  */
173  static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
174  const Scalar& krnSwMdc,
175  Params& params)
176  {
177  const double krwSw = 2.0; //Should not be used
178  params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
179  }
180 
181  /*
182  * Hysteresis parameters for gas-oil
183  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
184  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
185  * \param params Parameters
186  */
187  static void gasOilHysteresisParams(Scalar& pcSwMdc,
188  Scalar& krnSwMdc,
189  const Params& params)
190  {
191  pcSwMdc = params.gasOilParams().pcSwMdc();
192  krnSwMdc = params.gasOilParams().krnSwMdc();
193 
194  Valgrind::CheckDefined(pcSwMdc);
195  Valgrind::CheckDefined(krnSwMdc);
196  }
197 
198  /*
199  * Hysteresis parameters for gas-oil
200  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
201  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
202  * \param params Parameters
203  */
204  static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
205  const Scalar& krnSwMdc,
206  Params& params)
207  {
208  const double krwSw = 2.0; //Should not be used
209  params.gasOilParams().update(pcSwMdc, krwSw, krnSwMdc);
210  }
211 
221  template <class FluidState, class Evaluation = typename FluidState::Scalar>
222  static Evaluation pcgn(const Params& params,
223  const FluidState& fs)
224  {
225  const auto& Sw = 1.0 - Opm::decay<Evaluation>(fs.saturation(gasPhaseIdx));
226  return GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(), Sw);
227  }
228 
238  template <class FluidState, class Evaluation = typename FluidState::Scalar>
239  static Evaluation pcnw(const Params& params,
240  const FluidState& fs)
241  {
242  const auto& Sw = Opm::decay<Evaluation>(fs.saturation(waterPhaseIdx));
243  Valgrind::CheckDefined(Sw);
244  const auto& result = OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(), Sw);
245  Valgrind::CheckDefined(result);
246  return result;
247  }
248 
252  template <class ContainerT, class FluidState>
253  static void saturations(ContainerT& /*values */,
254  const Params& /* params */,
255  const FluidState& /* fluidState */)
256  {
257  OPM_THROW(std::logic_error, "Not implemented: saturations()");
258  }
259 
263  template <class FluidState, class Evaluation = typename FluidState::Scalar>
264  static Evaluation Sg(const Params& /* params */,
265  const FluidState& /* fluidState */)
266  {
267  OPM_THROW(std::logic_error, "Not implemented: Sg()");
268  }
269 
273  template <class FluidState, class Evaluation = typename FluidState::Scalar>
274  static Evaluation Sn(const Params& /* params */,
275  const FluidState& /* fluidState */)
276  {
277  OPM_THROW(std::logic_error, "Not implemented: Sn()");
278  }
279 
283  template <class FluidState, class Evaluation = typename FluidState::Scalar>
284  static Evaluation Sw(const Params& /* params */,
285  const FluidState& /* fluidState */)
286  {
287  OPM_THROW(std::logic_error, "Not implemented: Sw()");
288  }
289 
305  template <class ContainerT, class FluidState>
306  static void relativePermeabilities(ContainerT& values,
307  const Params& params,
308  const FluidState& fluidState)
309  {
310  typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
311 
312  values[waterPhaseIdx] = krw<FluidState, Evaluation>(params, fluidState);
313  values[oilPhaseIdx] = krn<FluidState, Evaluation>(params, fluidState);
314  values[gasPhaseIdx] = krg<FluidState, Evaluation>(params, fluidState);
315  }
316 
320  template <class FluidState, class Evaluation = typename FluidState::Scalar>
321  static Evaluation krg(const Params& params,
322  const FluidState& fluidState)
323  {
324  const Evaluation& Sw = 1 - Opm::decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
325  return GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), Sw);
326  }
327 
331  template <class FluidState, class Evaluation = typename FluidState::Scalar>
332  static Evaluation krw(const Params& params,
333  const FluidState& fluidState)
334  {
335  const Evaluation& Sw = Opm::decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
336  return OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), Sw);
337  }
338 
342  template <class FluidState, class Evaluation = typename FluidState::Scalar>
343  static Evaluation krn(const Params& params,
344  const FluidState& fluidState)
345  {
346  Scalar Swco = params.Swl();
347  const Evaluation& Sw = Opm::decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
348  const Evaluation& Sg = Opm::decay<Evaluation>(fluidState.saturation(gasPhaseIdx));
349 
350  Scalar krocw = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), Swco);
351  Evaluation krow = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), Sw);
352  Evaluation krw = OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), Sw);
353  Evaluation krg = GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), 1 - Sg);
354  Evaluation krog = GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), 1 - Sg);
355 
356  return krocw*((krow/krocw + krw)*(krog/krocw + krg) - krw - krg);
357  }
358 
366  template <class FluidState>
367  static void updateHysteresis(Params& params, const FluidState& fluidState)
368  {
369  Scalar Sw = Opm::scalarValue(fluidState.saturation(waterPhaseIdx));
370  Scalar Sg = Opm::scalarValue(fluidState.saturation(gasPhaseIdx));
371 
372  params.oilWaterParams().update(/*pcSw=*/Sw, /*krwSw=*/Sw, /*krnSw=*/Sw);
373  params.gasOilParams().update(/*pcSw=*/1 - Sg, /*krwSw=*/1 - Sg, /*krnSw=*/1 - Sg);
374  }
375 };
376 } // namespace Opm
377 
378 #endif
static Evaluation krw(const Params &params, const FluidState &fluidState)
The relative permeability of the wetting phase.
Definition: EclStone2Material.hpp:332
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclStone2Material.hpp:367
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclStone2Material.hpp:306
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclStone2Material.hpp:100
Definition: Air_Mesitylene.hpp:33
static void capillaryPressures(ContainerT &values, const Params &params, const FluidState &state)
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclStone2Material.hpp:137
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclStone2Material.hpp:108
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: EclStone2Material.hpp:253
Default implementation for the parameters required by the three-phase capillary pressure/relperm Ston...
static Evaluation krn(const Params &params, const FluidState &fluidState)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition: EclStone2Material.hpp:343
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition: EclStone2Material.hpp:61
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition: EclStone2Material.hpp:264
static Evaluation krg(const Params &params, const FluidState &fluidState)
The relative permeability of the gas phase.
Definition: EclStone2Material.hpp:321
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclStone2Material.hpp:116
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition: EclStone2Material.hpp:274
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition...
Definition: EclStone2Material.hpp:120
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclStone2Material.hpp:104
static Evaluation pcgn(const Params &params, const FluidState &fs)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition: EclStone2Material.hpp:222
static Evaluation pcnw(const Params &params, const FluidState &fs)
Capillary pressure between the non-wetting liquid (i.e., oil) and the wetting liquid (i...
Definition: EclStone2Material.hpp:239
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition: EclStone2Material.hpp:284
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure...
Definition: EclStone2Material.hpp:112