Cantera  3.2.0a2
Loading...
Searching...
No Matches
HighPressureGasTransport.h
Go to the documentation of this file.
1/**
2 * @file HighPressureGasTransport.h
3 * Interface for class HighPressureGasTransport
4 */
5
6// This file is part of Cantera. See License.txt in the top-level directory or
7// at https://cantera.org/license.txt for license and copyright information.
8
9#ifndef CT_HIGHPRESSUREGASTRAN_H
10#define CT_HIGHPRESSUREGASTRAN_H
11
12// Cantera includes
13#include "GasTransport.h"
15
16namespace Cantera
17{
18
19/**
20 * Base class for high pressure gas transport models. This class handles
21 * common functionality used by high pressure transport models, including
22 * critical property storage and binary diffusion coefficient correction.
23 *
24 * @ingroup tranprops
25 */
27{
28protected:
29 //! default constructor
31
32public:
33 /**
34 * Computes the matrix of binary diffusion coefficients [m²/s] using the Takahashi
35 * correction factor.
36 *
37 * @param[in] ld Leading dimension of the flattened array `d` used to store the
38 * diffusion coefficient matrix; usually equal to the number of
39 * species.
40 * @param[out] d Diffusion coefficient matrix stored in column-major (Fortran)
41 * order, such that @f$ \mathcal{D}_{ij} = \tt{d[ld*j + i]} @f$; must
42 * be at least `ld` times the number of species in length.
43 * @see GasTransport::fitDiffCoeffs()
44 */
45 void getBinaryDiffCoeffs(const size_t ld, double* const d) override;
46
47 /**
48 * Returns the mixture-averaged diffusion coefficients [m²/s].
49 *
50 * This method is the same as GasTransport::getMixDiffCoeffs() with the exception
51 * that the binary diffusion coefficients are multiplied by the Takahashi correction
52 * factor, which is described in takahashiCorrectionFactor() .
53 *
54 * @param[out] d Vector of mixture diffusion coefficients, @f$ D_{km}' @f$ ,
55 * for each species. length #m_nsp.
56 */
57 void getMixDiffCoeffs(double* const d) override;
58
59 /**
60 * Returns the mixture-averaged diffusion coefficients [m²/s].
61 *
62 * This method is the same as GasTransport::getMixDiffCoeffsMole() with the exception
63 * that the binary diffusion coefficients are multiplied by the Takahashi correction
64 * factor, which is described in takahashiCorrectionFactor() .
65 *
66 * @param[out] d vector of mixture-averaged diffusion coefficients for
67 * each species, length #m_nsp.
68 */
69 void getMixDiffCoeffsMole(double* const d) override;
70
71 /**
72 * Returns the mixture-averaged diffusion coefficients [m²/s].
73 *
74 * This method is the same as GasTransport::getMixDiffCoeffsMass() with the exception
75 * that the binary diffusion coefficients are multiplied by the Takahashi correction
76 * factor, which is described in takahashiCorrectionFactor() .
77 *
78 * @param[out] d vector of mixture-averaged diffusion coefficients for
79 * each species, length #m_nsp.
80 */
81 void getMixDiffCoeffsMass(double* const d) override;
82
83 /**
84 * Updates the matrix of species-pair Takahashi correction factors for use in
85 * computing the binary diffusion coefficients, see takahashiCorrectionFactor()
86 */
87 virtual void updateCorrectionFactors();
88
89protected:
90 /**
91 * Obtain required parameters from the 'critical-parameters' species input section,
92 * and checks the critical-properties.yaml file if an acentric factor is not
93 * specified.
94 *
95 * The way that GasTransport parses the input file is that if an acentric
96 * factor is not specified, it is quietly set to 0.0. A user may have the proper
97 * acentric factor in the critical-properties.yaml file, so that file is checked if
98 * a zero value is present.
99 */
100 void getTransportData() override;
101
102 /**
103 * Computes and stores the estimate of the critical properties for each species.
104 *
105 * This method sets the species composition vector to unity for species i and zero
106 * for all other species, and then queries the thermo object for the critical
107 * properties and stores them. It then resets the composition vector to the original
108 * state. This method only needs to be called once, as the critical properties for
109 * the pure species do not change.
110 *
111 * All species must have critical properties defined in the input file, either via
112 * critical properties or by specific values of the equation of state that are
113 * not zero.
114 */
116
117 //! Returns the stored value of the critical temperature for species 'i'.
118 double Tcrit_i(size_t i);
119
120 //! Returns the stored value of the critical pressure for species 'i'.
121 double Pcrit_i(size_t i);
122
123 //! Returns the stored value of the critical volume for species 'i'.
124 double Vcrit_i(size_t i);
125
126 //! Returns the stored value of the critical compressibility for species 'i'.
127 double Zcrit_i(size_t i);
128
129 vector<double> m_Tcrit; //!< Critical temperature [K] of each species
130 vector<double> m_Pcrit; //!< Critical pressure [Pa] of each species
131 vector<double> m_Vcrit; //!< Critical volume [m³/kmol] of each species
132 vector<double> m_Zcrit; //!< Critical compressibility [unitless] of each species
133
134 //! Matrix of Takahashi binary diffusion coefficient corrections. Size is
135 //! #m_nsp x #m_nsp.
137
138 friend class TransportFactory;
139};
140
141/**
142 * The implementation employs a method of corresponding states, using the Takahashi
143 * @cite takahashi1975 approach for binary diffusion coefficients (using mixture
144 * averaging rules for the mixture properties), the Lucas method for the viscosity, and
145 * a method from Ely and Hanley for the thermal conductivity. All methods are described
146 * in Poling et al. @cite poling2001 (viscosity in Ch. 9, thermal conductivity in
147 * Ch. 10, and diffusion coefficients in Ch. 11).
148 *
149 * @ingroup tranprops
150 */
152{
153protected:
154 //! default constructor
156
157public:
158 string transportModel() const override {
159 return "high-pressure";
160 }
161
162 void init(ThermoPhase* thermo, int mode=0) override;
163
164 /**
165 * Returns the mixture high-pressure thermal conductivity [W/m/K]
166 * using a method of corresponding states by Ely and Hanley.
167 *
168 * This method for computing the thermal conductivity is described in
169 * @cite ely-hanley1983 . There are references to earlier work in
170 * @cite ely-hanley1981 in the description of the model for thermal conductivity.
171 * The equations referenced in this description primarily are from
172 * @cite ely-hanley1983 , and equations that come from @cite ely-hanley1981 are
173 * marked as such. The equations referenced here are the same as the ones from
174 * the papers from Ely and Hanley.
175 *
176 * The thermal conductivity is assumed to have two components: a
177 * translational/collisional part and an internal part that is related to the
178 * transfer of energy to internal degrees of freedom.
179 *
180 * @f[
181 * \lambda_{mix}(\rho, T) = \lambda^{''}_{mix}(T) + \lambda^{'}_{mix}(\rho, T)
182 * @f]
183 *
184 * The first term on the right-hand side is the internal part and the second term
185 * is the translational part. The internal part is only a function of temperature,
186 * and the translational part is a function of both temperature and density.
187 *
188 * The internal component of the thermal conductivity is defined by Equation 2
189 * in @cite ely-hanley1983 :
190 *
191 * @f[
192 * \lambda^{''}_{mix}(T) = \sum_i \sum_j X_i X_j \lambda^{''}_{ij}(T)
193 * @f]
194 *
195 * The mixing rule for combining pure-species internal thermal conductivity
196 * components is given by Equation 3 in @cite ely-hanley1983 :
197 *
198 * @f[
199 * (\lambda^{''}_{ij})^{-1} = 2[(\lambda^{''}_{i})^{-1}
200 * + (\lambda^{''}_{j})^{-1}]
201 * @f]
202 *
203 * The pure species internal thermal conductivity is given by Equation 1 in
204 * @cite ely-hanley1983 :
205 *
206 * @f[
207 * \lambda^{''}_{i} = \frac{\eta_i^*}{MW_i}f_{int} \left(C_{p,i}^0 - 2.5R \right)
208 * @f]
209 *
210 * Where @f$ \eta_i^* @f$ is the referred to as the dilute gas viscosity and is the
211 * component that is temperature dependent described in @cite ely-hanley1981
212 * (see elyHanleyDilutePureSpeciesViscosity() ),
213 * @f$ MW_i @f$ is the molecular weight of species i, @f$ f_{int} @f$ is a constant
214 * and is 1.32, @f$ C_{p,i}^0 @f$ is the ideal gas heat capacity of species i,
215 * and R is the gas constant (J/kmol/K).
216 *
217 * For the translational component of the thermal conductivity, Equation 5 in
218 * @cite ely-hanley1983 is used:
219 *
220 * @f[
221 * \lambda^{'}_{mix}(\rho, T) = \lambda^{'}_{0}(\rho_0, T_0) F_{\lambda}
222 * @f]
223 *
224 * Where @f$ \lambda^{'}_{0}(\rho_0, T_0) @f$ is the internal component of the
225 * thermal conductivity of a reference fluid that is at a reference temperature
226 * and density, @f$ T_0 @f$ and @f$ \rho_0 @f$. These reference conditions are not
227 * the same as the conditions that the mixture is at (@f$ T @f$ and @f$ \rho @f$).
228 * The subscript 0 refers to the reference fluid. The term @f$ F_{\lambda} @f$ is
229 * defined by Equation 6 in @cite ely-hanley1983 :
230 *
231 * @f[
232 * F_{\lambda} = \left( \frac{MW_0}{MW_m} \right)^{1/2} f_{m,0}^{1/2}
233 * h_{m,0}^{-2/3}
234 * @f]
235 *
236 * Where @f$ MW_0 @f$ is the molecular weight of the reference fluid, @f$ MW_m @f$
237 * is the molecular weight of the mixture.
238 *
239 * @note: @f$ MW_m @f$ is an effective mixture molecular weight and should not be
240 * confused with a mole-weighted average of the molecular weights of the
241 * species in the mixture.
242 *
243 * The terms @f$ f_{m,0} @f$ and @f$ h_{m,0} @f$ are called reducing ratios. The
244 * @f$ h_{m,0} @f$ quantity is defined by Equation 9 in @cite ely-hanley1983 :
245 *
246 * @f[
247 * h_{m,0} = \sum_i \sum_j X_i X_j h_{ij}
248 * @f]
249 *
250 * with @f$ h_{ij} @f$ defined by Equation 11 in @cite ely-hanley1983 :
251 *
252 * @f[
253 * h_{ij} = \frac{1}{8} \left( h_i^{1/3} + h_j^{1/3} \right)^3 (1 - l_{ij})
254 * @f]
255 *
256 * The variable @f$ l_{ij} @f$ is a binary interaction parameter and is assumed to
257 * be zero as is done by Ely and Hanley. The pure species reducing ratio @f$ h_i @f$
258 * is defined by Equation 13 in @cite ely-hanley1983 :
259 *
260 * @f[
261 * h_i = \frac{V_c^i}{V_c^0} \phi_i(T_R^i, V_R^i, \omega_i)
262 * @f]
263 *
264 * Where @f$ \phi_i @f$ is a shape factor of Leach and Leland (see phiShapeFactor()),
265 * @f$ T_R^i @f$ is the reduced temperature of species i, @f$ V_R^i @f$ is the
266 * reduced volume of species i, and @f$ \omega_i @f$ is the acentric factor of
267 * species i.
268 *
269 * At this point, the value of @f$ h_{m,0} @f$ can be computed. This
270 * leaves the other reducing ratio, @f$ f_{m,0} @f$, to be computed. The value of
271 * that reducing ratio is defined by Equation 8 in @cite ely-hanley1983 :
272 *
273 * @f[
274 * f_{m,0} = \frac{1}{h_{m,0}} \sum_i \sum_j X_i X_j f_{ij} h_{ij}
275 * @f]
276 *
277 *
278 * The value of @f$ h_{ij} @f$ is the same as the one defined earlier. The
279 * combining rule for @f$ f_{ij} @f$ is given by Equation 10 in
280 * @cite ely-hanley1983 :
281 *
282 * @f[
283 * f_{ij} = (f_i f_j)^{1/2} (1-\kappa_{ij})
284 * @f]
285 *
286 * @f$ \kappa_{ij} @f$ is binary interaction parameters and is assumed to be zero
287 * as was done in the work of Ely and Hanley. The species-specific value of
288 * @f$ f_i @f$ is defined by Equation 12 in @cite ely-hanley1983 :
289 *
290 * @f[
291 * f_i = \frac{T_c^i}{T_c^0} \theta_i(T_R^i, V_R^i, \omega_i)
292 * @f]
293 *
294 * Where @f$ \theta_i @f$ is a shape factor of Leach and Leland
295 * (see thetaShapeFactor()), @f$ T_c^i @f$ is the critical temperature of species
296 * i, and @f$ T_c^0 @f$ is the critical temperature of the reference fluid.
297 *
298 * The value of @f$ h_{m,0} @f$ can be computed from
299 * the equation that defined the value of @f$ f_{m,0} h_{m,0} @f$. The value of
300 * @f$ h_{m,0} @f$ must be computed first and then the value of @f$ f_{m,0} @f$
301 * can be computed.
302 *
303 * The expression for @f$ MW_m @f$ is somewhat complex, but keep in mind that all
304 * terms except @f$ MW_m @f$ are known at this point and so simple algebra can be
305 * used to isolate the molecular weight of the mixture. The expression for the
306 * mixture molecular weight is given by Equation 14 in @cite ely-hanley1983 :
307 *
308 * @f[
309 * MW_m^{1/2} f_{m,0}^{1/2} h_{m,0}^{-4/3} = \sum_i \sum_j X_i X_j MW_{ij}^{-1/2}
310 * f_{ij}^{1/2} h_{ij}^{-4/3}
311 * @f]
312 *
313 * where the mixing rule for the molecular weight of the mixture is given by
314 * Equation 15 in @cite ely-hanley1983 :
315 *
316 * @f[
317 * MW_{ij} = \frac{1}{2} (MW_i^{-1} + MW_j^{-1})
318 * @f]
319 *
320 * For clarity, if we assume the right-hand-side of the molecular weight mixture
321 * equation is A, then the mixture molecular weight is given by:
322 *
323 * @f[
324 * MW_m = A^{-2} f_{m,0} h_{m,0}^{-8/3}
325 * @f]
326 *
327 * All of the above descriptions allow for the calculation of @f$ F_{\lambda} @f$
328 * in the expression for the mixture value of the translational component of the
329 * thermal conductivity. The reference fluid value of the translational component
330 * of the thermal conductivity ( @f$ \lambda^{'}_{0}(\rho_0, T_0) @f$ ) is still
331 * needed, which requires the the temperature and density of the reference fluid.
332 *
333 * The reference fluid density is computed using Equation 7 in
334 * @cite ely-hanley1983 :
335 *
336 * @f[
337 * \rho_0 = \rho h_{m,0}
338 * @f]
339 *
340 * The reference fluid temperature is computed using Equation 7 in
341 * @cite ely-hanley1983 :
342 *
343 * @f[
344 * T_0 = \frac{T}{f_{m,0}}
345 * @f]
346 *
347 * The reference fluid translational component of thermal conductivity can now be
348 * computed using Equation 18 in @cite ely-hanley1983.
349 * See elyHanleyReferenceThermalConductivity().
350 *
351 * @note The correlations for the reference fluid viscosity return a value of
352 * micro-grams/cm/s and the parts of the thermal conductivity that utilize
353 * the correlation fit parameters give values in units of mW/m/K. Appropriate
354 * conversions were made in the relevant functions:
355 * elyHanleyDiluteReferenceViscosity() and
356 * elyHanleyReferenceThermalConductivity()
357 */
358 double thermalConductivity() override;
359
360 /**
361 * Returns the mixture high-pressure viscosity [Pa·s] using the Lucas method.
362 *
363 * This uses the approach described in chapter 9-7. In this method, the mixture
364 * viscosity at high pressure is computed using the pure-fluid high pressure
365 * relation of Lucas with the difference being that mixture values of the
366 * model parameters are used. These mixture values are computed using the mixing
367 * rules described in see computeMixtureParameters() .
368 *
369 * The mixture pseudo-critical temperature and pressure are calculated using
370 * Equations 9-5.18 and 9-5.19. The mixture molecular weight is computed using
371 * Equation 9-5.20. The mixture values of the low-pressure polarity and quantum
372 * correction factors are computed using Equations 9-5.21 and 9-5.22.
373 */
374 double viscosity() override;
375
376 friend class TransportFactory;
377
378protected:
379
380 /**
381 * Get the viscosity [Pa·s] of a pure species using the method of Ely and Hanley.
382 *
383 * Uses the method outlined in @cite ely-hanley1981 to compute the viscosity of a
384 * pure species.
385 *
386 * The primary equation is given by Equation 1 in @cite ely-hanley1981 :
387 *
388 * @f[
389 * \eta_i(\rho,T) = \eta_0(\rho_0, T_0) F_{\eta}
390 * @f]
391 *
392 * Where @f$ F_{\eta} @f$ is defined by Equation 2 in @cite ely-hanley1981 :
393 *
394 * @f[
395 * F_{\eta} = \left( \frac{MW_i}{MW_0} \right)^{1/2} f_{i,0}^{1/2} h_{i,0}^{-2/3}
396 * @f]
397 *
398 * The `0` subscripts refer to the reference fluid, which is methane in this case,
399 * and the `i` subscripts refer to the species of interest. No mixing rules
400 * need to be used here because this is for a pure species. As such, the terms
401 * @f$ f_{i,0} @f$ and @f$ h_{i,0} @f$ have simpler expressions. The value of
402 * @f$ f_{i,0} @f$ is given by Equation 7 in @cite ely-hanley1981 :
403 *
404 * @f[
405 * f_{i,0} = \frac{T_c^i}{T_c^0} \theta_i(T_R^i, V_R^i, \omega_i)
406 * @f]
407 *
408 * and the value of @f$ h_{i,0} @f$ is given by Equation 8 in @cite ely-hanley1981 :
409 *
410 * @f[
411 * h_{i,0} = \frac{V_c^i}{V_c^0} \phi_i(T_R^i, V_R^i, Z_{c,i}, \omega_i)
412 * @f]
413 *
414 * Where @f$ \theta_i @f$ and @f$ \phi_i @f$ are shape factors of Leach and Leland
415 * ( see thetaShapeFactor() and phiShapeFactor() ), @f$ T_c^i @f$ is the critical
416 * temperature of species i, @f$ T_c^0 @f$ is the critical temperature of the
417 * reference fluid, @f$ V_c^i @f$ is the critical volume of species i,
418 * @f$ V_c^0 @f$ is the critical volume of the reference fluid, @f$ Z_{c,i} @f$
419 * is the critical compressibility of species i, and @f$ \omega_i @f$
420 * is the acentric factor of species i.
421 *
422 * @param V Molar volume of the species
423 * @param Tc Critical temperature of the species
424 * @param Vc Critical volume of the species
425 * @param Zc Critical compressibility of the species
426 * @param acentric_factor Acentric factor of the species
427 * @param mw Molecular weight of the species
428 */
429 double elyHanleyDilutePureSpeciesViscosity(double V, double Tc, double Vc,
430 double Zc, double acentric_factor,
431 double mw);
432
433 /**
434 * Returns the theta shape factor of Leach and Leland for a pure species.
435 *
436 * The shape factor @f$ \theta_i @f$ is defined by Equation 11 in @cite ely-hanley1981
437 * as follows:
438 *
439 * @f[
440 * \theta_i(T_R^i, V_R^i, \omega_i) = 1 + (\omega_i - \omega_0)F(T_R^i, V_R^i)]
441 * @f]
442 *
443 * Where @f$ \omega_i @f$ is the acentric factor of species i, @f$ \omega_0 @f$ is
444 * the acentric factor of the reference fluid, and @f$ F(T_R^i, V_R^i) @f$ is a
445 * function of the reduced temperature and reduced volume of species i. The
446 * function @f$ F(T_R^i, V_R^i) @f$ is defined by Equation 13 in
447 * @cite ely-hanley1981 :
448 *
449 * @f[
450 * F(T_R^i, V_R^i) = a_1 + b_1 ln(T_+^i) + (c_1 + d_1/T_+^i) (V_+^i - 0.5)
451 * @f]
452 *
453 * Where @f$ T_+^i @f$ and @f$ V_+^i @f$ are limited values of the reduced
454 * temperature and pressure. The limited temperature is defined by Equation 15 in
455 * @cite ely-hanley1981 :
456 *
457 * @f[
458 * T_+^i = \text{min}(2, \text{max}(T_R^i, 0.5))
459 * @f]
460 *
461 * and the limited pressure is defined by Equation 16 in @cite ely-hanley1981 :
462 *
463 * @f[
464 * V_i^+ = \text{min}(2, \text{max}(V_R^i, 0.5))
465 * @f]
466 *
467 * The values of @f$ a_1 @f$, @f$ b_1 @f$, @f$ c_1 @f$, and @f$ d_1 @f$ are
468 * given in Table II of @cite ely-hanley1981. They are:
469 *
470 * @f[
471 * a_1 = 0.090569, b_1 = -0.862762, c_1 = 0.316636, d_1 = -0.465684
472 * @f]
473 *
474 * @param Tr Reduced temperature
475 * @param Vr Reduced volume
476 * @param acentric_factor Acentric factor
477 */
478 double thetaShapeFactor(double Tr, double Vr, double acentric_factor);
479
480
481 /**
482 * Returns the phi shape factor of Leach and Leland for a pure species.
483 *
484 * The shape factor @f$ \phi_i @f$ is defined by Equation 12 in
485 * @cite ely-hanley1981 as follows:
486 *
487 * @f[
488 * \phi_i(T_R^i, V_R^i, \omega_i, Z_c^i) = \frac{Z_c^0}{Z_c^i} [1 + (\omega_i -
489 * \omega_0)G(T_R^i, V_R^i)]
490 * @f]
491 *
492 * Where @f$ Z_c^0 @f$ is the critical compressibility of the reference fluid,
493 * @f$ Z_c^i @f$ is the critical compressibility of species i, @f$ \omega_0 @f$
494 * is the acentric factor of the reference fluid, and @f$ G(T_R^i, V_R^i) @f$ is a
495 * function of the reduced temperature and reduced volume of species i. The
496 * function @f$ G(T_R^i, V_R^i) @f$ is defined by Equation 14 in
497 * @cite ely-hanley1981 :
498 *
499 * @f[
500 * G(T_R^i, V_R^i) = a_2(V_+^i + b_2) + c_2(V_+^i + d_2)ln(T_+^i)
501 * @f]
502 *
503 * Where @f$ T_+^i @f$ and @f$ V_+^i @f$ are limited values of the reduced
504 * temperature and pressure. The limited temperature is defined by Equation 15 in
505 * @cite ely-hanley1981 :
506 *
507 * @f[
508 * T_+^i = \text{min}(2, \text{max}(T_R^i, 0.5))
509 * @f]
510 *
511 * and the limited pressure is defined by Equation 16 in @cite ely-hanley1981 :
512 *
513 * @f[
514 * V_i^+ = \text{min}(2, \text{max}(V_R^i, 0.5))
515 * @f]
516 *
517 * The values of @f$ a_2 @f$, @f$ b_2 @f$, @f$ c_2 @f$, and @f$ d_2 @f$ are
518 * given in Table II of @cite ely-hanley1981. They are:
519 *
520 * @f[
521 * a_2 = 0.394901, b_2 = -1.023545, c_2 = -0.932813, d_2 = -0.754639
522 * @f]
523 *
524 * @param Tr Reduced temperature
525 * @param Vr Reduced volume
526 * @param Zc Critical compressibility
527 * @param acentric_factor Acentric factor
528 */
529 double phiShapeFactor(double Tr, double Vr, double Zc, double acentric_factor);
530
531 /**
532 * Returns the viscosity [Pa·s] for the reference fluid (methane) for low pressures.
533 *
534 * Computes the temperature dependent (referred to as the dilute viscosity in the
535 * reference) component only (eta_0) from the expression in Table III in
536 * @cite ely-hanley1981 . Prevents inputs larger than 10,000 Kelvin by just
537 * returning the value at 10,000 Kelvin.
538 *
539 * @param T0 Temperature of the reference fluid
540 */
541 double elyHanleyDiluteReferenceViscosity(double T0);
542
543 /**
544 * Returns the thermal conductivity [W/m/K] of the reference fluid of methane
545 *
546 * Computes the temperature and density dependent reference fluid thermal
547 * conductivity from the expression in Table I in @cite ely-hanley1983 .
548 *
549 * The reference fluid translational component of thermal conductivity is computed
550 * using Equation 18 in @cite ely-hanley1983 :
551 *
552 * @f[
553 * \lambda^{'}_{0}(\rho_0, T_0) = \frac{15R}{4MW_0} \eta_0^* +
554 * \lambda_0^{(1)}\rho_0 +
555 * \Delta\lambda_0(\rho_0, T_0)
556 * @f]
557 *
558 * Where @f$ \eta_0^* @f$ is the dilute reference gas viscosity,
559 * @f$ \lambda_0^{(1)} @f$ is the first density correction to the thermal
560 * conductivity, and @f$ \Delta\lambda_0 @f$ is the high-density contribution.
561 * Ely and Hanley provide a correlation equation to solve for this quantity for
562 * the methane reference fluid. The details of the correlation and the parameters
563 * are shown in Table I of @cite ely-hanley1983.
564 *
565 * For completeness, the relations are reproduced here.
566 *
567 * @f[
568 * \lambda_0^*(T_0) = \frac{15R}{4MW_0} \eta_0^*
569 * @f]
570 *
571 * @f[
572 * \eta_0^* = \sum_{n=1}^{9} C_n T_0^{(n-4)/3}
573 * @f]
574 *
575 * @f[
576 * \lambda_0^{(1)} = b_1 + b_2[b_3 - ln(T_0 / b_4)]^2
577 * @f]
578 *
579 * @f[
580 * \Delta\lambda_0 = \text{exp}[a_1 + a_2/T_0]
581 * (exp[(a_3 + a_4/T_0^{3/4})\rho_0^{0.1}
582 * + (\rho_0 / \rho_0,c - 1)\rho_0^{0.5}
583 * (a_5 + a_6/T_0 + a_7/T_0^2) ] - 1)
584 * @f]
585 *
586 * The constants a, b, and C are not related to any ones used in the previous
587 * equations, their values are defined in Table I of @cite ely-hanley1983.
588 *
589 * @param rho0 Density of the reference fluid
590 * @param T0 Temperature of the reference fluid
591 */
592 double elyHanleyReferenceThermalConductivity(double rho0, double T0);
593
594 /**
595 * Computes the composition-dependent values of parameters that are needed for the
596 * Lucas viscosity model.
597 *
598 * The equations for the mixing rules defined on page 9.23 of @cite poling2001 for
599 * the Lucas method's composition dependent parameters. The primary mixing rules
600 * are defined below, and the reduced properties are just the properties divided
601 * by the pseudo-critical mixture properties defined below.
602 *
603 * @note Equation numbers are from @cite poling2001
604 *
605 * @f[
606 * T_{\text{c,m}} = \sum_i X_i T_{\text{c,i}}
607 *
608 * \quad \text{( Equation 9-5.18)}
609 * @f]
610 *
611 * Where @f$ T_{\text{c,i}} @f$ is the critical temperature of species i,
612 * and @f$ X_i @f$ is the mole fraction of species i.
613 *
614 * @f[
615 * P_{\text{c,m}} = R T_{\text{c,m}}
616 * \frac{\sum_i X_i Z_{\text{c,i}}}{\sum_i X_i V_{\text{c,i}}}
617 *
618 * \quad \text{(Equation 9-5.19)}
619 * @f]
620 *
621 * Where @f$ Z_{\text{c,i}} @f$ is the critical compressibility of species i, and
622 * @f$ V_{\text{c,i}} @f$ is the critical volume of species i.
623 *
624 * @f[
625 * M_m = \sum X_i M_i
626 *
627 * \quad \text{(Equation 9-5.20)}
628 * @f]
629 *
630 * Where @f$ M_i @f$ is the molecular weight of species i.
631 *
632 * @f[
633 * F_{P,m}^{\text{o}} = \sum X_i F_{P,i}^{\text{o}}
634 *
635 * \quad \text{(Equation 9-5.21)}
636 * @f]
637 *
638 * Where @f$ F_{P,i}^{\text{o}} @f$ is the low-pressure polarity correction
639 * factor of species i from equation 9-4.18.
640 *
641 * @f[
642 * F_{Q,m}^{\text{o}} = \left ( \sum X_i F_{Q,i}^{\text{o}} \right ) A
643 *
644 * \quad \text{(Equation 9-5.22)}
645 * @f]
646 *
647 * Where @f$ F_{Q,i}^{\text{o}} @f$ is the low-pressure quantum correction factor
648 * of species i from equation 9-4.19, and A is defined below.
649 *
650 * @f[
651 * A = 1 - 0.01 \left ( \frac{M_H}{M_L} \right )^{0.87}
652 *
653 * \quad \text{(Equation 9-5.23)}
654 * @f]
655 *
656 * For @f$ \frac{M_H}{M_L} > 9 @f$ and @f$ 0.05 < X_H < 0.7 @f$, otherwise A = 1.
657 * In the above equation, $M_H$ and $M_L$ are the molecular weights of the
658 * heaviest and lightest components in the mixture, and @f$ X_H @f$ is the mole
659 * fraction of the heaviest component.
660 *
661 * While it isn't returned as a parameter, the species-specific reduced dipole
662 * moment (@f$ \mu_r @f$) is used to compute the mixture polarity correction factor.
663 * It is defined as:
664 *
665 * @f[
666 * \mu_r = 52.46 \frac{\mu^2 P_{\text{c,i}}}{T_{\text{c,i}}}
667 *
668 * \quad \text{(Equation 9-4.17)}
669 * @f]
670 */
672
673 /**
674 * Returns the non-dimensional low-pressure mixture viscosity in using the Lucas
675 * method.
676 *
677 * @f[
678 * \eta \xi = F_P^{\text{o}} F_Q^{\text{o}} [0.807 T_r^{0.618}
679 * - 0.357 e^{-0.449 T_r}
680 * + 0.340e^{-4.058 T_r} + 0.018]
681 *
682 * \quad \text{(Equation 9-4.16)}
683 * @f]
684 *
685 * This function is structured such that it can be used for pure species or
686 * mixtures, with the only difference being the values that are passed to the
687 * function (pure values versus mixture values).
688 *
689 * For the definition of the mixture rules, see computeMixtureParameters() .
690 *
691 * @param Tr Reduced temperature [unitless]
692 * @param FP Polarity correction factor [unitless]
693 * @param FQ Quantum correction factor [unitless]
694 */
695 double lowPressureNondimensionalViscosity(double Tr, double FP, double FQ);
696
697 /**
698 * Returns the non-dimensional high-pressure mixture viscosity in using the Lucas
699 * method.
700 *
701 * @f[
702 * \eta \xi = Z_2 F_P F_Q
703 *
704 * \quad \text{(Equation 9-6.12)}
705 * @f]
706 *
707 * This returns the value of η*ξ (by multiplying both sides of 9-6.12 by ξ and
708 * returning the right-side of the resulting equation).
709 *
710 * This function is structured such that it can be used for pure species or
711 * mixtures, with the only difference being the values that are passed to the
712 * function (pure values versus mixture values).
713 *
714 * For the definition of the mixture rules, see computeMixtureParameters() .
715 *
716 * @param Tr Reduced temperature [unitless]
717 * @param Pr Reduced pressure [unitless]
718 * @param FP_low Low-pressure polarity correction factor [unitless]
719 * @param FQ_low Low-pressure quantum correction factor [unitless]
720 * @param P_vap Vapor pressure [Pa]
721 * @param P_crit Critical pressure [Pa]
722 */
723 double highPressureNondimensionalViscosity(double Tr, double Pr, double FP_low,
724 double FQ_low, double P_vap,
725 double P_crit);
726
727 /**
728 * Calculates quantum correction term of the Lucas method for a species based
729 * on the reduced temperature(Tr) and molecular weight(MW), used in viscosity
730 * calculation.
731 *
732 * @f[
733 * F_{Q}^{\text{o}} = 1.22 Q^{0.15} {1 + 0.00385[ (T_r - 12)^2 ]^{\frac{1}{MW}}
734 * \text{sign} (T_r - 12 )}
735 *
736 * \quad \text{(Equation 9-4.19)}
737 * @f]
738 *
739 * @param Q Species-specific constant
740 * @param Tr Reduced temperature [unitless]
741 * @param MW Molecular weight [kg/kmol]
742 */
743 double quantumCorrectionFactor(double Q, double Tr, double MW);
744
745 /**
746 * Returns the polarity correction term for a species based on reduced temperature,
747 * reduced dipole moment, and critical compressibility. Used in the calculation of
748 * viscosity.
749 *
750 * Calculates polarity correction term of the Lucas method for a species based
751 * on the reduced temperature(Tr) and molecular weight(MW). Equation 9.4.18.
752 *
753 * @f[
754 * \begin{equation}
755 * F_P^0 =
756 * \begin{cases}
757 * 1 & 0 \leq \mu_r < 0.022 \\
758 * 1 + 30.55(0.292 - Z_c)^{1.72} & 0.022 \leq \mu_r < 0.075 \\
759 * 1 + 30.55(0.292 - Z_c)^{1.72} \times 0.96 + 0.1(T_r - 0.7) & 0.075 \leq \mu_r
760 * \end{cases}
761 * \end{equation}
762 * @f]
763 *
764 * @note The original description in Poling(2001) neglects to mention what happens
765 * when the quantity raised to the 1.72 power goes negative. That is an undefined
766 * operation that generates real and imaginary numbers. For now, only positive
767 * values are allowed.
768 *
769 * @param mu_r Species Reduced dipole moment
770 * @param Tr Reduced temperature
771 * @param Z_c Species Critical compressibility
772 */
773 double polarityCorrectionFactor(double mu_r, double Tr, double Z_c);
774
775
776private:
777
778 /**
779 * @name Reference fluid properties
780 * These are the properties of the reference fluid, which is methane in this case.
781 * These are used by the thermalConductivity() method.
782 * @{
783 */
784 const double m_ref_MW = 16.04; //!< Molecular weight [kg/kmol]
785 const double m_ref_Tc = 190.4; //!< Critical temperature [K]
786 const double m_ref_Vc = 0.0986; //!< Critical volume [m^3/kmol]
787 const double m_ref_Zc = 0.288; //!< Critical compressibility [unitless]
788 const double m_ref_rhoc = 0.1628; //!< Critical density [g/cm^3]
789 const double m_ref_acentric_factor = 0.011; //!< Acentric factor [unitless]
790 /** @} */
791
792 /**
793 * @name Lucas method viscosity parameters
794 * These are the parameters that are needed to calculate the viscosity using the
795 * Lucas method.
796 * @{
797 */
798 double m_FQ_mix_o; //!< Quantum correction factor
799 double m_FP_mix_o; //!< Polarity correction factor
800 double m_Tr_mix; //!< Reduced temperature
801 double m_Pr_mix; //!< Reduced pressure
802 double m_Pc_mix; //!< Critical pressure
803 double m_Tc_mix; //!< Critical temperature
804 double m_MW_mix; //!< Molecular weight
805 double m_P_vap_mix; //!< Vapor pressure
806 /** @} */
807};
808
809/**
810 * Transport properties for high pressure gas mixtures using the Chung method for
811 * viscosity and thermal conductivity.
812 *
813 * The implementation employs a method of corresponding states, using the Takahashi
814 * @cite takahashi1975 approach for binary diffusion coefficients (using mixture
815 * averaging rules for the mixture properties), and the Chung method for the viscosity
816 * and thermal conductivity of a high-pressure gas mixture. All methods are described
817 * in Poling et al. @cite poling2001 (viscosity in Ch. 9, thermal conductivity in
818 * Ch. 10, and diffusion coefficients in Ch. 11).
819 *
820 * @note All equations that are cited in this implementation are from the 5th edition
821 * of the book "The Properties of Gases and Liquids" by Poling, Prausnitz, and
822 * O'Connell.
823 *
824 * @ingroup tranprops
825 */
827{
828protected:
829 //! default constructor
831
832public:
833 string transportModel() const override {
834 return "high-pressure-Chung";
835 }
836
837 void init(ThermoPhase* thermo, int mode=0) override;
838
839 /**
840 * Returns the high-pressure mixture viscosity [Pa·s] using the Chung method.
841 *
842 * Based on the high-pressure gas mixture viscosity model of Chung described in
843 * chapter 9-7 of Poling. This method uses the pure species high-pressure viscosity
844 * relation of Chung with the difference being that mixture values of the model
845 * are computed using a set of mixing rules given by Chung
846 * (see computeMixtureParameters() ). The mixing rules are defined in section
847 * 9-5 of @cite poling2001.
848 *
849 * Because this method is using the high-pressure viscosity model with mixture
850 * parameters, see highPressureViscosity() for details on the model.
851 */
852 double viscosity() override;
853
854 /**
855 * Calculates the high-pressure mixture thermal conductivity using the Chung method.
856 *
857 * This method obtains appropriate mixture values of the parameters needed for the
858 * Chung model and then calls the highPressureThermalConductivity() method to
859 * obtain the mixture thermal conductivity.
860 *
861 * The mixture values of the pseudo-critical temperature and other model parameters
862 * are calculated using the Chung mixing rules defined on page 9.25
863 * (see computeMixtureParameters() ).
864 *
865 * The mixture value of the specific heat is computed using equation 10-6.6, which
866 * is the mole fraction weighted sum of the pure species specific heats. This value
867 * is not directly computed by the computeMixtureParameters() method.
868 *
869 * @f[
870 * C_{v,m} = \sum_i X_i C_{v,i}
871 *
872 * \quad \text{(Equation 10-6.6)}
873 * @f]
874 *
875 * Where @f$ C_{v,i} @f$ is the specific heat of species i, and @f$ X_i @f$ is the
876 * mole fraction of species i.
877 */
878 double thermalConductivity() override;
879
880 friend class TransportFactory;
881
882protected:
883
884 //! Computes and stores pure-fluid specific properties each species.
886
887 /**
888 * Computes the composition-dependent values of the parameters needed for
889 * the Chung viscosity model.
890 *
891 * The equations for the mixing rules defined on page 9.25 for the Chung method's
892 * composition dependent parameters. The primary mixing rules are defined below.
893 *
894 * @f[
895 * \sigma_m^3 = \sum_{i} \sum_{j} X_i X_j \sigma_{ij}^3
896 *
897 * \quad \text{(Equation 9-5.25)}
898 * @f]
899 *
900 * Where @f$ \sigma_{ij} @f$ is the molecular diameter
901 *
902 * @f[
903 * T_m^* = \frac{T}{\left( \frac{\epsilon}{k} \right )_m}
904 *
905 * \quad \text{(Equation 9-5.26)}
906 * @f]
907 *
908 * Where @f$ k @f$ is the Boltzmann constant and @f$ \epsilon @f$ is the minimum
909 * of the pair-potential energy. In these equations, we do not need to worry about
910 * what the values of @f$ \epsilon @f$ and @f$ k @f$ are.
911 *
912 * @f[
913 * \left( \frac{\epsilon}{k} \right)_m = \frac{\sum_{i} \sum_{j} X_i X_j
914 * \left( \frac{\epsilon_{ij}}{k} \right)
915 * \sigma_{ij}^3}{\sigma_m^3}
916 *
917 * \quad \text{(Equation 9-5.27)}
918 * @f]
919 *
920 * @f[
921 * MW_m = \left[ \frac{\sum_{i} \sum_{j} X_i X_j \left( \frac{\epsilon_{ij}}{k} \right)
922 * \sigma_{ij}^2 MW_{ij}^{\frac{1}{2}}}{\left( \frac{\epsilon}{k} \right)_m
923 * \sigma_m^2} \right]^2
924 *
925 * \quad \text{(Equation 9-5.28)}
926 * @f]
927 *
928 * Where MW is the molecular weight.
929 *
930 * @f[
931 * \omega_m = \frac{\sum_{i} \sum_{j} X_i X_j \omega_{ij} \sigma{ij}^3}{\sigma_m^3}
932 *
933 * \quad \text{(Equation 9-5.29)}
934 * @f]
935 *
936 * Where @f$ \omega @f$ is the acentric factor.
937 *
938 * @f[
939 * \mu_m^4 = \sigma_m^3 \sum_{i} \sum_{j} \left( \frac{X_i X_j \mu_i^2 \mu_j^2}
940 * {\sigma_{ij}^3} \right)
941 *
942 * \quad \text{(Equation 9-5.30)}
943 * @f]
944 *
945 * Where @f$ \mu @f$ is the dipole moment.
946 *
947 * @f[
948 * \kappa_m = \sum_{i} \sum_{j} X_i X_j \kappa_{ij}
949 *
950 * \quad \text{(Equation 9-5.31)}
951 * @f]
952 *
953 * Where @f$ \kappa @f$ is the association factor, which is used for highly polar
954 * molecules. In this work, the value is assumed to be zero for all species.
955 *
956 * The combining rules for species-species values (subscripted with ij in the above
957 * equations) are defined below.
958 *
959 * @f[
960 * \sigma_{i} = 0.809 V_{c,i}^{1/3}
961 *
962 * \quad \text{(Equation 9-5.32)}
963 * @f]
964 *
965 * Where @f$ V_{c,i} @f$ is the critical volume of species i.
966 *
967 * @f[
968 * \sigma_{ij} = \xi_{ij} \left( \sigma_{i} \sigma_{j} \right)^{1/2}
969 *
970 * \quad \text{(Equation 9-5.33)}
971 * @f]
972 *
973 * Where @f$ \xi @f$ is a binary interaction parameter.
974 *
975 * @f[
976 * \left( \frac{\epsilon_i}{k} \right) = \frac{T_{c,i}}{1.2593}
977 *
978 * \quad \text{(Equation 9-5.34)}
979 * @f]
980 *
981 * @f[
982 * \frac{\epsilon_{ij}}{k} = \zeta_{ij} \left( \right)^{\frac{1}{2}}
983 *
984 * \quad \text{(Equation 9-5.35)}
985 * @f]
986 *
987 * Where @f$ \zeta @f$ is a binary interaction parameter.
988 *
989 * @f[
990 * \omega_{ij} = \frac{\omega_i + \omega_j}{2}
991 *
992 * \quad \text{(Equation 9-5.37)}
993 * @f]
994 *
995 * Where @f$ \omega @f$ is the acentric factor.
996 *
997 * @f[
998 * \kappa_{ij} = \left( \kappa_i \kappa_j \right)^{1/2}
999 *
1000 * \quad \text{(Equation 9-5.39)}
1001 * @f]
1002 *
1003 * Where @f$ \kappa @f$ is the association factor.
1004 *
1005 * @f[
1006 * MW_{ij} = \frac{2 MW_i MW_j}{MW_i + MW_j}
1007 *
1008 * \quad \text{(Equation 9-5.40)}
1009 * @f]
1010 *
1011 * @f$ \xi @f$ and @f$ \zeta @f$ are the binary interaction parameters, and are
1012 * assumed to be unity in this implementation, in keeping with the Chung method.
1013 *
1014 * The Chung viscosity correction factor is defined as:
1015 *
1016 * @f[
1017 * F_{c,m} = 1 - 0.275 \omega_m + 0.059035 \mu_{r,m}^4 + \kappa_m
1018 *
1019 * \quad \text{(Equation 9-5.41)}
1020 * @f]
1021 *
1022 * Where @f$ \omega_m @f$ is the mixture acentric factor, @f$ \mu_{r,m} @f$ is the
1023 * mixture reduced dipole moment, and @f$ \kappa_m @f$ is the mixture association
1024 * factor.
1025 *
1026 * The mixture reduced dipole moment is computed using:
1027 *
1028 * @f[
1029 * \mu_{r,m} = \frac{131.3 \mu_m}{( V_{c,m} T_{c,m})^{\frac{1}{2}}}
1030 *
1031 * \quad \text{(Equation 9-5.42)}
1032 * @f]
1033 *
1034 * Where @f$ V_{c,m} @f$ and @f$ T_{c,m} @f$ are computed using the following
1035 * equations.
1036 *
1037 * @f[
1038 * V_{c,m} = \left( \frac{\sigma_m}{0.809} \right)^3
1039 *
1040 * \quad \text{(Equation 9-5.43)}
1041 * @f]
1042 *
1043 * @f[
1044 * T_{c,m} = 1.2593 \left( \frac{\epsilon}{k} \right)_m
1045 *
1046 * \quad \text{(Equation 9-5.44)}
1047 * @f]
1048 *
1049 * In the equations, @f$ T_c @f$ must be in units of K, @f$ V_c @f$ must be in
1050 * units of cm^3/mol, and @f$ \mu @f$ must be in units of Debye.
1051 */
1053
1054 /**
1055 * Returns the low-pressure mixture viscosity [Pa·s] using the Chung method.
1056 *
1057 * @f[
1058 * \eta = 26.69 F_c \frac{(MW*T)^(1/2)}{\sigma^2 \Omega}
1059 *
1060 * \quad \text{(Equation 9-4.10)}
1061 * @f]
1062 *
1063 * T must be in units of K, MW must be units of kg/kmol, and @f$ \sigma @f$ must
1064 * be in units of Angstroms. The viscosity is computed in micropoise, but the
1065 * return value is in standard SI units (Pa*s).
1066 *
1067 * This function is structured such that it can be used for pure species or
1068 * mixtures, with the only difference being the values that are passed to the
1069 * function (pure values versus mixture values).
1070 *
1071 * @param T Temperature [K]
1072 * @param T_star Reduced temperature [unitless]
1073 * @param MW Molecular weight [kg/kmol]
1074 * @param acentric_factor Acentric factor [unitless]
1075 * @param mu_r Dipole moment [Debye]
1076 * @param sigma Lennard-Jones collision diameter [Angstroms]
1077 * @param kappa Polar correction factor [unitless]
1078 */
1079 double lowPressureViscosity(double T, double T_star, double MW,
1080 double acentric_factor, double mu_r,
1081 double sigma, double kappa);
1082
1083 /**
1084 * Returns the high-pressure mixture viscosity [micropoise] using the Chung
1085 * method.
1086 *
1087 * @f[
1088 * \eta = \eta^* \frac{36.344 (M*T_c)^(1/2)}{V^{\frac{2}{3}}}
1089 *
1090 * \quad \text{(Equation 9-6.18)}
1091 * @f]
1092 *
1093 * where:
1094 *
1095 * @f[
1096 * \begin{align*}
1097 * \eta &= \text{viscosity, } \mu P \\
1098 * M &= \text{molecular weight, kg/kmol} \\
1099 * T_c &= \text{critical temperature, K} \\
1100 * V_c &= \text{critical molar volume, cm}^3 / \text{mol} \\
1101 * \end{align*}
1102 * @f]
1103 *
1104 * and,
1105 *
1106 * @f[
1107 * \eta^* = \frac{(T^*)^{\frac{1}{2}}}{\Omega_v} {F_c[(G_2)^{-1} + E_6 y]}
1108 * + \eta^{**}
1109 *
1110 * \quad \text{(Equation 9-6.19)}
1111 * @f]
1112 *
1113 * The values of @f$ T^* @f$ and @f$ F_c @f$ are defined as follows.
1114 *
1115 * @f[
1116 * T^* = 1.2593 T_r
1117 *
1118 * \quad \text{(Equation 9-4.9)}
1119 * @f]
1120 *
1121 * @f[
1122 * F_c = 1 - 0.275 \omega + 0.059035 \mu_r^4 + \kappa
1123 *
1124 * \quad \text{(Equation 9-4.11)}
1125 * @f]
1126 *
1127 * The value of @f$ \Omega_v @f$ is the viscosity collision integral evaluated at
1128 * the non-dimensional reduced temperature @f$ T^* @f$. For details on the
1129 * collision integral see neufeldCollisionIntegral() .
1130 *
1131 * This function is structured such that it can be used for pure species or
1132 * mixtures, with the only difference being the values that are passed to the
1133 * function (pure values versus mixture values).
1134 *
1135 * @param T_star Reduced temperature [unitless]
1136 * @param MW Molecular weight [kg/kmol]
1137 * @param rho Density [mol/cm³]
1138 * @param Vc Critical volume [cm³/mol]
1139 * @param Tc Critical temperature [K]
1140 * @param acentric_factor Acentric factor [unitless]
1141 * @param mu_r Dipole moment [Debye]
1142 * @param kappa Polar correction factor [unitless]
1143 */
1144 double highPressureViscosity(double T_star, double MW, double rho, double Vc,
1145 double Tc, double acentric_factor, double mu_r, double kappa);
1146
1147 /**
1148 * Computes the high-pressure thermal conductivity [W/m/K] using the Chung method.
1149 *
1150 * The Chung method for computing high-pressure thermal conductivity is described
1151 * on page 10.23 of @cite poling2001 .
1152 *
1153 * @f[
1154 * \lambda = \frac{31.2 \eta^0 \Psi}{M'} \left( G_2^{-1} + B_6 y \right)
1155 * + q B_7 y^2 T_r^{1/2} G_2
1156 *
1157 * \quad \text{(Equation 10-5.5)}
1158 * @f]
1159 *
1160 * where:
1161 *
1162 * @f[
1163 * \begin{align*}
1164 * \lambda &= \text{thermal conductivity, W/(m·K)} \\
1165 * \eta^0 &= \text{low-pressure gas viscosity, N·s/m}^2 \\
1166 * M' &= \text{molecular weight, kg/mol} \\
1167 * \Psi &= f(C_v, \omega, T_r) \text{ [as defined under Eq. (10-3.14)]} \\
1168 * q &= 3.586 \times 10^{-3} \left( \frac{T_c}{M'} \right)^{1/2} V_c^{2/3} \\
1169 * T &= \text{temperature, K} \\
1170 * T_c &= \text{critical temperature, K} \\
1171 * T_r &= \text{reduced temperature, } \frac{T}{T_c} \\
1172 * V_c &= \text{critical molar volume, cm}^3/\text{mol} \\
1173 * \end{align*}
1174 * @f]
1175 *
1176 * where,
1177 *
1178 * @f[
1179 * y = \frac{V_c}{6V}
1180 *
1181 * \quad \text{(Equation 10-5.6)}
1182 * @f]
1183 *
1184 * V is the molar volume of the fluid in cm^3/mol.
1185 *
1186 * @f[
1187 * G_1 = \frac{1 - 0.5y}{(1-y)^3}
1188 *
1189 * \quad \text{(Equation 10-5.7)}
1190 * @f]
1191 *
1192 * @f[
1193 * G_2 = \frac{(B_1 / y)[1 - \text{exp}(-B_4 y)]
1194 * + B_2 G_1 \text{exp}(B_5 y) + B_3 G_1}{B_1 B_4 + B_2 + B_3}
1195 *
1196 * \quad \text{(Equation 10-5.8)}
1197 * @f]
1198 *
1199 * The coefficients @f$ B_1 @f$, through @f$ B_7 @f$ are functions of the
1200 * acentric factor, reduced dipole moment and the association factor.
1201 *
1202 * @f[
1203 * B_i = a_i + b_i \omega + c_i \mu_r^4 + d_i \kappa
1204 *
1205 * \quad \text{(Equation 10-5.9)}
1206 * @f]
1207 *
1208 * The constants in the above equation are from table 10-3 on page 10.23 of
1209 * @cite poling2001.
1210 *
1211 * The definition of the @f$ \Psi @f$ function is given by:
1212 *
1213 * @f[
1214 * \Psi = 1 + \alpha {\frac{0.215 + 0.28288\alpha - 1.061\beta + 0.26665Z}{0.6366
1215 * + \beta Z + 1.061 \alpha \beta}}
1216 * @f]
1217 *
1218 * with,
1219 *
1220 * @f[
1221 * \alpha = \frac{C_v}{R} - \frac{3}{2}
1222 * @f]
1223 *
1224 * @f[
1225 * \beta = 0.7862 - 0.7109 \omega + 1.3168 \omega^2
1226 * @f]
1227 *
1228 * @f[
1229 * Z = 2.0 + 10.5 T_r^2
1230 * @f]
1231 *
1232 * These functions are from page 10.12 of @cite poling2001 .
1233 *
1234 * This method utilizes the low-pressure Chung viscosity as that is a required
1235 * parameter in the model, and thus calls the low pressure viscosity
1236 * implementation. This is why it requires parameters typically associated with
1237 * the viscosity calculation.
1238 *
1239 * This function is structured such that it can be used for pure species or
1240 * mixtures, with the only difference being the values that are passed to the
1241 * function (pure values versus mixture values).
1242 *
1243 * For mixtures, the mixture values of the input variables are computed using the
1244 * mixing rules of Chung, see computeMixtureParameters() .
1245 *
1246 * @param T Temperature [K]
1247 * @param T_star Reduced temperature [unitless]
1248 * @param MW Molecular weight [kg/kmol]
1249 * @param rho Density [mol/cm³]
1250 * @param Cv Specific heat [J/kg/K]
1251 * @param Vc Critical volume [cm³/mol]
1252 * @param Tc Critical temperature [K]
1253 * @param sigma Lennard-Jones collision diameter [Angstroms]
1254 * @param acentric_factor Acentric factor [unitless]
1255 * @param mu_r Dipole moment [Debye]
1256 * @param kappa Polar correction factor [unitless]
1257 * @return High pressure thermal conductivity [W/m/K]
1258 */
1259 double highPressureThermalConductivity(double T, double T_star, double MW,
1260 double rho, double Cv, double Vc, double Tc, double sigma,
1261 double acentric_factor, double mu_r, double kappa);
1262
1263private:
1264
1265 /**
1266 * @name Pure fluid properties
1267 * These are the pure fluid properties of each species that are needed to compute
1268 * the mixture properties for the Chung viscosity and thermal conductivity models.
1269 * @{
1270 */
1271 vector<double> m_sigma_i; //!< Effective molecular diameter [Angstroms]
1272 vector<double> m_epsilon_over_k_i; //!< Characteristic temperature [K]
1273 vector<double> m_MW_i; //!< Molecular weight [kg/kmol]
1274 vector<double> m_acentric_factor_i; //!< Acentric factor [unitless]
1275 vector<double> m_kappa_i; //!< Association factor [unitless]
1276 /** @} */
1277
1278
1279 /**
1280 * @name Chung mixture parameters
1281 * These are the parameters that are needed to calculate the viscosity and thermal
1282 * conductivity using the Chung method.
1283 * @{
1284 */
1285 double m_Vc_mix = 0; //!< Mixture critical volume [m³/kmol]
1286 double m_Tc_mix = 0; //!< Mixture critical temperature [K]
1287
1288 double m_sigma_mix = 0; //!< Effective mixture molecular diameter [Angstroms]
1289 double m_epsilon_over_k_mix = 0; //!< Mixture characteristic temperature [K]
1290 double m_MW_mix = 0; //!< Effective mixture molecular weight [kg/kmol]
1291
1292 // These are used to compute the Fc factor in the Chung viscosity model
1293 double m_mu_mix = 0; //!< Mixture dipole moment [Debye]
1294 double m_mu_r_mix = 0; //!< Mixture reduced dipole moment [unitless]
1295 double m_acentric_factor_mix = 0; //!< Mixture acentric factor [unitless]
1296 double m_kappa_mix = 0; //!< Mixture association factor [unitless]
1297
1298 /** @} */
1299
1300};
1301
1302}
1303#endif
Headers for the MixTransport object, which models transport properties in ideal gas solutions using a...
Transport properties for high pressure gas mixtures using the Chung method for viscosity and thermal ...
double m_kappa_mix
Mixture association factor [unitless].
vector< double > m_MW_i
Molecular weight [kg/kmol].
void initializePureFluidProperties()
Computes and stores pure-fluid specific properties each species.
double thermalConductivity() override
Calculates the high-pressure mixture thermal conductivity using the Chung method.
double m_MW_mix
Effective mixture molecular weight [kg/kmol].
double m_mu_r_mix
Mixture reduced dipole moment [unitless].
double m_Tc_mix
Mixture critical temperature [K].
double m_epsilon_over_k_mix
Mixture characteristic temperature [K].
vector< double > m_epsilon_over_k_i
Characteristic temperature [K].
vector< double > m_acentric_factor_i
Acentric factor [unitless].
double highPressureThermalConductivity(double T, double T_star, double MW, double rho, double Cv, double Vc, double Tc, double sigma, double acentric_factor, double mu_r, double kappa)
Computes the high-pressure thermal conductivity [W/m/K] using the Chung method.
void init(ThermoPhase *thermo, int mode=0) override
Initialize a transport manager.
double lowPressureViscosity(double T, double T_star, double MW, double acentric_factor, double mu_r, double sigma, double kappa)
Returns the low-pressure mixture viscosity [Pa·s] using the Chung method.
double m_Vc_mix
Mixture critical volume [m³/kmol].
vector< double > m_sigma_i
Effective molecular diameter [Angstroms].
ChungHighPressureGasTransport()=default
default constructor
double m_mu_mix
Mixture dipole moment [Debye].
double viscosity() override
Returns the high-pressure mixture viscosity [Pa·s] using the Chung method.
double m_sigma_mix
Effective mixture molecular diameter [Angstroms].
double m_acentric_factor_mix
Mixture acentric factor [unitless].
vector< double > m_kappa_i
Association factor [unitless].
double highPressureViscosity(double T_star, double MW, double rho, double Vc, double Tc, double acentric_factor, double mu_r, double kappa)
Returns the high-pressure mixture viscosity [micropoise] using the Chung method.
string transportModel() const override
Identifies the model represented by this Transport object.
void computeMixtureParameters()
Computes the composition-dependent values of the parameters needed for the Chung viscosity model.
A class for full (non-sparse) matrices with Fortran-compatible data storage, which adds matrix operat...
Definition DenseMatrix.h:55
Base class for high pressure gas transport models.
HighPressureGasTransportBase()=default
default constructor
void getBinaryDiffCoeffs(const size_t ld, double *const d) override
Computes the matrix of binary diffusion coefficients [m²/s] using the Takahashi correction factor.
vector< double > m_Pcrit
Critical pressure [Pa] of each species.
DenseMatrix m_P_corr_ij
Matrix of Takahashi binary diffusion coefficient corrections.
virtual void updateCorrectionFactors()
Updates the matrix of species-pair Takahashi correction factors for use in computing the binary diffu...
double Pcrit_i(size_t i)
Returns the stored value of the critical pressure for species 'i'.
void getMixDiffCoeffs(double *const d) override
Returns the mixture-averaged diffusion coefficients [m²/s].
double Vcrit_i(size_t i)
Returns the stored value of the critical volume for species 'i'.
void getMixDiffCoeffsMole(double *const d) override
Returns the mixture-averaged diffusion coefficients [m²/s].
vector< double > m_Zcrit
Critical compressibility [unitless] of each species.
double Tcrit_i(size_t i)
Returns the stored value of the critical temperature for species 'i'.
double Zcrit_i(size_t i)
Returns the stored value of the critical compressibility for species 'i'.
vector< double > m_Vcrit
Critical volume [m³/kmol] of each species.
void getTransportData() override
Obtain required parameters from the 'critical-parameters' species input section, and checks the criti...
vector< double > m_Tcrit
Critical temperature [K] of each species.
void getMixDiffCoeffsMass(double *const d) override
Returns the mixture-averaged diffusion coefficients [m²/s].
void initializeCriticalProperties()
Computes and stores the estimate of the critical properties for each species.
The implementation employs a method of corresponding states, using the Takahashi approach for binary...
double polarityCorrectionFactor(double mu_r, double Tr, double Z_c)
Returns the polarity correction term for a species based on reduced temperature, reduced dipole momen...
const double m_ref_Vc
Critical volume [m^3/kmol].
double elyHanleyDilutePureSpeciesViscosity(double V, double Tc, double Vc, double Zc, double acentric_factor, double mw)
Get the viscosity [Pa·s] of a pure species using the method of Ely and Hanley.
double elyHanleyDiluteReferenceViscosity(double T0)
Returns the viscosity [Pa·s] for the reference fluid (methane) for low pressures.
double thermalConductivity() override
Returns the mixture high-pressure thermal conductivity [W/m/K] using a method of corresponding states...
double highPressureNondimensionalViscosity(double Tr, double Pr, double FP_low, double FQ_low, double P_vap, double P_crit)
Returns the non-dimensional high-pressure mixture viscosity in using the Lucas method.
double quantumCorrectionFactor(double Q, double Tr, double MW)
Calculates quantum correction term of the Lucas method for a species based on the reduced temperature...
const double m_ref_MW
Molecular weight [kg/kmol].
double lowPressureNondimensionalViscosity(double Tr, double FP, double FQ)
Returns the non-dimensional low-pressure mixture viscosity in using the Lucas method.
double elyHanleyReferenceThermalConductivity(double rho0, double T0)
Returns the thermal conductivity [W/m/K] of the reference fluid of methane.
double m_FP_mix_o
Polarity correction factor.
void init(ThermoPhase *thermo, int mode=0) override
Initialize a transport manager.
const double m_ref_Zc
Critical compressibility [unitless].
double viscosity() override
Returns the mixture high-pressure viscosity [Pa·s] using the Lucas method.
string transportModel() const override
Identifies the model represented by this Transport object.
void computeMixtureParameters()
Computes the composition-dependent values of parameters that are needed for the Lucas viscosity model...
double m_FQ_mix_o
Quantum correction factor.
const double m_ref_Tc
Critical temperature [K].
double phiShapeFactor(double Tr, double Vr, double Zc, double acentric_factor)
Returns the phi shape factor of Leach and Leland for a pure species.
const double m_ref_acentric_factor
Acentric factor [unitless].
HighPressureGasTransport()=default
default constructor
double thetaShapeFactor(double Tr, double Vr, double acentric_factor)
Returns the theta shape factor of Leach and Leland for a pure species.
const double m_ref_rhoc
Critical density [g/cm^3].
Class MixTransport implements mixture-averaged transport properties for ideal gas mixtures.
Base class for a phase with thermodynamic properties.
Factory class for creating new instances of classes derived from Transport.
ThermoPhase & thermo()
Phase object.
Definition Transport.h:103
Namespace for the Cantera kernel.
Definition AnyMap.cpp:595