All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
H2ON2LiquidPhaseFluidSystem.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_H2O_N2_LIQUIDPHASE_FLUID_SYSTEM_HPP
28 #define OPM_H2O_N2_LIQUIDPHASE_FLUID_SYSTEM_HPP
29 
30 #include "BaseFluidSystem.hpp"
31 #include "NullParameterCache.hpp"
32 
39 #include <opm/common/Valgrind.hpp>
40 
41 #include <opm/common/Exceptions.hpp>
42 #include <opm/common/ErrorMacros.hpp>
43 
44 #include <iostream>
45 #include <cassert>
46 
47 namespace Opm {
48 namespace FluidSystems {
49 
56 template <class Scalar, bool useComplexRelations = true>
58  : public BaseFluidSystem<Scalar, H2ON2LiquidPhase<Scalar, useComplexRelations> >
59 {
62 
63  // convenience typedefs
64  typedef Opm::H2O<Scalar> IapwsH2O;
66  typedef Opm::N2<Scalar> SimpleN2;
67 
68 public:
70  template <class Evaluation>
71  struct ParameterCache : public Opm::NullParameterCache<Evaluation>
72  {};
73 
74  /****************************************
75  * Fluid phase related static parameters
76  ****************************************/
77 
79  static const int numPhases = 1;
80 
82  static const int liquidPhaseIdx = 0;
83 
85  static const char* phaseName(unsigned phaseIdx OPM_OPTIM_UNUSED)
86  {
87  assert(phaseIdx == liquidPhaseIdx);
88 
89  return "liquid";
90  }
91 
93  static bool isLiquid(unsigned /*phaseIdx*/)
94  {
95  //assert(phaseIdx == liquidPhaseIdx);
96  return true; //only water phase present
97  }
98 
100  static bool isCompressible(unsigned /*phaseIdx*/)
101  {
102  //assert(0 <= phaseIdx && phaseIdx < numPhases);
103  // the water component decides for the liquid phase...
104  return H2O::liquidIsCompressible();
105  }
106 
108  static bool isIdealGas(unsigned /*phaseIdx*/)
109  {
110  //assert(0 <= phaseIdx && phaseIdx < numPhases);
111  return false; // not a gas (only liquid phase present)
112  }
113 
115  static bool isIdealMixture(unsigned /*phaseIdx*/)
116  {
117  //assert(0 <= phaseIdx && phaseIdx < numPhases);
118  // we assume Henry's and Rault's laws for the water phase and
119  // and no interaction between gas molecules of different
120  // components, so all phases are ideal mixtures!
121  return true;
122  }
123 
124  /****************************************
125  * Component related static parameters
126  ****************************************/
127 
129  static const int numComponents = 2;
130 
132  static const int H2OIdx = 0;
134  static const int N2Idx = 1;
135 
137  typedef TabulatedH2O H2O;
138  //typedef SimpleH2O H2O;
139  //typedef IapwsH2O H2O;
140 
142  typedef SimpleN2 N2;
143 
145  static const char* componentName(unsigned compIdx)
146  {
147  static const char* name[] = {
148  H2O::name(),
149  N2::name()
150  };
151 
152  assert(0 <= compIdx && compIdx < numComponents);
153  return name[compIdx];
154  }
155 
157  static Scalar molarMass(unsigned compIdx)
158  {
159  //assert(0 <= compIdx && compIdx < numComponents);
160  return (compIdx == H2OIdx)
161  ? H2O::molarMass()
162  : (compIdx == N2Idx)
163  ? N2::molarMass()
164  : 1e30;
165  }
166 
172  static Scalar criticalTemperature(unsigned compIdx)
173  {
174  //assert(0 <= compIdx && compIdx < numComponents);
175  return (compIdx == H2OIdx)
177  : (compIdx == N2Idx)
179  : 1e30;
180  }
181 
187  static Scalar criticalPressure(unsigned compIdx)
188  {
189  //assert(0 <= compIdx && compIdx < numComponents);
190  return (compIdx == H2OIdx)
192  : (compIdx == N2Idx)
194  : 1e30;
195  }
196 
202  static Scalar acentricFactor(unsigned compIdx)
203  {
204  //assert(0 <= compIdx && compIdx < numComponents);
205  return (compIdx == H2OIdx)
207  : (compIdx == N2Idx)
208  ? N2::acentricFactor()
209  : 1e30;
210  }
211 
212  /****************************************
213  * thermodynamic relations
214  ****************************************/
215 
222  static void init()
223  {
224  init(/*tempMin=*/273.15,
225  /*tempMax=*/623.15,
226  /*numTemp=*/50,
227  /*pMin=*/0.0,
228  /*pMax=*/20e6,
229  /*numP=*/50);
230  }
231 
243  static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
244  Scalar pressMin, Scalar pressMax, unsigned nPress)
245  {
246  if (H2O::isTabulated) {
247  TabulatedH2O::init(tempMin, tempMax, nTemp,
248  pressMin, pressMax, nPress);
249  }
250  }
251 
253  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
254  static LhsEval density(const FluidState& fluidState,
255  const ParameterCache<ParamCacheEval>& /*paramCache*/,
256  unsigned phaseIdx)
257  {
258  assert(0 <= phaseIdx && phaseIdx < numPhases);
259 
260  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
261  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
262 
263  LhsEval sumMoleFrac = 0;
264  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
265  sumMoleFrac += Opm::decay<LhsEval>(fluidState.moleFraction(phaseIdx, compIdx));
266 
267  assert(phaseIdx == liquidPhaseIdx);
268 
269  if (!useComplexRelations)
270  // assume pure water
271  return H2O::liquidDensity(T, p);
272  else
273  {
274  // See: Ochs 2008
275  const LhsEval& rholH2O = H2O::liquidDensity(T, p);
276  const LhsEval& clH2O = rholH2O/H2O::molarMass();
277 
278  const auto& xlH2O = Opm::decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, H2OIdx));
279  const auto& xlN2 = Opm::decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, N2Idx));
280 
281  // this assumes each nitrogen molecule displaces exactly one
282  // water molecule in the liquid
283  return clH2O*(H2O::molarMass()*xlH2O + N2::molarMass()*xlN2)/sumMoleFrac;
284  }
285  }
286 
288  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
289  static LhsEval viscosity(const FluidState& fluidState,
290  const ParameterCache<ParamCacheEval>& /*paramCache*/,
291  unsigned phaseIdx)
292  {
293  assert(phaseIdx == liquidPhaseIdx);
294 
295  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
296  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
297 
298  // assume pure water for the liquid phase
299  return H2O::liquidViscosity(T, p);
300  }
301 
303  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
304  static LhsEval fugacityCoefficient(const FluidState& fluidState,
305  const ParameterCache<ParamCacheEval>& /*paramCache*/,
306  unsigned phaseIdx,
307  unsigned compIdx)
308  {
309  assert(phaseIdx == liquidPhaseIdx);
310  assert(0 <= compIdx && compIdx < numComponents);
311 
312  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
313  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
314 
315  if (compIdx == H2OIdx)
316  return H2O::vaporPressure(T)/p;
318  }
319 
321  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
322  static LhsEval diffusionCoefficient(const FluidState& fluidState,
323  const ParameterCache<ParamCacheEval>& /*paramCache*/,
324  unsigned phaseIdx,
325  unsigned /*compIdx*/)
326 
327  {
328  assert(phaseIdx == liquidPhaseIdx);
329 
330  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
331  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
332 
334  }
335 
337  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
338  static LhsEval enthalpy(const FluidState& fluidState,
339  const ParameterCache<ParamCacheEval>& /*paramCache*/,
340  unsigned phaseIdx)
341  {
342  assert (phaseIdx == liquidPhaseIdx);
343 
344  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
345  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
346  Valgrind::CheckDefined(T);
347  Valgrind::CheckDefined(p);
348 
349  // TODO: way to deal with the solutes???
350  return H2O::liquidEnthalpy(T, p);
351  }
352 
354  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
355  static LhsEval thermalConductivity(const FluidState& fluidState,
356  const ParameterCache<ParamCacheEval>& /*paramCache*/,
357  const unsigned phaseIdx)
358  {
359  assert(phaseIdx == liquidPhaseIdx);
360 
361  if(useComplexRelations){
362  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
363  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
364  return H2O::liquidThermalConductivity(T, p);
365  }
366  else
367  return 0.578078; // conductivity of water[W / (m K ) ] IAPWS evaluated at p=.1 MPa, T=8C
368  }
369 
371  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
372  static LhsEval heatCapacity(const FluidState& fluidState,
373  const ParameterCache<ParamCacheEval>& /*paramCache*/,
374  unsigned phaseIdx)
375  {
376  assert (phaseIdx == liquidPhaseIdx);
377 
378  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
379  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
380 
381  return H2O::liquidHeatCapacity(T, p);
382  }
383 };
384 
385 } // namespace FluidSystems
386 
387 } // namespace Opm
388 
389 #endif
Material properties of pure water .
Definition: H2O.hpp:61
static const char * name()
A human readable name for nitrogen.
Definition: N2.hpp:56
static Evaluation henry(const Evaluation &temperature)
Henry coefficent for molecular nitrogen in liquid water.
Definition: H2O_N2.hpp:52
static Evaluation liquidDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of liquid at a given pressure and temperature .
Definition: TabulatedComponent.hpp:435
static Evaluation liquidEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the liquid .
Definition: TabulatedComponent.hpp:290
static LhsEval heatCapacity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Specific isobaric heat capacity of a fluid phase [J/kg].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:372
static Evaluation liquidHeatCapacity(const Evaluation &temperature, const Evaluation &pressure)
Specific isobaric heat capacity of the liquid .
Definition: TabulatedComponent.hpp:324
static Scalar molarMass()
The molar mass in of the component.
Definition: TabulatedComponent.hpp:224
A simple version of pure water.
Relations valid for an ideal gas.
TabulatedH2O H2O
The type of the component for pure water.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:137
static void init()
Initialize the fluid system&#39;s static parameters.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:222
Properties of pure molecular nitrogen .
Definition: N2.hpp:48
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:254
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:145
static const Scalar acentricFactor()
The acentric factor of water.
Definition: H2O.hpp:86
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:172
static Evaluation liquidDiffCoeff(const Evaluation &temperature, const Evaluation &)
Diffusion coefficent for molecular nitrogen in liquid water.
Definition: H2O_N2.hpp:102
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:202
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:115
static Scalar criticalPressure()
Returns the critical pressure of molecular nitrogen.
Definition: N2.hpp:74
static LhsEval fugacityCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:304
static const int N2Idx
The index of the component for molecular nitrogen.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:134
static const int liquidPhaseIdx
Index of the liquid phase.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:82
A parameter cache which does nothing.
Material properties of pure water .
static LhsEval diffusionCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx, unsigned)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition: H2ON2LiquidPhaseFluidSystem.hpp:322
The type of the fluid system&#39;s parameter cache.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:71
static Evaluation liquidThermalConductivity(const Evaluation &temperature, const Evaluation &pressure)
The thermal conductivity of liquid water .
Definition: TabulatedComponent.hpp:503
static const char * name()
A human readable name for the component.
Definition: TabulatedComponent.hpp:218
static Evaluation vaporPressure(const Evaluation &temperature)
The vapor pressure in of the component at a given temperature.
Definition: TabulatedComponent.hpp:258
static Scalar molarMass()
The molar mass in of molecular nitrogen.
Definition: N2.hpp:62
static bool isLiquid(unsigned)
Return whether a phase is liquid.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:93
static Scalar criticalPressure()
Returns the critical pressure in of the component.
Definition: TabulatedComponent.hpp:236
static const char * phaseName(unsigned phaseIdx OPM_OPTIM_UNUSED)
Return the human readable name of a fluid phase.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:85
static Scalar criticalTemperature()
Returns the critical temperature in of the component.
Definition: TabulatedComponent.hpp:230
A liquid-phase-only fluid system with water and nitrogen as components.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:57
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, const unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:355
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the fluid system&#39;s static parameters using problem specific temperature and pressure range...
Definition: H2ON2LiquidPhaseFluidSystem.hpp:243
static const int numComponents
Number of chemical species in the fluid system.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:129
Properties of pure molecular nitrogen .
A generic class which tabulates all thermodynamic properties of a given component.
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: H2ON2LiquidPhaseFluidSystem.hpp:338
The base class for all fluid systems.
Definition: BaseFluidSystem.hpp:43
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:289
SimpleN2 N2
The type of the component for pure molecular nitrogen.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:142
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:100
A generic class which tabulates all thermodynamic properties of a given component.
Definition: TabulatedComponent.hpp:58
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:108
A parameter cache which does nothing.
Definition: NullParameterCache.hpp:39
static const int H2OIdx
The index of the water component.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:132
Scalar Scalar
The type used for scalar quantities.
Definition: BaseFluidSystem.hpp:49
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the tables.
Definition: TabulatedComponent.hpp:75
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:187
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:157
static Scalar criticalTemperature()
Returns the critical temperature of molecular nitrogen.
Definition: N2.hpp:68
static const int numPhases
Number of fluid phases in the fluid system.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:79
static Evaluation liquidViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of liquid.
Definition: TabulatedComponent.hpp:469
The base class for all fluid systems.
static bool liquidIsCompressible()
Returns true iff the liquid phase is assumed to be compressible.
Definition: TabulatedComponent.hpp:399
Binary coefficients for water and nitrogen.