My Project
H2OAirMesityleneFluidSystem.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_AIR_MESITYLENE_FLUID_SYSTEM_HPP
28#define OPM_H2O_AIR_MESITYLENE_FLUID_SYSTEM_HPP
29
30#include "BaseFluidSystem.hpp"
32
44
45namespace Opm {
46
52template <class Scalar>
54 : public BaseFluidSystem<Scalar, H2OAirMesityleneFluidSystem<Scalar> >
55{
58
59 typedef ::Opm::H2O<Scalar> IapwsH2O;
60 typedef TabulatedComponent<Scalar, IapwsH2O, /*alongVaporPressure=*/false> TabulatedH2O;
61
62public:
63 template <class Evaluation>
64 struct ParameterCache : public NullParameterCache<Evaluation>
65 {};
66
69
71 typedef ::Opm::Air<Scalar> Air;
72
74 //typedef SimpleH2O H2O;
76 //typedef IapwsH2O H2O;
77
79 static const int numPhases = 3;
81 static const int numComponents = 3;
82
84 static const int waterPhaseIdx = 0;
86 static const int naplPhaseIdx = 1;
88 static const int gasPhaseIdx = 2;
89
91 static const int H2OIdx = 0;
93 static const int NAPLIdx = 1;
95 static const int airIdx = 2;
96
98 static void init()
99 {
100 init(/*tempMin=*/273.15,
101 /*tempMax=*/623.15,
102 /*numTemp=*/50,
103 /*pMin=*/0.0,
104 /*pMax=*/20e6,
105 /*numP=*/50);
106 }
107
119 static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
120 Scalar pressMin, Scalar pressMax, unsigned nPress)
121 {
122 if (H2O::isTabulated) {
123 TabulatedH2O::init(tempMin, tempMax, nTemp,
124 pressMin, pressMax, nPress);
125 }
126 }
127
129 static bool isLiquid(unsigned phaseIdx)
130 {
131 //assert(0 <= phaseIdx && phaseIdx < numPhases);
132 return phaseIdx != gasPhaseIdx;
133 }
134
136 static bool isIdealGas(unsigned phaseIdx)
137 { return phaseIdx == gasPhaseIdx && H2O::gasIsIdeal() && Air::gasIsIdeal() && NAPL::gasIsIdeal(); }
138
140 static bool isCompressible(unsigned phaseIdx)
141 {
142 //assert(0 <= phaseIdx && phaseIdx < numPhases);
143 // gases are always compressible
144 return (phaseIdx == gasPhaseIdx)
145 ? true
146 : (phaseIdx == waterPhaseIdx)
149 }
150
152 static bool isIdealMixture(unsigned /*phaseIdx*/)
153 {
154 //assert(0 <= phaseIdx && phaseIdx < numPhases);
155 // we assume Henry's and Rault's laws for the water phase and
156 // and no interaction between gas molecules of different
157 // components, so all phases are ideal mixtures!
158 return true;
159 }
160
162 static const char* phaseName(unsigned phaseIdx)
163 {
164 switch (phaseIdx) {
165 case waterPhaseIdx: return "water";
166 case naplPhaseIdx: return "napl";
167 case gasPhaseIdx: return "gas";
168 };
169 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
170 }
171
173 static const char* componentName(unsigned compIdx)
174 {
175 switch (compIdx) {
176 case H2OIdx: return H2O::name();
177 case airIdx: return Air::name();
178 case NAPLIdx: return NAPL::name();
179 };
180 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
181 }
182
184 static Scalar molarMass(unsigned compIdx)
185 {
186 return
187 (compIdx == H2OIdx)
189 : (compIdx == airIdx)
191 : (compIdx == NAPLIdx)
193 : -1e10;
194 //throw std::logic_error("Invalid component index "+std::to_string(compIdx));
195 }
196
198 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
199 static LhsEval density(const FluidState& fluidState,
200 const ParameterCache<ParamCacheEval>& /*paramCache*/,
201 unsigned phaseIdx)
202 {
203 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
204
205 if (phaseIdx == waterPhaseIdx) {
206 // See: Ochs 2008
207 const LhsEval& p =
209 ? decay<LhsEval>(fluidState.pressure(phaseIdx))
210 : 1e30;
211
212 const LhsEval& rholH2O = H2O::liquidDensity(T, p);
213 const LhsEval& clH2O = rholH2O/H2O::molarMass();
214
215 // this assumes each dissolved molecule displaces exactly one
216 // water molecule in the liquid
217 return
218 clH2O*(H2O::molarMass()*decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, H2OIdx)) +
219 Air::molarMass()*decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, airIdx)) +
220 NAPL::molarMass()*decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, NAPLIdx)));
221 }
222 else if (phaseIdx == naplPhaseIdx) {
223 // assume pure NAPL for the NAPL phase
224 const LhsEval& p =
226 ? decay<LhsEval>(fluidState.pressure(phaseIdx))
227 : 1e30;
228 return NAPL::liquidDensity(T, p);
229 }
230
231 assert (phaseIdx == gasPhaseIdx);
232 const LhsEval& pg = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
233 const LhsEval& pH2O = decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, H2OIdx))*pg;
234 const LhsEval& pAir = decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, airIdx))*pg;
235 const LhsEval& pNAPL = decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, NAPLIdx))*pg;
236 return
237 H2O::gasDensity(T, pH2O) +
238 Air::gasDensity(T, pAir) +
239 NAPL::gasDensity(T, pNAPL);
240 }
241
243 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
244 static LhsEval viscosity(const FluidState& fluidState,
245 const ParameterCache<ParamCacheEval>& /*paramCache*/,
246 unsigned phaseIdx)
247 {
248 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
249 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
250
251 if (phaseIdx == waterPhaseIdx) {
252 // assume pure water viscosity
253
254 return H2O::liquidViscosity(T,
255 p);
256 }
257 else if (phaseIdx == naplPhaseIdx) {
258 // assume pure NAPL viscosity
259 return NAPL::liquidViscosity(T, p);
260 }
261
262 assert (phaseIdx == gasPhaseIdx);
263
264 /* Wilke method. See:
265 *
266 * See: R. Reid, et al.: The Properties of Gases and Liquids,
267 * 4th edition, McGraw-Hill, 1987, 407-410
268 * 5th edition, McGraw-Hill, 20001, p. 9.21/22
269 *
270 * in this case, we use a simplified version in order to avoid
271 * computationally costly evaluation of sqrt and pow functions and
272 * divisions
273 * -- compare e.g. with Promo Class p. 32/33
274 */
275 const LhsEval mu[numComponents] = {
277 Air::gasViscosity(T, p),
279 };
280 // molar masses
281 const Scalar M[numComponents] = {
285 };
286
287 const LhsEval& xgAir = decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, airIdx));
288 const LhsEval& xgH2O = decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, H2OIdx));
289 const LhsEval& xgNapl = decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, NAPLIdx));
290 const LhsEval& xgAW = xgAir + xgH2O;
291 const LhsEval& muAW = (mu[airIdx]*xgAir + mu[H2OIdx]*xgH2O)/xgAW;
292 const LhsEval& MAW = (xgAir*Air::molarMass() + xgH2O*H2O::molarMass())/xgAW;
293
294 Scalar phiCAW = 0.3; // simplification for this particular system
295 /* actually like this
296 * Scalar phiCAW = std::pow(1.+std::sqrt(mu[NAPLIdx]/muAW)*std::pow(MAW/M[NAPLIdx],0.25),2)
297 * / std::sqrt(8.*(1.+M[NAPLIdx]/MAW));
298 */
299 const LhsEval& phiAWC = phiCAW * muAW*M[NAPLIdx]/(mu[NAPLIdx]*MAW);
300
301 return (xgAW*muAW)/(xgAW + xgNapl*phiAWC) + (xgNapl*mu[NAPLIdx])/(xgNapl + xgAW*phiCAW);
302 }
303
305 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
306 static LhsEval diffusionCoefficient(const FluidState& /*fluidState*/,
307 const ParameterCache<ParamCacheEval>& /*paramCache*/,
308 unsigned /*phaseIdx*/,
309 unsigned /*compIdx*/)
310 {
311 return 0;
312#if 0
314
315 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
316 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
317 LhsEval diffCont;
318
319 if (phaseIdx==gasPhaseIdx) {
320 const LhsEval& diffAC = BinaryCoeff::Air_Mesitylene::gasDiffCoeff(T, p);
321 const LhsEval& diffWC = BinaryCoeff::H2O_Mesitylene::gasDiffCoeff(T, p);
322 const LhsEval& diffAW = BinaryCoeff::H2O_Air::gasDiffCoeff(T, p);
323
324 const LhsEval& xga = decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, airIdx));
325 const LhsEval& xgw = decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, H2OIdx));
326 const LhsEval& xgc = decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, NAPLIdx));
327
328 if (compIdx==NAPLIdx) return (1 - xgw)/(xga/diffAW + xgc/diffWC);
329 else if (compIdx==H2OIdx) return (1 - xgc)/(xgw/diffWC + xga/diffAC);
330 else if (compIdx==airIdx) throw std::logic_error("Diffusivity of air in the gas phase "
331 "is constraint by sum of diffusive fluxes = 0 !\n");
332 }
333 else if (phaseIdx==waterPhaseIdx){
334 const LhsEval& diffACl = 1.e-9; // BinaryCoeff::Air_Mesitylene::liquidDiffCoeff(temperature, pressure);
335 const LhsEval& diffWCl = 1.e-9; // BinaryCoeff::H2O_Mesitylene::liquidDiffCoeff(temperature, pressure);
336 const LhsEval& diffAWl = 1.e-9; // BinaryCoeff::H2O_Air::liquidDiffCoeff(temperature, pressure);
337
338 const LhsEval& xwa = decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, airIdx));
339 const LhsEval& xww = decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, H2OIdx));
340 const LhsEval& xwc = decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, NAPLIdx));
341
342 switch (compIdx) {
343 case NAPLIdx:
344 diffCont = (1.- xww)/(xwa/diffAWl + xwc/diffWCl);
345 return diffCont;
346 case airIdx:
347 diffCont = (1.- xwc)/(xww/diffWCl + xwa/diffACl);
348 return diffCont;
349 case H2OIdx:
350 throw std::logic_error("Diffusivity of water in the water phase "
351 "is constraint by sum of diffusive fluxes = 0 !\n");
352 };
353 }
354 else if (phaseIdx==naplPhaseIdx) {
355 throw std::logic_error("Diffusion coefficients of "
356 "substances in liquid phase are undefined!\n");
357 }
358 return 0;
359#endif
360 }
361
363 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
364 static LhsEval fugacityCoefficient(const FluidState& fluidState,
365 const ParameterCache<ParamCacheEval>& /*paramCache*/,
366 unsigned phaseIdx,
367 unsigned compIdx)
368 {
369 assert(phaseIdx < numPhases);
370 assert(compIdx < numComponents);
371
372 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
373 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
374 Valgrind::CheckDefined(T);
375 Valgrind::CheckDefined(p);
376
377 if (phaseIdx == waterPhaseIdx) {
378 if (compIdx == H2OIdx)
379 return H2O::vaporPressure(T)/p;
380 else if (compIdx == airIdx)
381 return BinaryCoeff::H2O_N2::henry(T)/p;
382 else if (compIdx == NAPLIdx)
384 assert(false);
385 }
386 // for the NAPL phase, we assume currently that nothing is
387 // dissolved. this means that the affinity of the NAPL
388 // component to the NAPL phase is much higher than for the
389 // other components, i.e. the fugacity cofficient is much
390 // smaller.
391 else if (phaseIdx == naplPhaseIdx) {
392 const LhsEval& phiNapl = NAPL::vaporPressure(T)/p;
393 if (compIdx == NAPLIdx)
394 return phiNapl;
395 else if (compIdx == airIdx)
396 return 1e6*phiNapl;
397 else if (compIdx == H2OIdx)
398 return 1e6*phiNapl;
399 assert(false);
400 }
401
402 // for the gas phase, assume an ideal gas when it comes to
403 // fugacity (-> fugacity == partial pressure)
404 assert(phaseIdx == gasPhaseIdx);
405 return 1.0;
406 }
407
408
410 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
411 static LhsEval enthalpy(const FluidState& fluidState,
412 const ParameterCache<ParamCacheEval>& /*paramCache*/,
413 unsigned phaseIdx)
414 {
415 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
416 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
417
418 if (phaseIdx == waterPhaseIdx) {
419 return H2O::liquidEnthalpy(T, p);
420 }
421 else if (phaseIdx == naplPhaseIdx) {
422 return NAPL::liquidEnthalpy(T, p);
423 }
424 else if (phaseIdx == gasPhaseIdx) {
425 // gas phase enthalpy depends strongly on composition
426 LhsEval result = 0;
427 result += H2O::gasEnthalpy(T, p) * decay<LhsEval>(fluidState.massFraction(gasPhaseIdx, H2OIdx));
428 result += NAPL::gasEnthalpy(T, p) * decay<LhsEval>(fluidState.massFraction(gasPhaseIdx, airIdx));
429 result += Air::gasEnthalpy(T, p) * decay<LhsEval>(fluidState.massFraction(gasPhaseIdx, NAPLIdx));
430
431 return result;
432 }
433 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
434 }
435
437 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
438 static LhsEval thermalConductivity(const FluidState& fluidState,
439 const ParameterCache<ParamCacheEval>& /*paramCache*/,
440 unsigned phaseIdx)
441 {
442 assert(phaseIdx < numPhases);
443
444 if (phaseIdx == waterPhaseIdx){ // water phase
445 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
446 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
447
449 }
450 else if (phaseIdx == gasPhaseIdx) { // gas phase
451 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
452 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
453
454 return Air::gasThermalConductivity(T, p);
455 }
456
457 assert(phaseIdx == naplPhaseIdx);
458
459 // Taken from:
460 //
461 // D. K. H. Briggs: "Thermal Conductivity of Liquids",
462 // Ind. Eng. Chem., 1957, 49 (3), pp 418–421
463 //
464 // Convertion to SI units:
465 // 344e-6 cal/(s cm K) = 0.0143964 J/(s m K)
466 return 0.0143964;
467 }
468};
469} // namespace Opm
470
471#endif
A generic class which tabulates all thermodynamic properties of a given component.
static Evaluation gasThermalConductivity(const Evaluation &, const Evaluation &)
Specific heat conductivity of steam .
Definition: Air.hpp:217
static Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of at a given pressure and temperature [kg/m^3].
Definition: Air.hpp:102
static Scalar molarMass()
The molar mass in of .
Definition: Air.hpp:80
static const char * name()
A human readable name for the .
Definition: Air.hpp:60
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition: Air.hpp:72
static Evaluation gasEnthalpy(const Evaluation &temperature, const Evaluation &)
Specific enthalpy of liquid water with 273.15 K as basis.
Definition: Air.hpp:180
static Evaluation gasViscosity(const Evaluation &temperature, const Evaluation &)
The dynamic viscosity of at a given pressure and temperature.
Definition: Air.hpp:137
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 gasDiffCoeff(Evaluation temperature, Evaluation pressure)
Binary diffusion coefficent [m^2/s] for air and mesitylene.
Definition: Air_Mesitylene.hpp:59
static Evaluation gasDiffCoeff(const Evaluation &temperature, const Evaluation &pressure)
Binary diffusion coefficent for molecular water and air.
Definition: H2O_Air.hpp:70
static Evaluation henry(const Evaluation &)
Henry coefficent for mesitylene in liquid water.
Definition: H2O_Mesitylene.hpp:52
static Evaluation gasDiffCoeff(Evaluation temperature, Evaluation pressure)
Binary diffusion coefficent [m^2/s] for molecular water and mesitylene.
Definition: H2O_Mesitylene.hpp:67
static Evaluation henry(const Evaluation &temperature)
Henry coefficent for molecular nitrogen in liquid water.
Definition: H2O_N2.hpp:52
A fluid system with water, gas and NAPL as phases and water, air and mesitylene (DNAPL) as components...
Definition: H2OAirMesityleneFluidSystem.hpp:55
static const int numComponents
Number of chemical species in the fluid system.
Definition: H2OAirMesityleneFluidSystem.hpp:81
static const int airIdx
The index of the air pseudo-component.
Definition: H2OAirMesityleneFluidSystem.hpp:95
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: H2OAirMesityleneFluidSystem.hpp:364
static const int H2OIdx
The index of the water component.
Definition: H2OAirMesityleneFluidSystem.hpp:91
static const int waterPhaseIdx
The index of the water phase.
Definition: H2OAirMesityleneFluidSystem.hpp:84
static const int naplPhaseIdx
The index of the NAPL phase.
Definition: H2OAirMesityleneFluidSystem.hpp:86
Mesitylene< Scalar > NAPL
The type of the mesithylene/napl component.
Definition: H2OAirMesityleneFluidSystem.hpp:68
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: H2OAirMesityleneFluidSystem.hpp:119
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition: H2OAirMesityleneFluidSystem.hpp:129
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition: H2OAirMesityleneFluidSystem.hpp:438
TabulatedH2O H2O
The type of the water component.
Definition: H2OAirMesityleneFluidSystem.hpp:75
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: H2OAirMesityleneFluidSystem.hpp:140
static const int NAPLIdx
The index of the NAPL component.
Definition: H2OAirMesityleneFluidSystem.hpp:93
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition: H2OAirMesityleneFluidSystem.hpp:184
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition: H2OAirMesityleneFluidSystem.hpp:173
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: H2OAirMesityleneFluidSystem.hpp:199
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition: H2OAirMesityleneFluidSystem.hpp:162
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: H2OAirMesityleneFluidSystem.hpp:244
static LhsEval diffusionCoefficient(const FluidState &, const ParameterCache< ParamCacheEval > &, unsigned, unsigned)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition: H2OAirMesityleneFluidSystem.hpp:306
static const int numPhases
Number of fluid phases in the fluid system.
Definition: H2OAirMesityleneFluidSystem.hpp:79
static const int gasPhaseIdx
The index of the gas phase.
Definition: H2OAirMesityleneFluidSystem.hpp:88
static void init()
Initialize the fluid system's static parameters.
Definition: H2OAirMesityleneFluidSystem.hpp:98
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: H2OAirMesityleneFluidSystem.hpp:152
::Opm::Air< Scalar > Air
The type of the air component.
Definition: H2OAirMesityleneFluidSystem.hpp:71
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: H2OAirMesityleneFluidSystem.hpp:136
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: H2OAirMesityleneFluidSystem.hpp:411
Material properties of pure water .
Definition: H2O.hpp:65
Component for Mesitylene.
Definition: Mesitylene.hpp:45
static bool liquidIsCompressible()
Returns true iff the liquid phase is assumed to be compressible.
Definition: Mesitylene.hpp:218
static Evaluation vaporPressure(const Evaluation &temperature)
The saturation vapor pressure in of pure mesitylene at a given temperature according to Antoine afte...
Definition: Mesitylene.hpp:99
static Evaluation liquidEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of liquid mesitylene .
Definition: Mesitylene.hpp:118
static Scalar molarMass()
The molar mass in of mesitylene.
Definition: Mesitylene.hpp:58
static Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of pure mesitylene vapor at a given pressure and temperature .
Definition: Mesitylene.hpp:190
static const char * name()
A human readable name for the mesitylene.
Definition: Mesitylene.hpp:52
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition: Mesitylene.hpp:212
static Evaluation gasEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of mesitylene vapor .
Definition: Mesitylene.hpp:178
static Evaluation liquidDensity(const Evaluation &temperature, const Evaluation &)
The density of pure mesitylene at a given pressure and temperature .
Definition: Mesitylene.hpp:200
static Evaluation liquidViscosity(Evaluation temperature, const Evaluation &)
The dynamic viscosity of pure mesitylene.
Definition: Mesitylene.hpp:255
static Evaluation gasViscosity(Evaluation temperature, const Evaluation &, bool=true)
The dynamic viscosity of mesitylene vapor.
Definition: Mesitylene.hpp:229
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 gasEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the gas .
Definition: TabulatedComponent.hpp:282
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the tables.
Definition: TabulatedComponent.hpp:72
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 Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of gas at a given pressure and temperature .
Definition: TabulatedComponent.hpp:426
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition: TabulatedComponent.hpp:414
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 gasViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of gas.
Definition: TabulatedComponent.hpp:461
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
Definition: H2OAirMesityleneFluidSystem.hpp:65
Definition: MathToolbox.hpp:50