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