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