My Project
Loading...
Searching...
No Matches
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"
32
40
41#include <cassert>
42
43namespace Opm {
44
51template <class Scalar>
53 : public BaseFluidSystem<Scalar, H2ON2LiquidPhaseFluidSystem<Scalar> >
54{
57
58 // convenience typedefs
59 typedef ::Opm::H2O<Scalar> IapwsH2O;
60 typedef ::Opm::TabulatedComponent<Scalar, IapwsH2O > TabulatedH2O;
61 typedef ::Opm::N2<Scalar> SimpleN2;
62
63public:
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 int liquidPhaseIdx = 0;
78
80 static const char* phaseName([[maybe_unused]] unsigned phaseIdx)
81 {
82 assert(phaseIdx == liquidPhaseIdx);
83
84 return "liquid";
85 }
86
88 static bool isLiquid(unsigned /*phaseIdx*/)
89 {
90 //assert(phaseIdx == liquidPhaseIdx);
91 return true; //only water phase present
92 }
93
95 static bool isCompressible(unsigned /*phaseIdx*/)
96 {
97 //assert(0 <= phaseIdx && phaseIdx < numPhases);
98 // the water component decides for the liquid phase...
100 }
101
103 static bool isIdealGas(unsigned /*phaseIdx*/)
104 {
105 //assert(0 <= phaseIdx && phaseIdx < numPhases);
106 return false; // not a gas (only liquid phase present)
107 }
108
110 static bool isIdealMixture(unsigned /*phaseIdx*/)
111 {
112 //assert(0 <= phaseIdx && phaseIdx < numPhases);
113 // we assume Henry's and Rault's laws for the water phase and
114 // and no interaction between gas molecules of different
115 // components, so all phases are ideal mixtures!
116 return true;
117 }
118
119 /****************************************
120 * Component related static parameters
121 ****************************************/
122
124 static const int numComponents = 2;
125
127 static const int H2OIdx = 0;
129 static const int N2Idx = 1;
130
133 //typedef SimpleH2O H2O;
134 //typedef IapwsH2O H2O;
135
137 typedef SimpleN2 N2;
138
140 static const char* componentName(unsigned compIdx)
141 {
142 static const char* name[] = {
143 H2O::name(),
144 N2::name()
145 };
146
147 assert(compIdx < numComponents);
148 return name[compIdx];
149 }
150
152 static Scalar molarMass(unsigned compIdx)
153 {
154 //assert(0 <= compIdx && compIdx < numComponents);
155 return (compIdx == H2OIdx)
157 : (compIdx == N2Idx)
158 ? N2::molarMass()
159 : 1e30;
160 }
161
167 static Scalar criticalTemperature(unsigned compIdx)
168 {
169 //assert(0 <= compIdx && compIdx < numComponents);
170 return (compIdx == H2OIdx)
172 : (compIdx == N2Idx)
174 : 1e30;
175 }
176
182 static Scalar criticalPressure(unsigned compIdx)
183 {
184 //assert(0 <= compIdx && compIdx < numComponents);
185 return (compIdx == H2OIdx)
187 : (compIdx == N2Idx)
189 : 1e30;
190 }
191
197 static Scalar acentricFactor(unsigned compIdx)
198 {
199 //assert(0 <= compIdx && compIdx < numComponents);
200 return (compIdx == H2OIdx)
202 : (compIdx == N2Idx)
204 : 1e30;
205 }
206
207 /****************************************
208 * thermodynamic relations
209 ****************************************/
210
217 static void init()
218 {
219 init(/*tempMin=*/273.15,
220 /*tempMax=*/623.15,
221 /*numTemp=*/50,
222 /*pMin=*/0.0,
223 /*pMax=*/20e6,
224 /*numP=*/50);
225 }
226
238 static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
239 Scalar pressMin, Scalar pressMax, unsigned nPress)
240 {
241 if (H2O::isTabulated) {
242 TabulatedH2O::init(tempMin, tempMax, nTemp,
243 pressMin, pressMax, nPress);
244 }
245 }
246
248 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
249 static LhsEval density(const FluidState& fluidState,
250 const ParameterCache<ParamCacheEval>& /*paramCache*/,
251 unsigned phaseIdx)
252 {
253 assert(phaseIdx < numPhases);
254
255 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
256 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
257
258 LhsEval sumMoleFrac = 0;
259 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
260 sumMoleFrac += decay<LhsEval>(fluidState.moleFraction(phaseIdx, compIdx));
261
262 assert(phaseIdx == liquidPhaseIdx);
263
264 // assume ideal mixture where each molecule occupies the same volume regardless
265 // of whether it is water or nitrogen.
266 const LhsEval& clH2O = H2O::liquidDensity(T, p)/H2O::molarMass();
267
268 const auto& xlH2O = decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, H2OIdx));
269 const auto& xlN2 = decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, N2Idx));
270
271 return clH2O*(H2O::molarMass()*xlH2O + N2::molarMass()*xlN2)/sumMoleFrac;
272 }
273
275 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
276 static LhsEval viscosity(const FluidState& fluidState,
277 const ParameterCache<ParamCacheEval>& /*paramCache*/,
278 unsigned phaseIdx)
279 {
280 assert(phaseIdx == liquidPhaseIdx);
281
282 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
283 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
284
285 // assume pure water for the liquid phase
286 return H2O::liquidViscosity(T, p);
287 }
288
290 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
291 static LhsEval fugacityCoefficient(const FluidState& fluidState,
292 const ParameterCache<ParamCacheEval>& /*paramCache*/,
293 unsigned phaseIdx,
294 unsigned compIdx)
295 {
296 assert(phaseIdx == liquidPhaseIdx);
297 assert(compIdx < numComponents);
298
299 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
300 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
301
302 if (compIdx == H2OIdx)
303 return H2O::vaporPressure(T)/p;
304 return BinaryCoeff::H2O_N2::henry(T)/p;
305 }
306
308 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
309 static LhsEval diffusionCoefficient(const FluidState& fluidState,
310 const ParameterCache<ParamCacheEval>& /*paramCache*/,
311 unsigned phaseIdx,
312 unsigned /*compIdx*/)
313
314 {
315 assert(phaseIdx == liquidPhaseIdx);
316
317 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
318 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
319
321 }
322
324 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
325 static LhsEval enthalpy(const FluidState& fluidState,
326 const ParameterCache<ParamCacheEval>& /*paramCache*/,
327 unsigned phaseIdx)
328 {
329 assert (phaseIdx == liquidPhaseIdx);
330
331 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
332 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
333 Valgrind::CheckDefined(T);
334 Valgrind::CheckDefined(p);
335
336 // TODO: way to deal with the solutes???
337 return H2O::liquidEnthalpy(T, p);
338 }
339
341 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
342 static LhsEval thermalConductivity(const FluidState& fluidState,
343 const ParameterCache<ParamCacheEval>& /*paramCache*/,
344 const unsigned phaseIdx)
345 {
346 assert(phaseIdx == liquidPhaseIdx);
347
348 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
349 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
351 }
352
354 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
355 static LhsEval heatCapacity(const FluidState& fluidState,
356 const ParameterCache<ParamCacheEval>& /*paramCache*/,
357 unsigned phaseIdx)
358 {
359 assert (phaseIdx == liquidPhaseIdx);
360
361 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
362 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
363
364 return H2O::liquidHeatCapacity(T, p);
365 }
366};
367
368} // namespace Opm
369
370#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:46
Scalar Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:51
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:54
static const int numComponents
Number of chemical species in the fluid system.
Definition H2ON2LiquidPhaseFluidSystem.hpp:124
static const int N2Idx
The index of the component for molecular nitrogen.
Definition H2ON2LiquidPhaseFluidSystem.hpp:129
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:325
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:309
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition H2ON2LiquidPhaseFluidSystem.hpp:167
static const int liquidPhaseIdx
Index of the liquid phase.
Definition H2ON2LiquidPhaseFluidSystem.hpp:77
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition H2ON2LiquidPhaseFluidSystem.hpp:276
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, const unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition H2ON2LiquidPhaseFluidSystem.hpp:342
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:291
static LhsEval heatCapacity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Specific isobaric heat capacity of a fluid phase [J/kg].
Definition H2ON2LiquidPhaseFluidSystem.hpp:355
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition H2ON2LiquidPhaseFluidSystem.hpp:197
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition H2ON2LiquidPhaseFluidSystem.hpp:182
TabulatedH2O H2O
The type of the component for pure water.
Definition H2ON2LiquidPhaseFluidSystem.hpp:132
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition H2ON2LiquidPhaseFluidSystem.hpp:95
static const int H2OIdx
The index of the water component.
Definition H2ON2LiquidPhaseFluidSystem.hpp:127
static const int numPhases
Number of fluid phases in the fluid system.
Definition H2ON2LiquidPhaseFluidSystem.hpp:74
static void init()
Initialize the fluid system's static parameters.
Definition H2ON2LiquidPhaseFluidSystem.hpp:217
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition H2ON2LiquidPhaseFluidSystem.hpp:103
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:238
SimpleN2 N2
The type of the component for pure molecular nitrogen.
Definition H2ON2LiquidPhaseFluidSystem.hpp:137
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition H2ON2LiquidPhaseFluidSystem.hpp:80
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition H2ON2LiquidPhaseFluidSystem.hpp:110
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition H2ON2LiquidPhaseFluidSystem.hpp:140
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition H2ON2LiquidPhaseFluidSystem.hpp:152
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition H2ON2LiquidPhaseFluidSystem.hpp:249
static bool isLiquid(unsigned)
Return whether a phase is liquid.
Definition H2ON2LiquidPhaseFluidSystem.hpp:88
Material properties of pure water .
Definition H2O.hpp:65
Properties of pure molecular nitrogen .
Definition N2.hpp:49
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 const char * name()
A human readable name for nitrogen.
Definition N2.hpp:56
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 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
static const char * name()
A human readable name for the component.
Definition TabulatedComponent.hpp:215
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
The type of the fluid system's parameter cache.
Definition H2ON2LiquidPhaseFluidSystem.hpp:67