My Project
SinglePhaseFluidSystem.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_SINGLE_PHASE_FLUIDSYSTEM_HPP
28 #define OPM_SINGLE_PHASE_FLUIDSYSTEM_HPP
29 
30 #include "BaseFluidSystem.hpp"
31 #include "NullParameterCache.hpp"
32 
39 
41 
42 #include <limits>
43 #include <cassert>
44 
45 namespace Opm {
46 
56 template <class Scalar, class Fluid>
58  : public BaseFluidSystem<Scalar, SinglePhaseFluidSystem<Scalar, Fluid> >
59 {
62 
63 public:
65  template <class Evaluation>
66  struct ParameterCache : public NullParameterCache<Evaluation>
67  {};
68 
69  /****************************************
70  * Fluid phase related static parameters
71  ****************************************/
72 
74  static const int numPhases = 1;
75 
77  static const char* phaseName(unsigned phaseIdx OPM_OPTIM_UNUSED)
78  {
79  assert(phaseIdx < numPhases);
80 
81  return Fluid::name();
82  }
83 
85  static bool isLiquid(unsigned /*phaseIdx*/)
86  {
87  //assert(0 <= phaseIdx && phaseIdx < numPhases);
88 
89  return Fluid::isLiquid();
90  }
91 
93  static bool isCompressible(unsigned /*phaseIdx*/)
94  {
95  //assert(0 <= phaseIdx && phaseIdx < numPhases);
96 
97  // let the fluid decide
98  return Fluid::isCompressible();
99  }
100 
102  static bool isIdealMixture(unsigned /*phaseIdx*/)
103  {
104  //assert(0 <= phaseIdx && phaseIdx < numPhases);
105 
106  // we assume immisibility
107  return true;
108  }
109 
111  static bool isIdealGas(unsigned /*phaseIdx*/)
112  {
113  //assert(0 <= phaseIdx && phaseIdx < numPhases);
114 
115  // let the fluid decide
116  return Fluid::isIdealGas();
117  }
118 
119  /****************************************
120  * Component related static parameters
121  ****************************************/
122 
124  static const int numComponents = 1;
125 
127  static const char* componentName(unsigned compIdx OPM_OPTIM_UNUSED)
128  {
129  assert(compIdx < numComponents);
130 
131  return Fluid::name();
132  }
133 
135  static Scalar molarMass(unsigned /*compIdx*/)
136  {
137  //assert(0 <= compIdx && compIdx < numComponents);
138 
139  return Fluid::molarMass();
140  }
141 
147  static Scalar criticalTemperature(unsigned /*compIdx*/)
148  {
149  //assert(0 <= compIdx && compIdx < numComponents);
150 
151  return Fluid::criticalTemperature();
152  }
153 
159  static Scalar criticalPressure(unsigned /*compIdx*/)
160  {
161  //assert(0 <= compIdx && compIdx < numComponents);
162 
163  return Fluid::criticalPressure();
164  }
165 
171  static Scalar acentricFactor(unsigned /*compIdx*/)
172  {
173  //assert(0 <= compIdx && compIdx < numComponents);
174 
175  return Fluid::acentricFactor();
176  }
177 
178  /****************************************
179  * thermodynamic relations
180  ****************************************/
181 
183  static void init()
184  { }
185 
187  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
188  static LhsEval density(const FluidState& fluidState,
189  const ParameterCache<ParamCacheEval>& /*paramCache*/,
190  unsigned phaseIdx)
191  {
192  assert(phaseIdx < numPhases);
193 
194  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
195  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
196  return Fluid::density(T, p);
197  }
198 
200  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
201  static LhsEval viscosity(const FluidState& fluidState,
202  const ParameterCache<ParamCacheEval>& /*paramCache*/,
203  unsigned phaseIdx)
204  {
205  assert(phaseIdx < numPhases);
206 
207  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
208  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
209  return Fluid::viscosity(T, p);
210  }
211 
213  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
214  static LhsEval fugacityCoefficient(const FluidState& /*fluidState*/,
215  const ParameterCache<ParamCacheEval>& /*paramCache*/,
216  unsigned phaseIdx,
217  unsigned compIdx)
218  {
219  assert(phaseIdx < numPhases);
220  assert(compIdx < numComponents);
221 
222  if (phaseIdx == compIdx)
223  // TODO (?): calculate the real fugacity coefficient of
224  // the component in the fluid. Probably that's not worth
225  // the effort, since the fugacity coefficient of the other
226  // component is infinite anyway...
227  return 1.0;
228  return std::numeric_limits<Scalar>::infinity();
229  }
230 
232  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
233  static LhsEval enthalpy(const FluidState& fluidState,
234  const ParameterCache<ParamCacheEval>& /*paramCache*/,
235  unsigned phaseIdx)
236  {
237  assert(phaseIdx < numPhases);
238 
239  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
240  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
241  return Fluid::enthalpy(T, p);
242  }
243 
245  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
246  static LhsEval thermalConductivity(const FluidState& fluidState,
247  const ParameterCache<ParamCacheEval>& /*paramCache*/,
248  unsigned phaseIdx)
249  {
250  assert(phaseIdx < numPhases);
251 
252  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
253  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
254  return Fluid::thermalConductivity(T, p);
255  }
256 
258  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
259  static LhsEval heatCapacity(const FluidState& fluidState,
260  const ParameterCache<ParamCacheEval>& /*paramCache*/,
261  unsigned phaseIdx)
262  {
263  assert(phaseIdx < numPhases);
264 
265  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
266  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
267  return Fluid::heatCapacity(T, p);
268  }
269 };
270 
271 } // namespace Opm
272 
273 #endif
The base class for all fluid systems.
Represents the gas phase of a single (pseudo-) component.
Material properties of pure water .
Represents the liquid phase of a single (pseudo-) component.
Properties of pure molecular nitrogen .
A parameter cache which does nothing.
A simple version of pure water.
A generic class which tabulates all thermodynamic properties of a given component.
Provides the OPM_UNUSED macro.
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 single phase models.
Definition: SinglePhaseFluidSystem.hpp:59
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition: SinglePhaseFluidSystem.hpp:246
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: SinglePhaseFluidSystem.hpp:201
static const int numComponents
Number of chemical species in the fluid system.
Definition: SinglePhaseFluidSystem.hpp:124
static Scalar criticalTemperature(unsigned)
Critical temperature of a component [K].
Definition: SinglePhaseFluidSystem.hpp:147
static Scalar acentricFactor(unsigned)
The acentric factor of a component [].
Definition: SinglePhaseFluidSystem.hpp:171
static const int numPhases
Number of fluid phases in the fluid system.
Definition: SinglePhaseFluidSystem.hpp:74
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: SinglePhaseFluidSystem.hpp:214
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: SinglePhaseFluidSystem.hpp:111
static void init()
Initialize the fluid system's static parameters.
Definition: SinglePhaseFluidSystem.hpp:183
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: SinglePhaseFluidSystem.hpp:188
static Scalar molarMass(unsigned)
Return the molar mass of a component in [kg/mol].
Definition: SinglePhaseFluidSystem.hpp:135
static const char * componentName(unsigned compIdx OPM_OPTIM_UNUSED)
Return the human readable name of a component.
Definition: SinglePhaseFluidSystem.hpp:127
static Scalar criticalPressure(unsigned)
Critical pressure of a component [Pa].
Definition: SinglePhaseFluidSystem.hpp:159
static LhsEval heatCapacity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Specific isobaric heat capacity of a fluid phase [J/kg].
Definition: SinglePhaseFluidSystem.hpp:259
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: SinglePhaseFluidSystem.hpp:93
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: SinglePhaseFluidSystem.hpp:102
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: SinglePhaseFluidSystem.hpp:233
static bool isLiquid(unsigned)
Return whether a phase is liquid.
Definition: SinglePhaseFluidSystem.hpp:85
static const char * phaseName(unsigned phaseIdx OPM_OPTIM_UNUSED)
Return the human readable name of a fluid phase.
Definition: SinglePhaseFluidSystem.hpp:77
The type of the fluid system's parameter cache.
Definition: SinglePhaseFluidSystem.hpp:67