All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
H2OAirXyleneFluidSystem.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_XYLENE_FLUID_SYSTEM_HPP
28 #define OPM_H2O_AIR_XYLENE_FLUID_SYSTEM_HPP
29 
35 
39 
40 #include "BaseFluidSystem.hpp"
41 #include "NullParameterCache.hpp"
42 
43 namespace Opm {
44 namespace FluidSystems {
45 
51 template <class Scalar>
53  : public BaseFluidSystem<Scalar, H2OAirXylene<Scalar> >
54 {
57 
58 public:
59  template <class Evaluation>
60  struct ParameterCache : public Opm::NullParameterCache<Evaluation>
61  {};
62 
69 
71  static const int numPhases = 3;
73  static const int numComponents = 3;
74 
76  static const int waterPhaseIdx = 0;
78  static const int naplPhaseIdx = 1;
80  static const int gasPhaseIdx = 2;
81 
83  static const int H2OIdx = 0;
85  static const int NAPLIdx = 1;
87  static const int airIdx = 2;
88 
90  static void init()
91  { }
92 
94  static bool isLiquid(unsigned phaseIdx)
95  {
96  //assert(0 <= phaseIdx && phaseIdx < numPhases);
97  return phaseIdx != gasPhaseIdx;
98  }
99 
101  static bool isIdealGas(unsigned phaseIdx)
102  { return phaseIdx == gasPhaseIdx && H2O::gasIsIdeal() && Air::gasIsIdeal() && NAPL::gasIsIdeal(); }
103 
105  static bool isIdealMixture(unsigned /*phaseIdx*/)
106  {
107  //assert(0 <= phaseIdx && phaseIdx < numPhases);
108 
109  // we assume Henry's and Rault's laws for the water phase and
110  // and no interaction between gas molecules of different
111  // components, so all phases are ideal mixtures!
112  return true;
113  }
114 
116  static bool isCompressible(unsigned phaseIdx)
117  {
118  return
119  (phaseIdx == gasPhaseIdx)
120  // gases are always compressible
121  ? true
122  : (phaseIdx == waterPhaseIdx)
123  // the water component decides for the water phase...
125  // the NAPL component decides for the napl phase...
127  }
128 
130  static const char* phaseName(unsigned phaseIdx)
131  {
132  switch (phaseIdx) {
133  case waterPhaseIdx: return "water";
134  case naplPhaseIdx: return "napl";
135  case gasPhaseIdx: return "gas";
136  };
137  OPM_THROW(std::logic_error, "Invalid phase index " << phaseIdx);
138  }
139 
141  static const char* componentName(unsigned compIdx)
142  {
143  switch (compIdx) {
144  case H2OIdx: return H2O::name();
145  case airIdx: return Air::name();
146  case NAPLIdx: return NAPL::name();
147  };
148  OPM_THROW(std::logic_error, "Invalid component index " << compIdx);
149  }
150 
152  static Scalar molarMass(unsigned compIdx)
153  {
154  return
155  (compIdx == H2OIdx)
156  // gases are always compressible
157  ? H2O::molarMass()
158  : (compIdx == airIdx)
159  // the water component decides for the water comp...
160  ? Air::molarMass()
161  // the NAPL component decides for the napl comp...
162  : (compIdx == NAPLIdx)
163  ? NAPL::molarMass()
164  : 1e30;
165  }
166 
168  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
169  static LhsEval density(const FluidState& fluidState,
170  const ParameterCache<ParamCacheEval>& /*paramCache*/,
171  unsigned phaseIdx)
172  {
173  if (phaseIdx == waterPhaseIdx) {
174  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
175  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
176 
177  // See: Ochs 2008
178  // \todo: proper citation
179  const LhsEval& rholH2O = H2O::liquidDensity(T, p);
180  const LhsEval& clH2O = rholH2O/H2O::molarMass();
181 
182  const auto& xwH2O = Opm::decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, H2OIdx));
183  const auto& xwAir = Opm::decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, airIdx));
184  const auto& xwNapl = Opm::decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, NAPLIdx));
185  // this assumes each dissolved molecule displaces exactly one
186  // water molecule in the liquid
187  return clH2O*(H2O::molarMass()*xwH2O + Air::molarMass()*xwAir + NAPL::molarMass()*xwNapl);
188  }
189  else if (phaseIdx == naplPhaseIdx) {
190  // assume pure NAPL for the NAPL phase
191  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
192  return NAPL::liquidDensity(T, LhsEval(1e30));
193  }
194 
195  assert (phaseIdx == gasPhaseIdx);
196  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
197  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
198 
199  const LhsEval& pH2O = Opm::decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, H2OIdx))*p;
200  const LhsEval& pAir = Opm::decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, airIdx))*p;
201  const LhsEval& pNAPL = Opm::decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, NAPLIdx))*p;
202  return
203  H2O::gasDensity(T, pH2O) +
204  Air::gasDensity(T, pAir) +
205  NAPL::gasDensity(T, pNAPL);
206  }
207 
209  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
210  static LhsEval viscosity(const FluidState& fluidState,
211  const ParameterCache<ParamCacheEval>& /*paramCache*/,
212  unsigned phaseIdx)
213  {
214  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
215  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
216 
217  if (phaseIdx == waterPhaseIdx) {
218  // assume pure water viscosity
219  return H2O::liquidViscosity(T, p);
220  }
221  else if (phaseIdx == naplPhaseIdx) {
222  // assume pure NAPL viscosity
223  return NAPL::liquidViscosity(T, p);
224  }
225 
226  assert (phaseIdx == gasPhaseIdx);
227 
228  /* Wilke method. See:
229  *
230  * See: R. Reid, et al.: The Properties of Gases and Liquids,
231  * 4th edition, McGraw-Hill, 1987, 407-410
232  * 5th edition, McGraw-Hill, 20001, p. 9.21/22
233  *
234  * in this case, we use a simplified version in order to avoid
235  * computationally costly evaluation of sqrt and pow functions and
236  * divisions
237  * -- compare e.g. with Promo Class p. 32/33
238  */
239  const LhsEval mu[numComponents] = {
241  Air::simpleGasViscosity(T, p),
243  };
244  // molar masses
245  const Scalar M[numComponents] = {
246  H2O::molarMass(),
247  Air::molarMass(),
249  };
250 
251  const auto& xgAir = Opm::decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, airIdx));
252  const auto& xgH2O = Opm::decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, H2OIdx));
253  const auto& xgNapl = Opm::decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, NAPLIdx));
254 
255  const LhsEval& xgAW = xgAir + xgH2O;
256  const LhsEval& muAW = (mu[airIdx]*xgAir + mu[H2OIdx]*xgH2O)/ xgAW;
257 
258  const LhsEval& MAW = (xgAir*Air::molarMass() + xgH2O*H2O::molarMass())/xgAW;
259 
260  Scalar phiCAW = 0.3; // simplification for this particular system
261  /* actually like this
262  * Scalar phiCAW = std::pow(1.+std::sqrt(mu[NAPLIdx]/muAW)*std::pow(MAW/M[NAPLIdx],0.25),2)
263  * / std::sqrt(8.*(1.+M[NAPLIdx]/MAW));
264  */
265  const LhsEval& phiAWC = phiCAW * muAW*M[NAPLIdx]/(mu[NAPLIdx]*MAW);
266 
267  return (xgAW*muAW)/(xgAW+xgNapl*phiAWC) + (xgNapl*mu[NAPLIdx])/(xgNapl + xgAW*phiCAW);
268  }
269 
271  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
272  static LhsEval diffusionCoefficient(const FluidState& fluidState,
273  const ParameterCache<ParamCacheEval>& /*paramCache*/,
274  unsigned phaseIdx,
275  unsigned compIdx)
276  {
277  if (phaseIdx==gasPhaseIdx) {
278  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
279  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
280 
281  const LhsEval& diffAC = Opm::BinaryCoeff::Air_Xylene::gasDiffCoeff(T, p);
282  const LhsEval& diffWC = Opm::BinaryCoeff::H2O_Xylene::gasDiffCoeff(T, p);
283  const LhsEval& diffAW = Opm::BinaryCoeff::H2O_Air::gasDiffCoeff(T, p);
284 
285  const LhsEval& xga = Opm::decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, airIdx));
286  const LhsEval& xgw = Opm::decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, H2OIdx));
287  const LhsEval& xgc = Opm::decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, NAPLIdx));
288 
289  if (compIdx==NAPLIdx) return (1.- xgw)/(xga/diffAW + xgc/diffWC);
290  else if (compIdx==H2OIdx) return (1.- xgc)/(xgw/diffWC + xga/diffAC);
291  else if (compIdx==airIdx) OPM_THROW(std::logic_error,
292  "Diffusivity of air in the gas phase "
293  "is constraint by sum of diffusive fluxes = 0 !\n");
294  } else if (phaseIdx==waterPhaseIdx){
295  Scalar diffACl = 1.e-9; // BinaryCoeff::Air_Xylene::liquidDiffCoeff(temperature, pressure);
296  Scalar diffWCl = 1.e-9; // BinaryCoeff::H2O_Xylene::liquidDiffCoeff(temperature, pressure);
297  Scalar diffAWl = 1.e-9; // BinaryCoeff::H2O_Air::liquidDiffCoeff(temperature, pressure);
298 
299  const LhsEval& xwa = Opm::decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, airIdx));
300  const LhsEval& xww = Opm::decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, H2OIdx));
301  const LhsEval& xwc = Opm::decay<LhsEval>(fluidState.moleFraction(waterPhaseIdx, NAPLIdx));
302 
303  switch (compIdx) {
304  case NAPLIdx:
305  return (1.- xww)/(xwa/diffAWl + xwc/diffWCl);
306  case airIdx:
307  return (1.- xwc)/(xww/diffWCl + xwa/diffACl);
308  case H2OIdx:
309  OPM_THROW(std::logic_error,
310  "Diffusivity of water in the water phase "
311  "is constraint by sum of diffusive fluxes = 0 !\n");
312  };
313  } else if (phaseIdx==naplPhaseIdx) {
314 
315  OPM_THROW(std::logic_error,
316  "Diffusion coefficients of "
317  "substances in liquid phase are undefined!\n");
318  }
319  return 0;
320  }
321 
323  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
324  static LhsEval fugacityCoefficient(const FluidState& fluidState,
325  const ParameterCache<ParamCacheEval>& /*paramCache*/,
326  unsigned phaseIdx,
327  unsigned compIdx)
328  {
329  assert(0 <= phaseIdx && phaseIdx < numPhases);
330  assert(0 <= compIdx && compIdx < numComponents);
331 
332  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
333  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
334 
335  if (phaseIdx == waterPhaseIdx) {
336  if (compIdx == H2OIdx)
337  return H2O::vaporPressure(T)/p;
338  else if (compIdx == airIdx)
340  else if (compIdx == NAPLIdx)
342  }
343 
344  // for the NAPL phase, we assume currently that nothing is
345  // dissolved. this means that the affinity of the NAPL
346  // component to the NAPL phase is much higher than for the
347  // other components, i.e. the fugacity cofficient is much
348  // smaller.
349  if (phaseIdx == naplPhaseIdx) {
350  const LhsEval& phiNapl = NAPL::vaporPressure(T)/p;
351  if (compIdx == NAPLIdx)
352  return phiNapl;
353  else if (compIdx == airIdx)
354  return 1e6*phiNapl;
355  else if (compIdx == H2OIdx)
356  return 1e6*phiNapl;
357  }
358 
359  // for the gas phase, assume an ideal gas when it comes to
360  // fugacity (-> fugacity == partial pressure)
361  assert(phaseIdx == gasPhaseIdx);
362  return 1.0;
363  }
364 
366  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
367  static LhsEval enthalpy(const FluidState& fluidState,
368  const ParameterCache<ParamCacheEval>& /*paramCache*/,
369  unsigned phaseIdx)
370  {
371  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
372  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
373 
374  if (phaseIdx == waterPhaseIdx) {
375  return H2O::liquidEnthalpy(T, p);
376  }
377  else if (phaseIdx == naplPhaseIdx) {
378  return NAPL::liquidEnthalpy(T, p);
379  }
380  else if (phaseIdx == gasPhaseIdx) { // gas phase enthalpy depends strongly on composition
381  const LhsEval& hgc = NAPL::gasEnthalpy(T, p);
382  const LhsEval& hgw = H2O::gasEnthalpy(T, p);
383  const LhsEval& hga = Air::gasEnthalpy(T, p);
384 
385  LhsEval result = 0;
386  result += hgw * Opm::decay<LhsEval>(fluidState.massFraction(gasPhaseIdx, H2OIdx));
387  result += hga * Opm::decay<LhsEval>(fluidState.massFraction(gasPhaseIdx, airIdx));
388  result += hgc * Opm::decay<LhsEval>(fluidState.massFraction(gasPhaseIdx, NAPLIdx));
389 
390  return result;
391  }
392  OPM_THROW(std::logic_error, "Invalid phase index " << phaseIdx);
393  }
394 
395 private:
396  template <class LhsEval>
397  static LhsEval waterPhaseDensity_(const LhsEval& T,
398  const LhsEval& pw,
399  const LhsEval& xww,
400  const LhsEval& xwa,
401  const LhsEval& xwc)
402  {
403  const LhsEval& rholH2O = H2O::liquidDensity(T, pw);
404  const LhsEval& clH2O = rholH2O/H2O::molarMass();
405 
406  // this assumes each dissolved molecule displaces exactly one
407  // water molecule in the liquid
408  return clH2O*(xww*H2O::molarMass() + xwa*Air::molarMass() + xwc*NAPL::molarMass());
409  }
410 
411  template <class LhsEval>
412  static LhsEval gasPhaseDensity_(const LhsEval& T,
413  const LhsEval& pg,
414  const LhsEval& xgw,
415  const LhsEval& xga,
416  const LhsEval& xgc)
417  { return H2O::gasDensity(T, pg*xgw) + Air::gasDensity(T, pg*xga) + NAPL::gasDensity(T, pg*xgc); }
418 
419  template <class LhsEval>
420  static LhsEval NAPLPhaseDensity_(const LhsEval& T, const LhsEval& pn)
421  { return NAPL::liquidDensity(T, pn); }
422 };
423 } // namespace FluidSystems
424 } // namespace Opm
425 
426 #endif
Material properties of pure water .
Definition: H2O.hpp:61
static Evaluation liquidDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of pure water in at a given pressure and temperature.
Definition: H2O.hpp:678
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition: Xylene.hpp:277
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: H2OAirXyleneFluidSystem.hpp:116
static const int waterPhaseIdx
The index of the water phase.
Definition: H2OAirXyleneFluidSystem.hpp:76
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition: H2O.hpp:618
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: H2OAirXyleneFluidSystem.hpp:210
static const int numComponents
Number of chemical species in the fluid system.
Definition: H2OAirXyleneFluidSystem.hpp:73
static const int numPhases
Number of fluid phases in the fluid system.
Definition: H2OAirXyleneFluidSystem.hpp:71
static void init()
Initialize the fluid system&#39;s static parameters.
Definition: H2OAirXyleneFluidSystem.hpp:90
Relations valid for an ideal gas.
static const int naplPhaseIdx
The index of the NAPL phase.
Definition: H2OAirXyleneFluidSystem.hpp:78
static Evaluation gasEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the pure component in gas.
Definition: Xylene.hpp:211
static bool liquidIsCompressible()
Returns true iff the liquid phase is assumed to be compressible.
Definition: H2O.hpp:537
static Evaluation vaporPressure(const Evaluation &temperature)
The saturation vapor pressure in of pure xylene at a given temperature according to Antoine after Be...
Definition: Xylene.hpp:95
static Evaluation gasViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of steam.
Definition: H2O.hpp:779
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: H2OAirXyleneFluidSystem.hpp:367
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition: H2OAirXyleneFluidSystem.hpp:152
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition: H2OAirXyleneFluidSystem.hpp:130
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition: Air.hpp:73
static Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of steam in at a given pressure and temperature.
Definition: H2O.hpp:553
static Evaluation liquidEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of liquid water .
Definition: H2O.hpp:228
Binary coefficients for water and xylene.
Opm::Xylene< Scalar > NAPL
The type of the xylene/napl component.
Definition: H2OAirXyleneFluidSystem.hpp:66
static Scalar molarMass()
The molar mass in of .
Definition: Air.hpp:81
static const char * name()
A human readable name for the xylene.
Definition: Xylene.hpp:55
static const int airIdx
The index of the air pseudo-component.
Definition: H2OAirXyleneFluidSystem.hpp:87
Opm::Air< Scalar > Air
The type of the air component.
Definition: H2OAirXyleneFluidSystem.hpp:68
static Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density in of the component at a given pressure in and temperature in .
Definition: Xylene.hpp:220
A parameter cache which does nothing.
A simple class implementing the fluid properties of air.
Definition: Air.hpp:47
static Evaluation henry(const Evaluation &temperature)
Henry coefficent for air in liquid water.
Definition: H2O_Air.hpp:55
Material properties of pure water .
static Evaluation gasViscosity(Evaluation temperature, const Evaluation &)
The dynamic viscosity of the pure component at a given pressure in and temperature in ...
Definition: Xylene.hpp:290
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 Evaluation liquidEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the pure component in liquid.
Definition: Xylene.hpp:149
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: H2OAirXyleneFluidSystem.hpp:169
Component for Xylene.
static const char * name()
A human readable name for the water.
Definition: H2O.hpp:74
A fluid system with water, gas and NAPL as phases and water, air and NAPL (contaminant) as components...
Definition: H2OAirXyleneFluidSystem.hpp:52
static bool liquidIsCompressible()
Returns true iff the liquid phase is assumed to be compressible.
Definition: Xylene.hpp:283
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 LhsEval diffusionCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx, unsigned compIdx)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition: H2OAirXyleneFluidSystem.hpp:272
static Evaluation liquidViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of pure water.
Definition: H2O.hpp:804
A simple class implementing the fluid properties of air.
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: H2OAirXyleneFluidSystem.hpp:105
static Evaluation vaporPressure(Evaluation temperature)
The vapor pressure in of pure water at a given temperature.
Definition: H2O.hpp:132
static Evaluation gasDiffCoeff(Evaluation temperature, Evaluation pressure)
Binary diffusion coefficent [m^2/s] for molecular water and xylene.
Definition: H2O_Xylene.hpp:66
static Evaluation gasEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of water steam .
Definition: H2O.hpp:177
static Evaluation liquidViscosity(Evaluation temperature, const Evaluation &)
The dynamic liquid viscosity of the pure component.
Definition: Xylene.hpp:309
static const int H2OIdx
The index of the water component.
Definition: H2OAirXyleneFluidSystem.hpp:83
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: H2OAirXyleneFluidSystem.hpp:101
static Scalar molarMass()
The molar mass in of xylene.
Definition: Xylene.hpp:61
A generic class which tabulates all thermodynamic properties of a given component.
The base class for all fluid systems.
Definition: BaseFluidSystem.hpp:43
static const int NAPLIdx
The index of the NAPL component.
Definition: H2OAirXyleneFluidSystem.hpp:85
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition: H2OAirXyleneFluidSystem.hpp:94
static Evaluation gasDiffCoeff(Evaluation temperature, Evaluation pressure)
Binary diffusion coefficent [m^2/s] for air and xylene.
Definition: Air_Xylene.hpp:57
static const char * name()
A human readable name for the .
Definition: Air.hpp:61
Component for Xylene.
Definition: Xylene.hpp:46
static Evaluation liquidDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of the liquid component at a given pressure in and temperature in . ...
Definition: Xylene.hpp:265
A parameter cache which does nothing.
Definition: NullParameterCache.hpp:39
static const Scalar molarMass()
The molar mass in of water.
Definition: H2O.hpp:80
Scalar Scalar
The type used for scalar quantities.
Definition: BaseFluidSystem.hpp:49
Binary coefficients for water and xylene.
Definition: H2OAirXyleneFluidSystem.hpp:60
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: H2OAirXyleneFluidSystem.hpp:324
Opm::H2O< Scalar > H2O
The type of the water component.
Definition: H2OAirXyleneFluidSystem.hpp:64
static Evaluation henry(const Evaluation &)
Henry coefficent for xylene in liquid water.
Definition: H2O_Xylene.hpp:51
Binary coefficients for water and nitrogen.
The base class for all fluid systems.
static const int gasPhaseIdx
The index of the gas phase.
Definition: H2OAirXyleneFluidSystem.hpp:80
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition: H2OAirXyleneFluidSystem.hpp:141