My Project
Xylene.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_XYLENE_HPP
28 #define OPM_XYLENE_HPP
29 
33 
35 
36 #include <cmath>
37 
38 namespace Opm {
45 template <class Scalar>
46 class Xylene : public Component<Scalar, Xylene<Scalar> >
47 {
48  typedef Constants<Scalar> Consts;
49  typedef ::Opm::IdealGas<Scalar> IdealGas;
50 
51 public:
55  static const char* name()
56  { return "xylene"; }
57 
61  static Scalar molarMass()
62  { return 0.106; }
63 
67  static Scalar criticalTemperature()
68  { return 617.1; }
69 
73  static Scalar criticalPressure()
74  { return 35.4e5; }
75 
79  static Scalar tripleTemperature()
80  { throw std::runtime_error("Not implemented: tripleTemperature for xylene"); }
81 
85  static Scalar triplePressure()
86  { throw std::runtime_error("Not implemented: triplePressure for xylene"); }
87 
94  template <class Evaluation>
95  static Evaluation vaporPressure(const Evaluation& temperature)
96  {
97  const Scalar A = 7.00909;;
98  const Scalar B = 1462.266;;
99  const Scalar C = 215.110;;
100 
101  return 100*1.334*pow(10.0, (A - (B/(temperature - 273.15 + C))));
102  }
103 
112  template <class Evaluation>
113  static Evaluation spHeatCapLiquidPhase(const Evaluation& temperature, const Evaluation& /*pressure*/)
114  {
115  Evaluation CH3,C6H5,H;
116  // after Reid et al. : Missenard group contrib. method (s. example 5-8)
117  // Xylene: C9H12 : 3* CH3 ; 1* C6H5 (phenyl-ring) ; -2* H (this was too much!)
118  // linear interpolation between table values [J/(mol K)]
119 
120  if(temperature < 298.0){ // take care: extrapolation for Temperature<273
121  H = 13.4 + 1.2*(temperature - 273.0)/25.0; // 13.4 + 1.2 = 14.6 = H(T=298K) i.e. interpolation of table values 273<T<298
122  CH3 = 40.0 + 1.6*(temperature - 273.0)/25.0; // 40 + 1.6 = 41.6 = CH3(T=298K)
123  C6H5 = 113.0 + 4.2*(temperature - 273.0)/25.0; // 113 + 4.2 = 117.2 = C6H5(T=298K)
124  }
125  else if(temperature < 323.0){
126  H = 14.6 + 0.9*(temperature - 298.0)/25.0; // i.e. interpolation of table values 298<T<323
127  CH3 = 41.6 + 1.9*(temperature - 298.0)/25.0;
128  C6H5 = 117.2 + 6.2*(temperature - 298.0)/25.0;
129  }
130  else if(temperature < 348.0){
131  H = 15.5 + 1.2*(temperature - 323.0)/25.0; // i.e. interpolation of table values 323<T<348
132  CH3 = 43.5 + 2.3*(temperature - 323.0)/25.0;
133  C6H5 = 123.4 + 6.3*(temperature - 323.0)/25.0;
134  }
135  else {
136  H = 16.7 + 2.1*(temperature - 348.0)/25.0; // i.e. interpolation of table values 348<T<373
137  CH3 = 45.8 + 2.5*(temperature - 348.0)/25.0; // take care: extrapolation for Temperature>373
138  C6H5 = 129.7 + 6.3*(temperature - 348.0)/25.0; // most likely leads to underestimation
139  }
140 
141  return (C6H5 + 2*CH3 - H)/molarMass();// J/(mol K) -> J/(kg K)
142  }
143 
144 
148  template <class Evaluation>
149  static Evaluation liquidEnthalpy(const Evaluation& temperature, const Evaluation& pressure)
150  {
151  // Gauss quadrature rule:
152  // Interval: [0K; temperature (K)]
153  // Gauss-Legendre-Integration with variable transformation:
154  // \int_a^b f(T) dT \approx (b-a)/2 \sum_i=1^n \alpha_i f( (b-a)/2 x_i + (a+b)/2 )
155  // with: n=2, legendre -> x_i = +/- \sqrt(1/3), \apha_i=1
156  // here: a=0, b=actual temperature in Kelvin
157  // \leadsto h(T) = \int_0^T c_p(T) dT
158  // \approx 0.5 T * (cp( (0.5-0.5*\sqrt(1/3)) T) + cp((0.5+0.5*\sqrt(1/3)) T))
159  // = 0.5 T * (cp(0.2113 T) + cp(0.7887 T) )
160 
161  // enthalpy may have arbitrary reference state, but the empirical/fitted heatCapacity function needs Kelvin as input
162  return 0.5*temperature*(spHeatCapLiquidPhase(Evaluation(0.2113*temperature),pressure)
163  + spHeatCapLiquidPhase(Evaluation(0.7887*temperature),pressure));
164  }
165 
169  static Scalar boilingTemperature()
170  { return 412.3; }
171 
180  template <class Evaluation>
181  static Evaluation heatVap(Evaluation temperature,
182  const Evaluation& /*pressure*/)
183  {
184  temperature = min(temperature, criticalTemperature()); // regularization
185  temperature = max(temperature, 0.0); // regularization
186 
187  const Scalar T_crit = criticalTemperature();
188  const Scalar Tr1 = boilingTemperature()/criticalTemperature();
189  const Scalar p_crit = criticalPressure();
190 
191  // Chen method, eq. 7-11.4 (at boiling)
192  const Scalar DH_v_boil = Consts::R * T_crit * Tr1
193  * (3.978 * Tr1 - 3.958 + 1.555*std::log(p_crit * 1e-5 /*Pa->bar*/ ) )
194  / (1.07 - Tr1); /* [J/mol] */
195 
196  /* Variation with temperature according to Watson relation eq 7-12.1*/
197  const Evaluation& Tr2 = temperature/criticalTemperature();
198  const Scalar n = 0.375;
199  const Evaluation& DH_vap = DH_v_boil * pow(((1.0 - Tr2)/(1.0 - Tr1)), n);
200 
201  return (DH_vap/molarMass()); // we need [J/kg]
202  }
203 
210  template <class Evaluation>
211  static Evaluation gasEnthalpy(const Evaluation& temperature, const Evaluation& pressure)
212  {
213  return liquidEnthalpy(temperature, pressure) + heatVap(temperature, pressure);
214  }
215 
219  template <class Evaluation>
220  static Evaluation gasDensity(const Evaluation& temperature, const Evaluation& pressure)
221  {
222  return IdealGas::density(Evaluation(molarMass()), temperature, pressure);
223  }
224 
231  template <class Evaluation>
232  static Evaluation molarGasDensity(const Evaluation& temperature, const Evaluation& pressure)
233  {
234  return gasDensity(temperature, pressure) / molarMass();
235  }
236 
246  template <class Evaluation>
247  static Evaluation molarLiquidDensity(Evaluation temperature, const Evaluation& /*pressure*/)
248  {
249  // saturated molar volume according to Lide, CRC Handbook of
250  // Thermophysical and Thermochemical Data, CRC Press, 1994
251  // valid for 245 < Temperature < 600
252  temperature = min(temperature, 500.0); // regularization
253  temperature = max(temperature, 250.0); // regularization
254 
255  const Scalar A1 = 0.25919; // from table
256  const Scalar A2 = 0.0014569; // from table
257  const Evaluation& expo = 1.0 + pow((1.0 - temperature/criticalTemperature()), (2.0/7.0));
258  return 1.0/(A2*pow(A1, expo));
259  }
260 
264  template <class Evaluation>
265  static Evaluation liquidDensity(const Evaluation& temperature, const Evaluation& pressure)
266  { return molarLiquidDensity(temperature, pressure)*molarMass(); }
267 
271  static bool gasIsCompressible()
272  { return true; }
273 
277  static bool gasIsIdeal()
278  { return true; }
279 
283  static bool liquidIsCompressible()
284  { return false; }
285 
289  template <class Evaluation>
290  static Evaluation gasViscosity(Evaluation temperature, const Evaluation& /*pressure*/)
291  {
292  temperature = min(temperature, 500.0); // regularization
293  temperature = max(temperature, 250.0); // regularization
294 
295  const Evaluation& Tr = max(temperature/criticalTemperature(), 1e-10);
296  const Scalar Fp0 = 1.0;
297  const Scalar xi = 0.004623;
298  const Evaluation& eta_xi = Fp0*(0.807*pow(Tr, 0.618)
299  - 0.357*exp(-0.449*Tr)
300  + 0.34*exp(-4.058*Tr)
301  + 0.018);
302  return eta_xi/xi / 1e7; // [Pa s]
303  }
304 
308  template <class Evaluation>
309  static Evaluation liquidViscosity(Evaluation temperature, const Evaluation& /*pressure*/)
310  {
311  temperature = min(temperature, 500.0); // regularization
312  temperature = max(temperature, 250.0); // regularization
313 
314  const Scalar A = -3.82;
315  const Scalar B = 1027.0;
316  const Scalar C = -6.38e-4;
317  const Scalar D = 4.52e-7;
318 
319  return 1e-3*exp(A
320  + B/temperature
321  + C*temperature
322  + D*temperature*temperature); // in [cP]
323  }
324 };
325 
326 } // namespace Opm
327 
328 #endif
Abstract base class of a pure chemical species.
A central place for various physical constants occuring in some equations.
Relations valid for an ideal gas.
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
Abstract base class of a pure chemical species.
Definition: Component.hpp:42
A central place for various physical constants occuring in some equations.
Definition: Constants.hpp:41
static const Scalar R
The ideal gas constant [J/(mol K)].
Definition: Constants.hpp:45
Relations valid for an ideal gas.
Definition: IdealGas.hpp:38
static Evaluation density(const Evaluation &avgMolarMass, const Evaluation &temperature, const Evaluation &pressure)
The density of the gas in , depending on pressure, temperature and average molar mass of the gas.
Definition: IdealGas.hpp:48
Component for Xylene.
Definition: Xylene.hpp:47
static Evaluation liquidViscosity(Evaluation temperature, const Evaluation &)
The dynamic liquid viscosity of the pure component.
Definition: Xylene.hpp:309
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
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 gasEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the pure component in gas.
Definition: Xylene.hpp:211
static const char * name()
A human readable name for the xylene.
Definition: Xylene.hpp:55
static bool liquidIsCompressible()
Returns true iff the liquid phase is assumed to be compressible.
Definition: Xylene.hpp:283
static Evaluation molarLiquidDensity(Evaluation temperature, const Evaluation &)
The molar density of pure xylene at a given pressure and temperature .
Definition: Xylene.hpp:247
static Evaluation spHeatCapLiquidPhase(const Evaluation &temperature, const Evaluation &)
Specific heat cap of liquid xylene .
Definition: Xylene.hpp:113
static Scalar criticalTemperature()
Returns the critical temperature of xylene.
Definition: Xylene.hpp:67
static Evaluation liquidEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the pure component in liquid.
Definition: Xylene.hpp:149
static Scalar tripleTemperature()
Returns the temperature at xylene's triple point.
Definition: Xylene.hpp:79
static Scalar molarMass()
The molar mass in of xylene.
Definition: Xylene.hpp:61
static Evaluation heatVap(Evaluation temperature, const Evaluation &)
Latent heat of vaporization for xylene .
Definition: Xylene.hpp:181
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 molarGasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of xylene gas at a given pressure and temperature.
Definition: Xylene.hpp:232
static Scalar triplePressure()
Returns the pressure at xylene's triple point.
Definition: Xylene.hpp:85
static Scalar boilingTemperature()
Returns the temperature at xylene's boiling point (1 atm).
Definition: Xylene.hpp:169
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
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition: Xylene.hpp:277
static bool gasIsCompressible()
Returns true iff the gas phase is assumed to be compressible.
Definition: Xylene.hpp:271
static Scalar criticalPressure()
Returns the critical pressure of xylene.
Definition: Xylene.hpp:73