My Project
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 
40 
41 #include <iostream>
42 #include <cassert>
43 
44 namespace Opm {
45 
52 template <class Scalar>
54  : public BaseFluidSystem<Scalar, H2ON2LiquidPhaseFluidSystem<Scalar> >
55 {
58 
59  // convenience typedefs
60  typedef ::Opm::H2O<Scalar> IapwsH2O;
61  typedef ::Opm::TabulatedComponent<Scalar, IapwsH2O > TabulatedH2O;
62  typedef ::Opm::N2<Scalar> SimpleN2;
63 
64 public:
66  template <class Evaluation>
67  struct ParameterCache : public NullParameterCache<Evaluation>
68  {};
69 
70  /****************************************
71  * Fluid phase related static parameters
72  ****************************************/
73 
75  static const int numPhases = 1;
76 
78  static const int liquidPhaseIdx = 0;
79 
81  static const char* phaseName([[maybe_unused]] unsigned phaseIdx)
82  {
83  assert(phaseIdx == liquidPhaseIdx);
84 
85  return "liquid";
86  }
87 
89  static bool isLiquid(unsigned /*phaseIdx*/)
90  {
91  //assert(phaseIdx == liquidPhaseIdx);
92  return true; //only water phase present
93  }
94 
96  static bool isCompressible(unsigned /*phaseIdx*/)
97  {
98  //assert(0 <= phaseIdx && phaseIdx < numPhases);
99  // the water component decides for the liquid phase...
100  return H2O::liquidIsCompressible();
101  }
102 
104  static bool isIdealGas(unsigned /*phaseIdx*/)
105  {
106  //assert(0 <= phaseIdx && phaseIdx < numPhases);
107  return false; // not a gas (only liquid phase present)
108  }
109 
111  static bool isIdealMixture(unsigned /*phaseIdx*/)
112  {
113  //assert(0 <= phaseIdx && phaseIdx < numPhases);
114  // we assume Henry's and Rault's laws for the water phase and
115  // and no interaction between gas molecules of different
116  // components, so all phases are ideal mixtures!
117  return true;
118  }
119 
120  /****************************************
121  * Component related static parameters
122  ****************************************/
123 
125  static const int numComponents = 2;
126 
128  static const int H2OIdx = 0;
130  static const int N2Idx = 1;
131 
133  typedef TabulatedH2O H2O;
134  //typedef SimpleH2O H2O;
135  //typedef IapwsH2O H2O;
136 
138  typedef SimpleN2 N2;
139 
141  static const char* componentName(unsigned compIdx)
142  {
143  static const char* name[] = {
144  H2O::name(),
145  N2::name()
146  };
147 
148  assert(compIdx < numComponents);
149  return name[compIdx];
150  }
151 
153  static Scalar molarMass(unsigned compIdx)
154  {
155  //assert(0 <= compIdx && compIdx < numComponents);
156  return (compIdx == H2OIdx)
157  ? H2O::molarMass()
158  : (compIdx == N2Idx)
159  ? N2::molarMass()
160  : 1e30;
161  }
162 
168  static Scalar criticalTemperature(unsigned compIdx)
169  {
170  //assert(0 <= compIdx && compIdx < numComponents);
171  return (compIdx == H2OIdx)
173  : (compIdx == N2Idx)
175  : 1e30;
176  }
177 
183  static Scalar criticalPressure(unsigned compIdx)
184  {
185  //assert(0 <= compIdx && compIdx < numComponents);
186  return (compIdx == H2OIdx)
188  : (compIdx == N2Idx)
190  : 1e30;
191  }
192 
198  static Scalar acentricFactor(unsigned compIdx)
199  {
200  //assert(0 <= compIdx && compIdx < numComponents);
201  return (compIdx == H2OIdx)
203  : (compIdx == N2Idx)
205  : 1e30;
206  }
207 
208  /****************************************
209  * thermodynamic relations
210  ****************************************/
211 
218  static void init()
219  {
220  init(/*tempMin=*/273.15,
221  /*tempMax=*/623.15,
222  /*numTemp=*/50,
223  /*pMin=*/0.0,
224  /*pMax=*/20e6,
225  /*numP=*/50);
226  }
227 
239  static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
240  Scalar pressMin, Scalar pressMax, unsigned nPress)
241  {
242  if (H2O::isTabulated) {
243  TabulatedH2O::init(tempMin, tempMax, nTemp,
244  pressMin, pressMax, nPress);
245  }
246  }
247 
249  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
250  static LhsEval density(const FluidState& fluidState,
251  const ParameterCache<ParamCacheEval>& /*paramCache*/,
252  unsigned phaseIdx)
253  {
254  assert(phaseIdx < numPhases);
255 
256  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
257  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
258 
259  LhsEval sumMoleFrac = 0;
260  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
261  sumMoleFrac += decay<LhsEval>(fluidState.moleFraction(phaseIdx, compIdx));
262 
263  assert(phaseIdx == liquidPhaseIdx);
264 
265  // assume ideal mixture where each molecule occupies the same volume regardless
266  // of whether it is water or nitrogen.
267  const LhsEval& clH2O = H2O::liquidDensity(T, p)/H2O::molarMass();
268 
269  const auto& xlH2O = decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, H2OIdx));
270  const auto& xlN2 = decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, N2Idx));
271 
272  return clH2O*(H2O::molarMass()*xlH2O + N2::molarMass()*xlN2)/sumMoleFrac;
273  }
274 
276  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
277  static LhsEval viscosity(const FluidState& fluidState,
278  const ParameterCache<ParamCacheEval>& /*paramCache*/,
279  unsigned phaseIdx)
280  {
281  assert(phaseIdx == liquidPhaseIdx);
282 
283  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
284  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
285 
286  // assume pure water for the liquid phase
287  return H2O::liquidViscosity(T, p);
288  }
289 
291  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
292  static LhsEval fugacityCoefficient(const FluidState& fluidState,
293  const ParameterCache<ParamCacheEval>& /*paramCache*/,
294  unsigned phaseIdx,
295  unsigned compIdx)
296  {
297  assert(phaseIdx == liquidPhaseIdx);
298  assert(compIdx < numComponents);
299 
300  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
301  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
302 
303  if (compIdx == H2OIdx)
304  return H2O::vaporPressure(T)/p;
305  return BinaryCoeff::H2O_N2::henry(T)/p;
306  }
307 
309  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
310  static LhsEval diffusionCoefficient(const FluidState& fluidState,
311  const ParameterCache<ParamCacheEval>& /*paramCache*/,
312  unsigned phaseIdx,
313  unsigned /*compIdx*/)
314 
315  {
316  assert(phaseIdx == liquidPhaseIdx);
317 
318  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
319  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
320 
322  }
323 
325  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
326  static LhsEval enthalpy(const FluidState& fluidState,
327  const ParameterCache<ParamCacheEval>& /*paramCache*/,
328  unsigned phaseIdx)
329  {
330  assert (phaseIdx == liquidPhaseIdx);
331 
332  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
333  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
334  Valgrind::CheckDefined(T);
335  Valgrind::CheckDefined(p);
336 
337  // TODO: way to deal with the solutes???
338  return H2O::liquidEnthalpy(T, p);
339  }
340 
342  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
343  static LhsEval thermalConductivity(const FluidState& fluidState,
344  const ParameterCache<ParamCacheEval>& /*paramCache*/,
345  const unsigned phaseIdx)
346  {
347  assert(phaseIdx == liquidPhaseIdx);
348 
349  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
350  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
351  return H2O::liquidThermalConductivity(T, p);
352  }
353 
355  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
356  static LhsEval heatCapacity(const FluidState& fluidState,
357  const ParameterCache<ParamCacheEval>& /*paramCache*/,
358  unsigned phaseIdx)
359  {
360  assert (phaseIdx == liquidPhaseIdx);
361 
362  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
363  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
364 
365  return H2O::liquidHeatCapacity(T, p);
366  }
367 };
368 
369 } // namespace Opm
370 
371 #endif
The base class for all fluid systems.
Material properties of pure water .
Binary coefficients for water and nitrogen.
Relations valid for an ideal gas.
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.
Some templates to wrap the valgrind client request macros.
The base class for all fluid systems.
Definition: BaseFluidSystem.hpp:44
Scalar Scalar
The type used for scalar quantities.
Definition: BaseFluidSystem.hpp:49
static Evaluation liquidDiffCoeff(const Evaluation &temperature, const Evaluation &)
Diffusion coefficent for molecular nitrogen in liquid water.
Definition: H2O_N2.hpp:102
static Evaluation henry(const Evaluation &temperature)
Henry coefficent for molecular nitrogen in liquid water.
Definition: H2O_N2.hpp:52
A liquid-phase-only fluid system with water and nitrogen as components.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:55
static const char * phaseName([[maybe_unused]] unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:81
static const int numComponents
Number of chemical species in the fluid system.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:125
static const int N2Idx
The index of the component for molecular nitrogen.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:130
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: H2ON2LiquidPhaseFluidSystem.hpp:326
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:310
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:168
static const int liquidPhaseIdx
Index of the liquid phase.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:78
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:277
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, const unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:343
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:292
static LhsEval heatCapacity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Specific isobaric heat capacity of a fluid phase [J/kg].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:356
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:198
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:183
TabulatedH2O H2O
The type of the component for pure water.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:133
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:96
static const int H2OIdx
The index of the water component.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:128
static const int numPhases
Number of fluid phases in the fluid system.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:75
static void init()
Initialize the fluid system's static parameters.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:218
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:141
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:104
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the fluid system's static parameters using problem specific temperature and pressure range...
Definition: H2ON2LiquidPhaseFluidSystem.hpp:239
SimpleN2 N2
The type of the component for pure molecular nitrogen.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:138
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:111
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition: H2ON2LiquidPhaseFluidSystem.hpp:153
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:250
static bool isLiquid(unsigned)
Return whether a phase is liquid.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:89
Material properties of pure water .
Definition: H2O.hpp:62
Properties of pure molecular nitrogen .
Definition: N2.hpp:49
static const char * name()
A human readable name for nitrogen.
Definition: N2.hpp:56
static Scalar criticalPressure()
Returns the critical pressure of molecular nitrogen.
Definition: N2.hpp:74
static Scalar criticalTemperature()
Returns the critical temperature of molecular nitrogen.
Definition: N2.hpp:68
static Scalar acentricFactor()
Acentric factor of .
Definition: N2.hpp:85
static Scalar molarMass()
The molar mass in of molecular nitrogen.
Definition: N2.hpp:62
A parameter cache which does nothing.
Definition: NullParameterCache.hpp:40
A generic class which tabulates all thermodynamic properties of a given component.
Definition: TabulatedComponent.hpp:56
static Evaluation liquidHeatCapacity(const Evaluation &temperature, const Evaluation &pressure)
Specific isobaric heat capacity of the liquid .
Definition: TabulatedComponent.hpp:333
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the tables.
Definition: TabulatedComponent.hpp:72
static Scalar criticalTemperature()
Returns the critical temperature in of the component.
Definition: TabulatedComponent.hpp:227
static Scalar criticalPressure()
Returns the critical pressure in of the component.
Definition: TabulatedComponent.hpp:233
static Scalar molarMass()
The molar mass in of the component.
Definition: TabulatedComponent.hpp:221
static bool liquidIsCompressible()
Returns true iff the liquid phase is assumed to be compressible.
Definition: TabulatedComponent.hpp:408
static Evaluation liquidViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of liquid.
Definition: TabulatedComponent.hpp:478
static const char * name()
A human readable name for the component.
Definition: TabulatedComponent.hpp:215
static Scalar acentricFactor()
Returns the acentric factor of the component.
Definition: TabulatedComponent.hpp:239
static Evaluation liquidDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of liquid at a given pressure and temperature .
Definition: TabulatedComponent.hpp:444
static Evaluation vaporPressure(const Evaluation &temperature)
The vapor pressure in of the component at a given temperature.
Definition: TabulatedComponent.hpp:267
static Evaluation liquidThermalConductivity(const Evaluation &temperature, const Evaluation &pressure)
The thermal conductivity of liquid water .
Definition: TabulatedComponent.hpp:512
static Evaluation liquidEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the liquid .
Definition: TabulatedComponent.hpp:299
The type of the fluid system's parameter cache.
Definition: H2ON2LiquidPhaseFluidSystem.hpp:68