VanGenuchten.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_VAN_GENUCHTEN_HPP
28 #define OPM_VAN_GENUCHTEN_HPP
29 
30 #include "VanGenuchtenParams.hpp"
31 
33 
34 #include <algorithm>
35 #include <cmath>
36 #include <cassert>
37 
38 namespace Opm {
54 template <class TraitsT, class ParamsT = VanGenuchtenParams<TraitsT> >
55 class VanGenuchten : public TraitsT
56 {
57 public:
59  typedef TraitsT Traits;
60 
62  typedef ParamsT Params;
63 
65  typedef typename Traits::Scalar Scalar;
66 
68  static const int numPhases = Traits::numPhases;
69  static_assert(numPhases == 2,
70  "The van Genuchten capillary pressure law only "
71  "applies to the case of two fluid phases");
72 
75  static const bool implementsTwoPhaseApi = true;
76 
79  static const bool implementsTwoPhaseSatApi = true;
80 
83  static const bool isSaturationDependent = true;
84 
87  static const bool isPressureDependent = false;
88 
91  static const bool isTemperatureDependent = false;
92 
95  static const bool isCompositionDependent = false;
96 
113  template <class Container, class FluidState>
114  static void capillaryPressures(Container& values, const Params& params, const FluidState& fs)
115  {
116  typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
117 
118  values[Traits::wettingPhaseIdx] = 0.0; // reference phase
119  values[Traits::nonWettingPhaseIdx] = pcnw<FluidState, Evaluation>(params, fs);
120  }
121 
126  template <class Container, class FluidState>
127  static void saturations(Container& values, const Params& params, const FluidState& fs)
128  {
129  typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
130 
131  values[Traits::wettingPhaseIdx] = Sw<FluidState, Evaluation>(params, fs);
132  values[Traits::nonWettingPhaseIdx] = 1 - values[Traits::wettingPhaseIdx];
133  }
134 
145  template <class Container, class FluidState>
146  static void relativePermeabilities(Container& values, const Params& params, const FluidState& fs)
147  {
148  typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
149 
150  values[Traits::wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fs);
151  values[Traits::nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fs);
152  }
153 
168  template <class FluidState, class Evaluation = typename FluidState::Scalar>
169  static Evaluation pcnw(const Params& params, const FluidState& fs)
170  {
171  const Evaluation& Sw =
172  Opm::decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
173 
174  assert(0 <= Sw && Sw <= 1);
175 
176  return twoPhaseSatPcnw(params, Sw);
177  }
178 
193  template <class Evaluation>
194  static Evaluation twoPhaseSatPcnw(const Params& params, const Evaluation& Sw)
195  {
196  return Opm::pow(Opm::pow(Sw, -1.0/params.vgM()) - 1, 1.0/params.vgN())/params.vgAlpha();
197  }
198 
211  template <class FluidState, class Evaluation = typename FluidState::Scalar>
212  static Evaluation Sw(const Params& params, const FluidState& fs)
213  {
214  Evaluation pC =
215  Opm::decay<Evaluation>(fs.pressure(Traits::nonWettingPhaseIdx))
216  - Opm::decay<Evaluation>(fs.pressure(Traits::wettingPhaseIdx));
217  return twoPhaseSatSw(params, pC);
218  }
219 
220  template <class Evaluation>
221  static Evaluation twoPhaseSatSw(const Params& params, const Evaluation& pC)
222  {
223  assert(pC >= 0);
224 
225  return Opm::pow(Opm::pow(params.vgAlpha()*pC, params.vgN()) + 1, -params.vgM());
226  }
227 
232  template <class FluidState, class Evaluation = typename FluidState::Scalar>
233  static Evaluation Sn(const Params& params, const FluidState& fs)
234  { return 1 - Sw<FluidState, Evaluation>(params, fs); }
235 
236  template <class Evaluation>
237  static Evaluation twoPhaseSatSn(const Params& params, const Evaluation& pC)
238  { return 1 - twoPhaseSatSw(params, pC); }
239 
250  template <class FluidState, class Evaluation = typename FluidState::Scalar>
251  static Evaluation krw(const Params& params, const FluidState& fs)
252  {
253  const Evaluation& Sw =
254  Opm::decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
255 
256  return twoPhaseSatKrw(params, Sw);
257  }
258 
259  template <class Evaluation>
260  static Evaluation twoPhaseSatKrw(const Params& params, const Evaluation& Sw)
261  {
262  assert(0.0 <= Sw && Sw <= 1.0);
263 
264  Evaluation r = 1.0 - Opm::pow(1.0 - Opm::pow(Sw, 1/params.vgM()), params.vgM());
265  return Opm::sqrt(Sw)*r*r;
266  }
267 
277  template <class FluidState, class Evaluation = typename FluidState::Scalar>
278  static Evaluation krn(const Params& params, const FluidState& fs)
279  {
280  const Evaluation& Sw =
281  1.0 - Opm::decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
282 
283  return twoPhaseSatKrn(params, Sw);
284  }
285 
286  template <class Evaluation>
287  static Evaluation twoPhaseSatKrn(const Params& params, Evaluation Sw)
288  {
289  assert(0 <= Sw && Sw <= 1);
290 
291  return
292  Opm::pow(1 - Sw, 1.0/3) *
293  Opm::pow(1 - Opm::pow(Sw, 1/params.vgM()), 2*params.vgM());
294  }
295 };
296 } // namespace Opm
297 
298 #endif
static Evaluation krn(const Params &params, const FluidState &fs)
The relative permeability for the non-wetting phase of the medium according to van Genuchten...
Definition: VanGenuchten.hpp:278
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: VanGenuchten.hpp:91
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: VanGenuchten.hpp:83
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
TraitsT Traits
The traits class for this material law.
Definition: VanGenuchten.hpp:59
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition...
Definition: VanGenuchten.hpp:95
ParamsT Params
The type of the parameter objects for this law.
Definition: VanGenuchten.hpp:62
static Evaluation krw(const Params &params, const FluidState &fs)
The relative permeability for the wetting phase of the medium according to van Genuchten&#39;s curve with...
Definition: VanGenuchten.hpp:251
static Evaluation Sw(const Params &params, const FluidState &fs)
The saturation-capillary pressure curve according to van Genuchten.
Definition: VanGenuchten.hpp:212
static void relativePermeabilities(Container &values, const Params &params, const FluidState &fs)
The relative permeability-saturation curves according to van Genuchten.
Definition: VanGenuchten.hpp:146
static void saturations(Container &values, const Params &params, const FluidState &fs)
Calculate the saturations of the phases starting from their pressure differences. ...
Definition: VanGenuchten.hpp:127
static Evaluation pcnw(const Params &params, const FluidState &fs)
The capillary pressure-saturation curve according to van Genuchten.
Definition: VanGenuchten.hpp:169
Implementation of the van Genuchten capillary pressure - saturation relation.
Definition: VanGenuchten.hpp:55
Definition: Air_Mesitylene.hpp:33
Specification of the material parameters for the van Genuchten constitutive relations.
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: VanGenuchten.hpp:79
static Evaluation Sn(const Params &params, const FluidState &fs)
Calculate the non-wetting phase saturations depending on the phase pressures.
Definition: VanGenuchten.hpp:233
static Evaluation twoPhaseSatPcnw(const Params &params, const Evaluation &Sw)
The saturation-capillary pressure curve according to van Genuchten using a material law specific API...
Definition: VanGenuchten.hpp:194
static void capillaryPressures(Container &values, const Params &params, const FluidState &fs)
The capillary pressure-saturation curves according to van Genuchten.
Definition: VanGenuchten.hpp:114
static const int numPhases
The number of fluid phases.
Definition: VanGenuchten.hpp:68
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: VanGenuchten.hpp:75
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure...
Definition: VanGenuchten.hpp:87
Traits::Scalar Scalar
The type of the scalar values for this law.
Definition: VanGenuchten.hpp:65