My Project
TwoPhaseImmiscibleFluidSystem.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_TWO_PHASE_IMMISCIBLE_FLUID_SYSTEM_HPP
28 #define OPM_TWO_PHASE_IMMISCIBLE_FLUID_SYSTEM_HPP
29 
30 #include <limits>
31 #include <cassert>
32 
36 
37 #include "BaseFluidSystem.hpp"
38 #include "NullParameterCache.hpp"
39 
40 namespace Opm {
41 
55 template <class Scalar, class WettingPhase, class NonwettingPhase>
57  : public BaseFluidSystem<Scalar, TwoPhaseImmiscibleFluidSystem<Scalar, WettingPhase, NonwettingPhase> >
58 {
59  // do not try to instanciate this class, it has only static members!
61  {}
62 
65 
66 public:
67  template <class Evaluation>
68  struct ParameterCache : public NullParameterCache<Evaluation>
69  {};
70 
71  /****************************************
72  * Fluid phase related static parameters
73  ****************************************/
74 
76  static const int numPhases = 2;
77 
79  static const int wettingPhaseIdx = 0;
81  static const int nonWettingPhaseIdx = 1;
82 
84  static const char* phaseName(unsigned phaseIdx)
85  {
86  assert(phaseIdx < numPhases);
87 
88  static const char* name[] = {
89  "wetting",
90  "nonwetting"
91  };
92  return name[phaseIdx];
93  }
94 
96  static bool isLiquid(unsigned phaseIdx)
97  {
98  //assert(0 <= phaseIdx && phaseIdx < numPhases);
99  return
100  (phaseIdx == wettingPhaseIdx)
101  ? WettingPhase::isLiquid()
102  : NonwettingPhase::isLiquid();
103  }
104 
106  static bool isCompressible(unsigned phaseIdx)
107  {
108  //assert(0 <= phaseIdx && phaseIdx < numPhases);
109 
110  return
111  (phaseIdx == wettingPhaseIdx)
112  ? WettingPhase::isCompressible()
113  : NonwettingPhase::isCompressible();
114  }
115 
117  static bool isIdealGas(unsigned phaseIdx)
118  {
119  //assert(0 <= phaseIdx && phaseIdx < numPhases);
120 
121  // let the fluids decide
122  return
123  (phaseIdx == wettingPhaseIdx)
124  ? WettingPhase::isIdealGas()
125  : NonwettingPhase::isIdealGas();
126  }
127 
129  static bool isIdealMixture(unsigned /*phaseIdx*/)
130  {
131  //assert(0 <= phaseIdx && phaseIdx < numPhases);
132 
133  // we assume immisibility
134  return true;
135  }
136 
137  /****************************************
138  * Component related static parameters
139  ****************************************/
140 
142  static const int numComponents = 2;
143 
145  static const int wettingCompIdx = 0;
147  static const int nonWettingCompIdx = 1;
148 
150  static const char* componentName(unsigned compIdx)
151  {
152  assert(compIdx < numComponents);
153 
154  if (compIdx == wettingCompIdx)
155  return WettingPhase::name();
156  return NonwettingPhase::name();
157  }
158 
160  static Scalar molarMass(unsigned compIdx)
161  {
162  //assert(0 <= compIdx && compIdx < numComponents);
163 
164  // let the fluids decide
165  return
166  (compIdx == wettingCompIdx)
167  ? WettingPhase::molarMass()
168  : NonwettingPhase::molarMass();
169  }
170 
174  static Scalar criticalTemperature(unsigned compIdx)
175  {
176  //assert(0 <= compIdx && compIdx < numComponents);
177  // let the fluids decide
178  return
179  (compIdx == wettingCompIdx)
180  ? WettingPhase::criticalTemperature()
181  : NonwettingPhase::criticalTemperature();
182  }
183 
187  static Scalar criticalPressure(unsigned compIdx)
188  {
189  //assert(0 <= compIdx && compIdx < numComponents);
190  // let the fluids decide
191  return
192  (compIdx == wettingCompIdx)
193  ? WettingPhase::criticalPressure()
194  : NonwettingPhase::criticalPressure();
195  }
196 
200  static Scalar acentricFactor(unsigned compIdx)
201  {
202  //assert(0 <= compIdx && compIdx < numComponents);
203  // let the fluids decide
204  return
205  (compIdx == wettingCompIdx)
206  ? WettingPhase::acentricFactor()
207  : NonwettingPhase::acentricFactor();
208  }
209 
210  /****************************************
211  * thermodynamic relations
212  ****************************************/
213 
215  static void init()
216  {
217  // two gaseous phases at once do not make sense physically!
218  // (But two liquids are fine)
219  assert(WettingPhase::isLiquid() || NonwettingPhase::isLiquid());
220  }
221 
223  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
224  static LhsEval density(const FluidState& fluidState,
225  const ParameterCache<ParamCacheEval>& /*paramCache*/,
226  unsigned phaseIdx)
227  {
228  assert(phaseIdx < numPhases);
229 
230  const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
231  const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
232  if (phaseIdx == wettingPhaseIdx)
233  return WettingPhase::density(temperature, pressure);
234  return NonwettingPhase::density(temperature, pressure);
235  }
236 
238  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
239  static LhsEval viscosity(const FluidState& fluidState,
240  const ParameterCache<ParamCacheEval>& /*paramCache*/,
241  unsigned phaseIdx)
242  {
243  assert(phaseIdx < numPhases);
244 
245  const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
246  const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
247  if (phaseIdx == wettingPhaseIdx)
248  return WettingPhase::viscosity(temperature, pressure);
249  return NonwettingPhase::viscosity(temperature, pressure);
250  }
251 
253  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
254  static LhsEval fugacityCoefficient(const FluidState& /*fluidState*/,
255  const ParameterCache<ParamCacheEval>& /*paramCache*/,
256  unsigned phaseIdx,
257  unsigned compIdx)
258  {
259  assert(phaseIdx < numPhases);
260  assert(compIdx < numComponents);
261 
262  if (phaseIdx == compIdx)
263  // TODO (?): calculate the real fugacity coefficient of
264  // the component in the fluid. Probably that's not worth
265  // the effort, since the fugacity coefficient of the other
266  // component is infinite anyway...
267  return 1.0;
268  return std::numeric_limits<Scalar>::infinity();
269  }
270 
272  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
273  static LhsEval enthalpy(const FluidState& fluidState,
274  const ParameterCache<ParamCacheEval>& /*paramCache*/,
275  unsigned phaseIdx)
276  {
277  assert(phaseIdx < numPhases);
278 
279  const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
280  const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
281  if (phaseIdx == wettingPhaseIdx)
282  return WettingPhase::enthalpy(temperature, pressure);
283  return NonwettingPhase::enthalpy(temperature, pressure);
284  }
285 
287  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
288  static LhsEval thermalConductivity(const FluidState& fluidState,
289  const ParameterCache<ParamCacheEval>& /*paramCache*/,
290  unsigned phaseIdx)
291  {
292  assert(phaseIdx < numPhases);
293 
294  const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
295  const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
296  if (phaseIdx == wettingPhaseIdx)
297  return WettingPhase::thermalConductivity(temperature, pressure);
298  return NonwettingPhase::thermalConductivity(temperature, pressure);
299  }
300 
302  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
303  static LhsEval heatCapacity(const FluidState& fluidState,
304  const ParameterCache<ParamCacheEval>& /*paramCache*/,
305  unsigned phaseIdx)
306  {
307  assert(phaseIdx < numPhases);
308 
309  const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
310  const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
311  if (phaseIdx == wettingPhaseIdx)
312  return WettingPhase::heatCapacity(temperature, pressure);
313  return NonwettingPhase::heatCapacity(temperature, pressure);
314  }
315 };
316 
317 } // namespace Opm
318 
319 #endif
The base class for all fluid systems.
Represents the gas phase of a single (pseudo-) component.
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system assum...
Represents the liquid phase of a single (pseudo-) component.
A parameter cache which does nothing.
The base class for all fluid systems.
Definition: BaseFluidSystem.hpp:44
Scalar Scalar
The type used for scalar quantities.
Definition: BaseFluidSystem.hpp:49
A parameter cache which does nothing.
Definition: NullParameterCache.hpp:40
A fluid system for two-phase models assuming immiscibility and thermodynamic equilibrium.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:58
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:129
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition: TwoPhaseImmiscibleFluidSystem.hpp:200
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition: TwoPhaseImmiscibleFluidSystem.hpp:288
static void init()
Initialize the fluid system's static parameters.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:215
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:96
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:150
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:84
static LhsEval fugacityCoefficient(const FluidState &, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:254
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition: TwoPhaseImmiscibleFluidSystem.hpp:160
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:117
static const int numPhases
Number of fluid phases in the fluid system.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:76
static const int wettingCompIdx
Index of the wetting phase's component.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:145
static const int nonWettingCompIdx
Index of the non-wetting phase's component.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:147
static LhsEval enthalpy(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy [J/kg...
Definition: TwoPhaseImmiscibleFluidSystem.hpp:273
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: TwoPhaseImmiscibleFluidSystem.hpp:239
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition: TwoPhaseImmiscibleFluidSystem.hpp:174
static LhsEval heatCapacity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Specific isobaric heat capacity of a fluid phase [J/kg].
Definition: TwoPhaseImmiscibleFluidSystem.hpp:303
static const int nonWettingPhaseIdx
Index of the non-wetting phase.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:81
static const int numComponents
Number of chemical species in the fluid system.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:142
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:224
static const int wettingPhaseIdx
Index of the wetting phase.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:79
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition: TwoPhaseImmiscibleFluidSystem.hpp:187
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: TwoPhaseImmiscibleFluidSystem.hpp:106
Definition: TwoPhaseImmiscibleFluidSystem.hpp:69