The implementation employs a method of corresponding states, using the Takahashi [48] approach for binary diffusion coefficients (using mixture averaging rules for the mixture properties), the Lucas method for the viscosity, and a method from Ely and Hanley for the thermal conductivity. More...
#include <HighPressureGasTransport.h>
The implementation employs a method of corresponding states, using the Takahashi [48] approach for binary diffusion coefficients (using mixture averaging rules for the mixture properties), the Lucas method for the viscosity, and a method from Ely and Hanley for the thermal conductivity.
All methods are described in Poling et al. [38] (viscosity in Ch. 9, thermal conductivity in Ch. 10, and diffusion coefficients in Ch. 11).
Definition at line 151 of file HighPressureGasTransport.h.
Public Member Functions | |
string | transportModel () const override |
Identifies the model represented by this Transport object. | |
void | init (ThermoPhase *thermo, int mode=0) override |
Initialize a transport manager. | |
double | thermalConductivity () override |
Returns the mixture high-pressure thermal conductivity [W/m/K] using a method of corresponding states by Ely and Hanley. | |
double | viscosity () override |
Returns the mixture high-pressure viscosity [Pa·s] using the Lucas method. | |
![]() | |
void | getBinaryDiffCoeffs (const size_t ld, double *const d) override |
Computes the matrix of binary diffusion coefficients [m²/s] using the Takahashi correction factor. | |
void | getMixDiffCoeffs (double *const d) override |
Returns the mixture-averaged diffusion coefficients [m²/s]. | |
void | getMixDiffCoeffsMole (double *const d) override |
Returns the mixture-averaged diffusion coefficients [m²/s]. | |
void | getMixDiffCoeffsMass (double *const d) override |
Returns the mixture-averaged diffusion coefficients [m²/s]. | |
virtual void | updateCorrectionFactors () |
Updates the matrix of species-pair Takahashi correction factors for use in computing the binary diffusion coefficients, see takahashiCorrectionFactor() | |
![]() | |
MixTransport ()=default | |
Default constructor. | |
string | transportModel () const override |
Identifies the model represented by this Transport object. | |
void | getThermalDiffCoeffs (double *const dt) override |
Return the thermal diffusion coefficients [kg/m/s]. | |
double | thermalConductivity () override |
Returns the mixture thermal conductivity [W/m/K]. | |
void | getMobilities (double *const mobil) override |
Get the electrical mobilities [m²/V/s]. | |
void | update_T () override |
Update the internal parameters whenever the temperature has changed. | |
void | update_C () override |
Update the internal parameters whenever the concentrations have changed. | |
void | getSpeciesFluxes (size_t ndim, const double *const grad_T, size_t ldx, const double *const grad_X, size_t ldf, double *const fluxes) override |
Get the species diffusive mass fluxes [kg/m²/s] with respect to the mass averaged velocity, given the gradients in mole fraction and temperature. | |
void | init (ThermoPhase *thermo, int mode=0) override |
Initialize a transport manager. | |
![]() | |
double | viscosity () override |
Get the viscosity [Pa·s] of the mixture. | |
void | getSpeciesViscosities (double *const visc) override |
Get the pure-species viscosities [Pa·s]. | |
void | getBinaryDiffCoeffs (const size_t ld, double *const d) override |
Returns the matrix of binary diffusion coefficients [m²/s]. | |
void | getMixDiffCoeffs (double *const d) override |
Returns the Mixture-averaged diffusion coefficients [m²/s]. | |
void | getMixDiffCoeffsMole (double *const d) override |
Returns the mixture-averaged diffusion coefficients [m²/s]. | |
void | getMixDiffCoeffsMass (double *const d) override |
Returns the mixture-averaged diffusion coefficients [m²/s]. | |
void | getViscosityPolynomial (size_t i, double *coeffs) const override |
Return the polynomial fits to the viscosity of species i . | |
void | getConductivityPolynomial (size_t i, double *coeffs) const override |
Return the temperature fits of the heat conductivity of species i . | |
void | getBinDiffusivityPolynomial (size_t i, size_t j, double *coeffs) const override |
Return the polynomial fits to the binary diffusivity of species pair (i, j) | |
void | getCollisionIntegralPolynomial (size_t i, size_t j, double *astar_coeffs, double *bstar_coeffs, double *cstar_coeffs) const override |
Return the polynomial fits to the collision integral of species pair (i, j) | |
void | setViscosityPolynomial (size_t i, double *coeffs) override |
Modify the polynomial fits to the viscosity of species i | |
void | setConductivityPolynomial (size_t i, double *coeffs) override |
Modify the temperature fits of the heat conductivity of species i | |
void | setBinDiffusivityPolynomial (size_t i, size_t j, double *coeffs) override |
Modify the polynomial fits to the binary diffusivity of species pair (i, j) | |
void | setCollisionIntegralPolynomial (size_t i, size_t j, double *astar_coeffs, double *bstar_coeffs, double *cstar_coeffs, bool actualT) override |
Modify the polynomial fits to the collision integral of species pair (i, j) | |
void | init (ThermoPhase *thermo, int mode=0) override |
Initialize a transport manager. | |
bool | CKMode () const override |
Boolean indicating the form of the transport properties polynomial fits. | |
void | invalidateCache () override |
Invalidate any cached values which are normally updated only when a change in state is detected. | |
![]() | |
Transport ()=default | |
Constructor. | |
Transport (const Transport &)=delete | |
Transport & | operator= (const Transport &)=delete |
virtual string | transportModel () const |
Identifies the model represented by this Transport object. | |
ThermoPhase & | thermo () |
Phase object. | |
void | checkSpeciesIndex (size_t k) const |
Check that the specified species index is in range. | |
void | checkSpeciesArraySize (size_t kk) const |
Check that an array size is at least m_nsp. | |
virtual void | getSpeciesFluxes (size_t ndim, const double *const grad_T, size_t ldx, const double *const grad_X, size_t ldf, double *const fluxes) |
Get the species diffusive mass fluxes [kg/m²/s] with respect to the specified solution averaged velocity, given the mole fraction and temperature gradients. | |
virtual void | getMolarFluxes (const double *const state1, const double *const state2, const double delta, double *const cfluxes) |
Get the molar fluxes [kmol/m²/s], given the thermodynamic state at two nearby points. | |
virtual void | getMassFluxes (const double *state1, const double *state2, double delta, double *mfluxes) |
Get the mass fluxes [kg/m²/s], given the thermodynamic state at two nearby points. | |
virtual void | getThermalDiffCoeffs (double *const dt) |
Return a vector of thermal diffusion coefficients [kg/m/s]. | |
virtual void | getBinaryDiffCoeffs (const size_t ld, double *const d) |
Returns the matrix of binary diffusion coefficients [m²/s]. | |
virtual void | getMultiDiffCoeffs (const size_t ld, double *const d) |
Return the multicomponent diffusion coefficients [m²/s]. | |
virtual void | getMixDiffCoeffs (double *const d) |
Return a vector of mixture averaged diffusion coefficients [m²/s]. | |
virtual void | getMixDiffCoeffsMole (double *const d) |
Returns a vector of mixture averaged diffusion coefficients [m²/s]. | |
virtual void | getMixDiffCoeffsMass (double *const d) |
Returns a vector of mixture averaged diffusion coefficients [m²/s]. | |
virtual void | getViscosityPolynomial (size_t i, double *coeffs) const |
Return the polynomial fits to the viscosity of species i . | |
virtual void | getConductivityPolynomial (size_t i, double *coeffs) const |
Return the temperature fits of the heat conductivity of species i . | |
virtual void | getBinDiffusivityPolynomial (size_t i, size_t j, double *coeffs) const |
Return the polynomial fits to the binary diffusivity of species pair (i, j) | |
virtual void | getCollisionIntegralPolynomial (size_t i, size_t j, double *astar_coeffs, double *bstar_coeffs, double *cstar_coeffs) const |
Return the polynomial fits to the collision integral of species pair (i, j) | |
virtual void | setViscosityPolynomial (size_t i, double *coeffs) |
Modify the polynomial fits to the viscosity of species i | |
virtual void | setConductivityPolynomial (size_t i, double *coeffs) |
Modify the temperature fits of the heat conductivity of species i | |
virtual void | setBinDiffusivityPolynomial (size_t i, size_t j, double *coeffs) |
Modify the polynomial fits to the binary diffusivity of species pair (i, j) | |
virtual void | setCollisionIntegralPolynomial (size_t i, size_t j, double *astar_coeffs, double *bstar_coeffs, double *cstar_coeffs, bool flag) |
Modify the polynomial fits to the collision integral of species pair (i, j) | |
AnyMap | parameters () const |
Return the parameters for a phase definition which are needed to reconstruct an identical object using the newTransport() function. | |
AnyMap | fittingErrors () const |
Get error metrics about any functional fits calculated for pure species transport properties. | |
virtual void | invalidateCache () |
Invalidate any cached values which are normally updated only when a change in state is detected. | |
virtual double | bulkViscosity () |
The bulk viscosity [Pa·s]. | |
virtual double | electricalConductivity () |
Get the electrical conductivity [siemens/m]. | |
Protected Member Functions | |
HighPressureGasTransport ()=default | |
default constructor | |
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 | thetaShapeFactor (double Tr, double Vr, double acentric_factor) |
Returns the theta shape factor of Leach and Leland for a pure species. | |
double | phiShapeFactor (double Tr, double Vr, double Zc, double acentric_factor) |
Returns the phi shape factor of Leach and Leland for a pure species. | |
double | elyHanleyDiluteReferenceViscosity (double T0) |
Returns the viscosity [Pa·s] for the reference fluid (methane) for low pressures. | |
double | elyHanleyReferenceThermalConductivity (double rho0, double T0) |
Returns the thermal conductivity [W/m/K] of the reference fluid of methane. | |
void | computeMixtureParameters () |
Computes the composition-dependent values of parameters that are needed for the Lucas viscosity model. | |
double | lowPressureNondimensionalViscosity (double Tr, double FP, double FQ) |
Returns the non-dimensional low-pressure mixture viscosity in using the Lucas method. | |
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(Tr) and molecular weight(MW), used in viscosity calculation. | |
double | polarityCorrectionFactor (double mu_r, double Tr, double Z_c) |
Returns the polarity correction term for a species based on reduced temperature, reduced dipole moment, and critical compressibility. | |
![]() | |
HighPressureGasTransportBase ()=default | |
default constructor | |
void | getTransportData () override |
Obtain required parameters from the 'critical-parameters' species input section, and checks the critical-properties.yaml file if an acentric factor is not specified. | |
void | initializeCriticalProperties () |
Computes and stores the estimate of the critical properties for each species. | |
double | Tcrit_i (size_t i) |
Returns the stored value of the critical temperature for species 'i'. | |
double | Pcrit_i (size_t i) |
Returns the stored value of the critical pressure for species 'i'. | |
double | Vcrit_i (size_t i) |
Returns the stored value of the critical volume for species 'i'. | |
double | Zcrit_i (size_t i) |
Returns the stored value of the critical compressibility for species 'i'. | |
![]() | |
void | updateCond_T () |
Update the temperature dependent parts of the species thermal conductivities. | |
![]() | |
virtual void | update_T () |
virtual void | update_C ()=0 |
virtual void | updateViscosity_T () |
Update the temperature-dependent viscosity terms. | |
virtual void | updateSpeciesViscosities () |
Update the pure-species viscosities. | |
virtual void | updateDiff_T () |
Update the binary diffusion coefficients. | |
virtual void | setupCollisionParameters () |
Setup parameters for a new kinetic-theory-based transport manager for low-density gases. | |
void | setupCollisionIntegral () |
Setup range for polynomial fits to collision integrals of Monchick & Mason [31]. | |
void | makePolarCorrections (size_t i, size_t j, double &f_eps, double &f_sigma) |
Corrections for polar-nonpolar binary diffusion coefficients. | |
void | fitCollisionIntegrals (MMCollisionInt &integrals) |
Generate polynomial fits to collision integrals. | |
virtual void | fitProperties (MMCollisionInt &integrals) |
Generate polynomial fits to the viscosity \( \eta \) and conductivity \( \lambda \). | |
virtual void | fitDiffCoeffs (MMCollisionInt &integrals) |
Generate polynomial fits to the binary diffusion coefficients. | |
void | getBinDiffCorrection (double t, MMCollisionInt &integrals, size_t k, size_t j, double xk, double xj, double &fkj, double &fjk) |
Second-order correction to the binary diffusion coefficients. | |
Private Attributes | |
Reference fluid properties | |
These are the properties of the reference fluid, which is methane in this case. These are used by the thermalConductivity() method. | |
const double | m_ref_MW = 16.04 |
Molecular weight [kg/kmol]. | |
const double | m_ref_Tc = 190.4 |
Critical temperature [K]. | |
const double | m_ref_Vc = 0.0986 |
Critical volume [m^3/kmol]. | |
const double | m_ref_Zc = 0.288 |
Critical compressibility [unitless]. | |
const double | m_ref_rhoc = 0.1628 |
Critical density [g/cm^3]. | |
const double | m_ref_acentric_factor = 0.011 |
Acentric factor [unitless]. | |
Lucas method viscosity parameters | |
These are the parameters that are needed to calculate the viscosity using the Lucas method. | |
double | m_FQ_mix_o |
Quantum correction factor. | |
double | m_FP_mix_o |
Polarity correction factor. | |
double | m_Tr_mix |
Reduced temperature. | |
double | m_Pr_mix |
Reduced pressure. | |
double | m_Pc_mix |
Critical pressure. | |
double | m_Tc_mix |
Critical temperature. | |
double | m_MW_mix |
Molecular weight. | |
double | m_P_vap_mix |
Vapor pressure. | |
Friends | |
class | TransportFactory |
Additional Inherited Members | |
![]() | |
vector< double > | m_Tcrit |
Critical temperature [K] of each species. | |
vector< double > | m_Pcrit |
Critical pressure [Pa] of each species. | |
vector< double > | m_Vcrit |
Critical volume [m³/kmol] of each species. | |
vector< double > | m_Zcrit |
Critical compressibility [unitless] of each species. | |
DenseMatrix | m_P_corr_ij |
Matrix of Takahashi binary diffusion coefficient corrections. | |
![]() | |
vector< double > | m_cond |
vector of species thermal conductivities [W/m/K] | |
double | m_lambda = 0.0 |
Internal storage for the calculated mixture thermal conductivity [W/m/K]. | |
bool | m_spcond_ok = false |
Update boolean for the species thermal conductivities. | |
bool | m_condmix_ok = false |
Update boolean for the mixture rule for the mixture thermal conductivity. | |
![]() | |
vector< double > | m_molefracs |
Vector of species mole fractions. | |
double | m_viscmix = 0.0 |
Internal storage for the viscosity of the mixture [Pa·s]. | |
bool | m_visc_ok = false |
Update boolean for mixture rule for the mixture viscosity. | |
bool | m_viscwt_ok = false |
Update boolean for the weighting factors for the mixture viscosity. | |
bool | m_spvisc_ok = false |
Update boolean for the species viscosities. | |
bool | m_bindiff_ok = false |
Update boolean for the binary diffusivities at unit pressure. | |
int | m_mode = 0 |
Type of the polynomial fits to temperature. | |
DenseMatrix | m_phi |
Viscosity weighting function. size = m_nsp * m_nsp. | |
vector< double > | m_spwork |
work space length = m_nsp | |
vector< double > | m_visc |
vector of species viscosities [Pa·s]. | |
vector< vector< double > > | m_visccoeffs |
Polynomial fits to the viscosity of each species. | |
vector< double > | m_mw |
Local copy of the species molecular weights. | |
DenseMatrix | m_wratjk |
Holds square roots of molecular weight ratios. | |
DenseMatrix | m_wratkj1 |
Holds square roots of molecular weight ratios. | |
vector< double > | m_sqvisc |
vector of square root of species viscosities. | |
vector< double > | m_polytempvec |
Powers of the ln temperature, up to fourth order. | |
double | m_temp = -1.0 |
Current value of the temperature [K] at which the properties in this object are calculated. | |
double | m_kbt = 0.0 |
Current value of Boltzmann constant times the temperature [J]. | |
double | m_sqrt_t = 0.0 |
current value of temperature to 1/2 power | |
double | m_logt = 0.0 |
Current value of the log of the temperature. | |
double | m_t14 = 0.0 |
Current value of temperature to 1/4 power. | |
vector< vector< double > > | m_diffcoeffs |
Polynomial fits to the binary diffusivity of each species. | |
DenseMatrix | m_bdiff |
Matrix of binary diffusion coefficients at the reference pressure and the current temperature Size is m_nsp x m_nsp. | |
vector< vector< double > > | m_condcoeffs |
temperature fits of the heat conduction | |
vector< vector< int > > | m_poly |
Indices for the (i,j) interaction in collision integral fits. | |
vector< vector< double > > | m_omega22_poly |
Fit for omega22 collision integral. | |
vector< vector< int > > | m_star_poly_uses_actualT |
Flag to indicate for which (i,j) interaction pairs the actual temperature is used instead of the reduced temperature. | |
vector< vector< double > > | m_astar_poly |
Fit for astar collision integral. | |
vector< vector< double > > | m_bstar_poly |
Fit for bstar collision integral. | |
vector< vector< double > > | m_cstar_poly |
Fit for cstar collision integral. | |
vector< double > | m_zrot |
Rotational relaxation number for each species. | |
vector< double > | m_crot |
Dimensionless rotational heat capacity of each species. | |
vector< bool > | m_polar |
Vector of booleans indicating whether a species is a polar molecule. | |
vector< double > | m_alpha |
Polarizability [m³] of each species in the phase. | |
vector< double > | m_eps |
Lennard-Jones well-depth [J] of the species in the current phase. | |
vector< double > | m_sigma |
Lennard-Jones diameter [m] of the species in the current phase. | |
DenseMatrix | m_reducedMass |
This is the reduced mass [kg] of the interaction between species i and j. | |
DenseMatrix | m_diam |
hard-sphere diameter [m] for (i,j) collision | |
DenseMatrix | m_epsilon |
The effective well depth [J] for (i,j) collisions. | |
DenseMatrix | m_dipole |
The effective dipole moment [Coulomb·m] for (i,j) collisions. | |
DenseMatrix | m_delta |
Reduced dipole moment of the interaction between two species. | |
vector< double > | m_w_ac |
Pitzer acentric factor [dimensionless]. | |
vector< double > | m_disp |
Dispersion coefficient normalized by the square of the elementary charge [m⁵]. | |
vector< double > | m_quad_polar |
Quadrupole polarizability. | |
![]() | |
ThermoPhase * | m_thermo |
pointer to the object representing the phase | |
size_t | m_nsp = 0 |
Number of species in the phase. | |
AnyMap | m_fittingErrors |
Maximum errors associated with fitting pure species transport properties. | |
|
protecteddefault |
default constructor
|
inlineoverridevirtual |
Identifies the model represented by this Transport object.
Each derived class should override this method to return a meaningful identifier.
Reimplemented from Transport.
Definition at line 158 of file HighPressureGasTransport.h.
|
overridevirtual |
Initialize a transport manager.
This routine sets up a transport manager. It calculates the collision integrals and populates species-dependent data structures.
thermo | Pointer to the ThermoPhase object |
mode | Chemkin compatible mode or not. This alters the specification of the collision integrals. defaults to no. |
Reimplemented from GasTransport.
Definition at line 343 of file HighPressureGasTransport.cpp.
|
overridevirtual |
Returns the mixture high-pressure thermal conductivity [W/m/K] using a method of corresponding states by Ely and Hanley.
This method for computing the thermal conductivity is described in [9] . There are references to earlier work in [8] in the description of the model for thermal conductivity. The equations referenced in this description primarily are from [9] , and equations that come from [8] are marked as such. The equations referenced here are the same as the ones from the papers from Ely and Hanley.
The thermal conductivity is assumed to have two components: a translational/collisional part and an internal part that is related to the transfer of energy to internal degrees of freedom.
\[ \lambda_{mix}(\rho, T) = \lambda^{''}_{mix}(T) + \lambda^{'}_{mix}(\rho, T) \]
The first term on the right-hand side is the internal part and the second term is the translational part. The internal part is only a function of temperature, and the translational part is a function of both temperature and density.
The internal component of the thermal conductivity is defined by Equation 2 in [9] :
\[ \lambda^{''}_{mix}(T) = \sum_i \sum_j X_i X_j \lambda^{''}_{ij}(T) \]
The mixing rule for combining pure-species internal thermal conductivity components is given by Equation 3 in [9] :
\[ (\lambda^{''}_{ij})^{-1} = 2[(\lambda^{''}_{i})^{-1} + (\lambda^{''}_{j})^{-1}] \]
The pure species internal thermal conductivity is given by Equation 1 in [9] :
\[ \lambda^{''}_{i} = \frac{\eta_i^*}{MW_i}f_{int} \left(C_{p,i}^0 - 2.5R \right) \]
Where \( \eta_i^* \) is the referred to as the dilute gas viscosity and is the component that is temperature dependent described in [8] (see elyHanleyDilutePureSpeciesViscosity() ), \( MW_i \) is the molecular weight of species i, \( f_{int} \) is a constant and is 1.32, \( C_{p,i}^0 \) is the ideal gas heat capacity of species i, and R is the gas constant (J/kmol/K).
For the translational component of the thermal conductivity, Equation 5 in [9] is used:
\[ \lambda^{'}_{mix}(\rho, T) = \lambda^{'}_{0}(\rho_0, T_0) F_{\lambda} \]
Where \( \lambda^{'}_{0}(\rho_0, T_0) \) is the internal component of the thermal conductivity of a reference fluid that is at a reference temperature and density, \( T_0 \) and \( \rho_0 \). These reference conditions are not the same as the conditions that the mixture is at ( \( T \) and \( \rho \)). The subscript 0 refers to the reference fluid. The term \( F_{\lambda} \) is defined by Equation 6 in [9] :
\[ F_{\lambda} = \left( \frac{MW_0}{MW_m} \right)^{1/2} f_{m,0}^{1/2} h_{m,0}^{-2/3} \]
Where \( MW_0 \) is the molecular weight of the reference fluid, \( MW_m \) is the molecular weight of the mixture.
The terms \( f_{m,0} \) and \( h_{m,0} \) are called reducing ratios. The \( h_{m,0} \) quantity is defined by Equation 9 in [9] :
\[ h_{m,0} = \sum_i \sum_j X_i X_j h_{ij} \]
with \( h_{ij} \) defined by Equation 11 in [9] :
\[ h_{ij} = \frac{1}{8} \left( h_i^{1/3} + h_j^{1/3} \right)^3 (1 - l_{ij}) \]
The variable \( l_{ij} \) is a binary interaction parameter and is assumed to be zero as is done by Ely and Hanley. The pure species reducing ratio \( h_i \) is defined by Equation 13 in [9] :
\[ h_i = \frac{V_c^i}{V_c^0} \phi_i(T_R^i, V_R^i, \omega_i) \]
Where \( \phi_i \) is a shape factor of Leach and Leland (see phiShapeFactor()), \( T_R^i \) is the reduced temperature of species i, \( V_R^i \) is the reduced volume of species i, and \( \omega_i \) is the acentric factor of species i.
At this point, the value of \( h_{m,0} \) can be computed. This leaves the other reducing ratio, \( f_{m,0} \), to be computed. The value of that reducing ratio is defined by Equation 8 in [9] :
\[ f_{m,0} = \frac{1}{h_{m,0}} \sum_i \sum_j X_i X_j f_{ij} h_{ij} \]
The value of \( h_{ij} \) is the same as the one defined earlier. The combining rule for \( f_{ij} \) is given by Equation 10 in [9] :
\[ f_{ij} = (f_i f_j)^{1/2} (1-\kappa_{ij}) \]
\( \kappa_{ij} \) is binary interaction parameters and is assumed to be zero as was done in the work of Ely and Hanley. The species-specific value of \( f_i \) is defined by Equation 12 in [9] :
\[ f_i = \frac{T_c^i}{T_c^0} \theta_i(T_R^i, V_R^i, \omega_i) \]
Where \( \theta_i \) is a shape factor of Leach and Leland (see thetaShapeFactor()), \( T_c^i \) is the critical temperature of species i, and \( T_c^0 \) is the critical temperature of the reference fluid.
The value of \( h_{m,0} \) can be computed from the equation that defined the value of \( f_{m,0} h_{m,0} \). The value of \( h_{m,0} \) must be computed first and then the value of \( f_{m,0} \) can be computed.
The expression for \( MW_m \) is somewhat complex, but keep in mind that all terms except \( MW_m \) are known at this point and so simple algebra can be used to isolate the molecular weight of the mixture. The expression for the mixture molecular weight is given by Equation 14 in [9] :
\[ 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} f_{ij}^{1/2} h_{ij}^{-4/3} \]
where the mixing rule for the molecular weight of the mixture is given by Equation 15 in [9] :
\[ MW_{ij} = \frac{1}{2} (MW_i^{-1} + MW_j^{-1}) \]
For clarity, if we assume the right-hand-side of the molecular weight mixture equation is A, then the mixture molecular weight is given by:
\[ MW_m = A^{-2} f_{m,0} h_{m,0}^{-8/3} \]
All of the above descriptions allow for the calculation of \( F_{\lambda} \) in the expression for the mixture value of the translational component of the thermal conductivity. The reference fluid value of the translational component of the thermal conductivity ( \( \lambda^{'}_{0}(\rho_0, T_0) \) ) is still needed, which requires the the temperature and density of the reference fluid.
The reference fluid density is computed using Equation 7 in [9] :
\[ \rho_0 = \rho h_{m,0} \]
The reference fluid temperature is computed using Equation 7 in [9] :
\[ T_0 = \frac{T}{f_{m,0}} \]
The reference fluid translational component of thermal conductivity can now be computed using Equation 18 in [9]. See elyHanleyReferenceThermalConductivity().
Reimplemented from Transport.
Definition at line 350 of file HighPressureGasTransport.cpp.
|
overridevirtual |
Returns the mixture high-pressure viscosity [Pa·s] using the Lucas method.
This uses the approach described in chapter 9-7. In this method, the mixture viscosity at high pressure is computed using the pure-fluid high pressure relation of Lucas with the difference being that mixture values of the model parameters are used. These mixture values are computed using the mixing rules described in see computeMixtureParameters() .
The mixture pseudo-critical temperature and pressure are calculated using Equations 9-5.18 and 9-5.19. The mixture molecular weight is computed using Equation 9-5.20. The mixture values of the low-pressure polarity and quantum correction factors are computed using Equations 9-5.21 and 9-5.22.
Reimplemented from GasTransport.
Definition at line 555 of file HighPressureGasTransport.cpp.
|
protected |
Get the viscosity [Pa·s] of a pure species using the method of Ely and Hanley.
Uses the method outlined in [8] to compute the viscosity of a pure species.
The primary equation is given by Equation 1 in [8] :
\[ \eta_i(\rho,T) = \eta_0(\rho_0, T_0) F_{\eta} \]
Where \( F_{\eta} \) is defined by Equation 2 in [8] :
\[ F_{\eta} = \left( \frac{MW_i}{MW_0} \right)^{1/2} f_{i,0}^{1/2} h_{i,0}^{-2/3} \]
The 0
subscripts refer to the reference fluid, which is methane in this case, and the i
subscripts refer to the species of interest. No mixing rules need to be used here because this is for a pure species. As such, the terms \( f_{i,0} \) and \( h_{i,0} \) have simpler expressions. The value of \( f_{i,0} \) is given by Equation 7 in [8] :
\[ f_{i,0} = \frac{T_c^i}{T_c^0} \theta_i(T_R^i, V_R^i, \omega_i) \]
and the value of \( h_{i,0} \) is given by Equation 8 in [8] :
\[ h_{i,0} = \frac{V_c^i}{V_c^0} \phi_i(T_R^i, V_R^i, Z_{c,i}, \omega_i) \]
Where \( \theta_i \) and \( \phi_i \) are shape factors of Leach and Leland ( see thetaShapeFactor() and phiShapeFactor() ), \( T_c^i \) is the critical temperature of species i, \( T_c^0 \) is the critical temperature of the reference fluid, \( V_c^i \) is the critical volume of species i, \( V_c^0 \) is the critical volume of the reference fluid, \( Z_{c,i} \) is the critical compressibility of species i, and \( \omega_i \) is the acentric factor of species i.
V | Molar volume of the species |
Tc | Critical temperature of the species |
Vc | Critical volume of the species |
Zc | Critical compressibility of the species |
acentric_factor | Acentric factor of the species |
mw | Molecular weight of the species |
Definition at line 456 of file HighPressureGasTransport.cpp.
|
protected |
Returns the theta shape factor of Leach and Leland for a pure species.
The shape factor \( \theta_i \) is defined by Equation 11 in [8] as follows:
\[ \theta_i(T_R^i, V_R^i, \omega_i) = 1 + (\omega_i - \omega_0)F(T_R^i, V_R^i)] \]
Where \( \omega_i \) is the acentric factor of species i, \( \omega_0 \) is the acentric factor of the reference fluid, and \( F(T_R^i, V_R^i) \) is a function of the reduced temperature and reduced volume of species i. The function \( F(T_R^i, V_R^i) \) is defined by Equation 13 in [8] :
\[ F(T_R^i, V_R^i) = a_1 + b_1 ln(T_+^i) + (c_1 + d_1/T_+^i) (V_+^i - 0.5) \]
Where \( T_+^i \) and \( V_+^i \) are limited values of the reduced temperature and pressure. The limited temperature is defined by Equation 15 in [8] :
\[ T_+^i = \text{min}(2, \text{max}(T_R^i, 0.5)) \]
and the limited pressure is defined by Equation 16 in [8] :
\[ V_i^+ = \text{min}(2, \text{max}(V_R^i, 0.5)) \]
The values of \( a_1 \), \( b_1 \), \( c_1 \), and \( d_1 \) are given in Table II of [8]. They are:
\[ a_1 = 0.090569, b_1 = -0.862762, c_1 = 0.316636, d_1 = -0.465684 \]
Tr | Reduced temperature |
Vr | Reduced volume |
acentric_factor | Acentric factor |
Definition at line 478 of file HighPressureGasTransport.cpp.
|
protected |
Returns the phi shape factor of Leach and Leland for a pure species.
The shape factor \( \phi_i \) is defined by Equation 12 in [8] as follows:
\[ \phi_i(T_R^i, V_R^i, \omega_i, Z_c^i) = \frac{Z_c^0}{Z_c^i} [1 + (\omega_i - \omega_0)G(T_R^i, V_R^i)] \]
Where \( Z_c^0 \) is the critical compressibility of the reference fluid, \( Z_c^i \) is the critical compressibility of species i, \( \omega_0 \) is the acentric factor of the reference fluid, and \( G(T_R^i, V_R^i) \) is a function of the reduced temperature and reduced volume of species i. The function \( G(T_R^i, V_R^i) \) is defined by Equation 14 in [8] :
\[ G(T_R^i, V_R^i) = a_2(V_+^i + b_2) + c_2(V_+^i + d_2)ln(T_+^i) \]
Where \( T_+^i \) and \( V_+^i \) are limited values of the reduced temperature and pressure. The limited temperature is defined by Equation 15 in [8] :
\[ T_+^i = \text{min}(2, \text{max}(T_R^i, 0.5)) \]
and the limited pressure is defined by Equation 16 in [8] :
\[ V_i^+ = \text{min}(2, \text{max}(V_R^i, 0.5)) \]
The values of \( a_2 \), \( b_2 \), \( c_2 \), and \( d_2 \) are given in Table II of [8]. They are:
\[ a_2 = 0.394901, b_2 = -1.023545, c_2 = -0.932813, d_2 = -0.754639 \]
Tr | Reduced temperature |
Vr | Reduced volume |
Zc | Critical compressibility |
acentric_factor | Acentric factor |
Definition at line 489 of file HighPressureGasTransport.cpp.
|
protected |
Returns the viscosity [Pa·s] for the reference fluid (methane) for low pressures.
Computes the temperature dependent (referred to as the dilute viscosity in the reference) component only (eta_0) from the expression in Table III in [8] . Prevents inputs larger than 10,000 Kelvin by just returning the value at 10,000 Kelvin.
T0 | Temperature of the reference fluid |
Definition at line 500 of file HighPressureGasTransport.cpp.
|
protected |
Returns the thermal conductivity [W/m/K] of the reference fluid of methane.
Computes the temperature and density dependent reference fluid thermal conductivity from the expression in Table I in [9] .
The reference fluid translational component of thermal conductivity is computed using Equation 18 in [9] :
\[ \lambda^{'}_{0}(\rho_0, T_0) = \frac{15R}{4MW_0} \eta_0^* + \lambda_0^{(1)}\rho_0 + \Delta\lambda_0(\rho_0, T_0) \]
Where \( \eta_0^* \) is the dilute reference gas viscosity, \( \lambda_0^{(1)} \) is the first density correction to the thermal conductivity, and \( \Delta\lambda_0 \) is the high-density contribution. Ely and Hanley provide a correlation equation to solve for this quantity for the methane reference fluid. The details of the correlation and the parameters are shown in Table I of [9].
For completeness, the relations are reproduced here.
\[ \lambda_0^*(T_0) = \frac{15R}{4MW_0} \eta_0^* \]
\[ \eta_0^* = \sum_{n=1}^{9} C_n T_0^{(n-4)/3} \]
\[ \lambda_0^{(1)} = b_1 + b_2[b_3 - ln(T_0 / b_4)]^2 \]
\[ \Delta\lambda_0 = \text{exp}[a_1 + a_2/T_0] (exp[(a_3 + a_4/T_0^{3/4})\rho_0^{0.1} + (\rho_0 / \rho_0,c - 1)\rho_0^{0.5} (a_5 + a_6/T_0 + a_7/T_0^2) ] - 1) \]
The constants a, b, and C are not related to any ones used in the previous equations, their values are defined in Table I of [9].
rho0 | Density of the reference fluid |
T0 | Temperature of the reference fluid |
Definition at line 521 of file HighPressureGasTransport.cpp.
|
protected |
Computes the composition-dependent values of parameters that are needed for the Lucas viscosity model.
The equations for the mixing rules defined on page 9.23 of [38] for the Lucas method's composition dependent parameters. The primary mixing rules are defined below, and the reduced properties are just the properties divided by the pseudo-critical mixture properties defined below.
\[ T_{\text{c,m}} = \sum_i X_i T_{\text{c,i}} \quad \text{( Equation 9-5.18)} \]
Where \( T_{\text{c,i}} \) is the critical temperature of species i, and \( X_i \) is the mole fraction of species i.
\[ P_{\text{c,m}} = R T_{\text{c,m}} \frac{\sum_i X_i Z_{\text{c,i}}}{\sum_i X_i V_{\text{c,i}}} \quad \text{(Equation 9-5.19)} \]
Where \( Z_{\text{c,i}} \) is the critical compressibility of species i, and \( V_{\text{c,i}} \) is the critical volume of species i.
\[ M_m = \sum X_i M_i \quad \text{(Equation 9-5.20)} \]
Where \( M_i \) is the molecular weight of species i.
\[ F_{P,m}^{\text{o}} = \sum X_i F_{P,i}^{\text{o}} \quad \text{(Equation 9-5.21)} \]
Where \( F_{P,i}^{\text{o}} \) is the low-pressure polarity correction factor of species i from equation 9-4.18.
\[ F_{Q,m}^{\text{o}} = \left ( \sum X_i F_{Q,i}^{\text{o}} \right ) A \quad \text{(Equation 9-5.22)} \]
Where \( F_{Q,i}^{\text{o}} \) is the low-pressure quantum correction factor of species i from equation 9-4.19, and A is defined below.
\[ A = 1 - 0.01 \left ( \frac{M_H}{M_L} \right )^{0.87} \quad \text{(Equation 9-5.23)} \]
For \( \frac{M_H}{M_L} > 9 \) and \( 0.05 < X_H < 0.7 \), otherwise A = 1. In the above equation, $M_H$ and $M_L$ are the molecular weights of the heaviest and lightest components in the mixture, and \( X_H \) is the mole fraction of the heaviest component.
While it isn't returned as a parameter, the species-specific reduced dipole moment ( \( \mu_r \)) is used to compute the mixture polarity correction factor. It is defined as:
\[ \mu_r = 52.46 \frac{\mu^2 P_{\text{c,i}}}{T_{\text{c,i}}} \quad \text{(Equation 9-4.17)} \]
Definition at line 573 of file HighPressureGasTransport.cpp.
|
protected |
Returns the non-dimensional low-pressure mixture viscosity in using the Lucas method.
\[ \eta \xi = F_P^{\text{o}} F_Q^{\text{o}} [0.807 T_r^{0.618} - 0.357 e^{-0.449 T_r} + 0.340e^{-4.058 T_r} + 0.018] \quad \text{(Equation 9-4.16)} \]
This function is structured such that it can be used for pure species or mixtures, with the only difference being the values that are passed to the function (pure values versus mixture values).
For the definition of the mixture rules, see computeMixtureParameters() .
Tr | Reduced temperature [unitless] |
FP | Polarity correction factor [unitless] |
FQ | Quantum correction factor [unitless] |
Definition at line 675 of file HighPressureGasTransport.cpp.
|
protected |
Returns the non-dimensional high-pressure mixture viscosity in using the Lucas method.
\[ \eta \xi = Z_2 F_P F_Q \quad \text{(Equation 9-6.12)} \]
This returns the value of η*ξ (by multiplying both sides of 9-6.12 by ξ and returning the right-side of the resulting equation).
This function is structured such that it can be used for pure species or mixtures, with the only difference being the values that are passed to the function (pure values versus mixture values).
For the definition of the mixture rules, see computeMixtureParameters() .
Tr | Reduced temperature [unitless] |
Pr | Reduced pressure [unitless] |
FP_low | Low-pressure polarity correction factor [unitless] |
FQ_low | Low-pressure quantum correction factor [unitless] |
P_vap | Vapor pressure [Pa] |
P_crit | Critical pressure [Pa] |
Definition at line 683 of file HighPressureGasTransport.cpp.
|
protected |
Calculates quantum correction term of the Lucas method for a species based on the reduced temperature(Tr) and molecular weight(MW), used in viscosity calculation.
\[ F_{Q}^{\text{o}} = 1.22 Q^{0.15} {1 + 0.00385[ (T_r - 12)^2 ]^{\frac{1}{MW}} \text{sign} (T_r - 12 )} \quad \text{(Equation 9-4.19)} \]
Q | Species-specific constant |
Tr | Reduced temperature [unitless] |
MW | Molecular weight [kg/kmol] |
Definition at line 753 of file HighPressureGasTransport.cpp.
|
protected |
Returns the polarity correction term for a species based on reduced temperature, reduced dipole moment, and critical compressibility.
Used in the calculation of viscosity.
Calculates polarity correction term of the Lucas method for a species based on the reduced temperature(Tr) and molecular weight(MW). Equation 9.4.18.
\[ \begin{equation} F_P^0 = \begin{cases} 1 & 0 \leq \mu_r < 0.022 \\ 1 + 30.55(0.292 - Z_c)^{1.72} & 0.022 \leq \mu_r < 0.075 \\ 1 + 30.55(0.292 - Z_c)^{1.72} \times 0.96 + 0.1(T_r - 0.7) & 0.075 \leq \mu_r \end{cases} \end{equation} \]
mu_r | Species Reduced dipole moment |
Tr | Reduced temperature |
Z_c | Species Critical compressibility |
Definition at line 760 of file HighPressureGasTransport.cpp.
|
friend |
Definition at line 376 of file HighPressureGasTransport.h.
|
private |
Molecular weight [kg/kmol].
Definition at line 784 of file HighPressureGasTransport.h.
|
private |
Critical temperature [K].
Definition at line 785 of file HighPressureGasTransport.h.
|
private |
Critical volume [m^3/kmol].
Definition at line 786 of file HighPressureGasTransport.h.
|
private |
Critical compressibility [unitless].
Definition at line 787 of file HighPressureGasTransport.h.
|
private |
Critical density [g/cm^3].
Definition at line 788 of file HighPressureGasTransport.h.
|
private |
Acentric factor [unitless].
Definition at line 789 of file HighPressureGasTransport.h.
|
private |
Quantum correction factor.
Definition at line 798 of file HighPressureGasTransport.h.
|
private |
Polarity correction factor.
Definition at line 799 of file HighPressureGasTransport.h.
|
private |
Reduced temperature.
Definition at line 800 of file HighPressureGasTransport.h.
|
private |
Reduced pressure.
Definition at line 801 of file HighPressureGasTransport.h.
|
private |
Critical pressure.
Definition at line 802 of file HighPressureGasTransport.h.
|
private |
Critical temperature.
Definition at line 803 of file HighPressureGasTransport.h.
|
private |
Molecular weight.
Definition at line 804 of file HighPressureGasTransport.h.
|
private |
Vapor pressure.
Definition at line 805 of file HighPressureGasTransport.h.