Class MixTransport implements mixture-averaged transport properties for ideal gas mixtures. More...
#include <MixTransport.h>
Class MixTransport implements mixture-averaged transport properties for ideal gas mixtures.
The model is based on that described in Kee, et al. [20].
Specific mixture-averaged formulas are implemented by:
Definition at line 35 of file MixTransport.h.
Public Member Functions | |
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 | |
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]. | |
virtual void | getTransportData () |
Read the transport database. | |
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. | |
Protected Attributes | |
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. | |
|
default |
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.
Reimplemented in UnityLewisTransport.
Definition at line 41 of file MixTransport.h.
|
overridevirtual |
Return the thermal diffusion coefficients [kg/m/s].
Model by S. Chapman and T.G. Cowling [4]. For more information about this implementation and its validation, see T. Zirwes and A. Kronenburg [56].
The thermal diffusion coefficient of species \( k \) is computed from
\[ D_k^{T}= \frac{1}{2}\rho\frac{M_k}{\bar{M}}D_{mk}'\Theta_k \]
with
\[ \Theta_k=\frac{15}{2}\frac{\bar{M}^2}{\rho}\sum_i\left(\frac{1.2C_{ki}^*-1}{D_{ki}}\right)\left(\frac{Y_k\frac{\eta_i}{M_i}a_i-Y_i\frac{\eta_k}{M_k}a_k}{M_k+M_i}\right) \]
where \( C_{k,i}^* \) is a reduced collision integral and
\[ a_k=\left(1+\frac{1.065}{2\sqrt{2}X_k}\sum_{i\ne k}X_i\Phi_{k,i}\right)^{-1}, \]
with \( \Phi_{k,i} \) the Wilke mixing operator. The thermodiffusion coefficients are then normalized with
\[ \hat{D}^T_k=D^T_k-Y_k\sum_i D^T_i. \]
This ensures that the sum of all thermodiffusion coefficients and thus the sum of all thermodiffusion fluxes are zero.
[out] | dt | Vector of thermal diffusion coefficients |
Reimplemented from Transport.
Reimplemented in UnityLewisTransport.
Definition at line 51 of file MixTransport.cpp.
|
overridevirtual |
Returns the mixture thermal conductivity [W/m/K].
The thermal conductivity is computed from the following mixture rule:
\[ \lambda = 0.5 \left( \sum_k X_k \lambda_k + \frac{1}{\sum_k X_k/\lambda_k} \right) \]
It's used to compute the flux of energy due to a thermal gradient
\[ \mathbf{q} = - \lambda \nabla T \]
Reimplemented from Transport.
Definition at line 32 of file MixTransport.cpp.
|
overridevirtual |
Get the electrical mobilities [m²/V/s].
This function returns the mobilities. Here, the mobility is calculated from the diffusion coefficient using the Einstein relation:
\[ \mu^e_k = \frac{F D_{km}'}{R T} \]
mobil | Returns the mobilities of the species in array mobil . The array must be dimensioned at least as large as the number of species. |
Reimplemented from Transport.
Definition at line 23 of file MixTransport.cpp.
|
overridevirtual |
Update the internal parameters whenever the temperature has changed.
This is called whenever a transport property is requested if the temperature has changed since the last call to update_T().
Reimplemented from GasTransport.
Definition at line 151 of file MixTransport.cpp.
|
overridevirtual |
Update the internal parameters whenever the concentrations have changed.
This is called whenever a transport property is requested if the concentrations have changed since the last call to update_C().
Implements GasTransport.
Definition at line 164 of file MixTransport.cpp.
|
overridevirtual |
Get the species diffusive mass fluxes [kg/m²/s] with respect to the mass averaged velocity, given the gradients in mole fraction and temperature.
The diffusive mass flux of species k is computed from
\[ \mathbf{j}_k = -\rho \frac{M_k}{\overline{M}} D_{km}' \nabla X_k. \]
ndim | Number of dimensions in the flux expressions | |
[in] | grad_T | Gradient of the temperature (length ndim ) |
ldx | Leading dimension of the grad_X array (usually equal to the number of species) | |
[in] | grad_X | Gradients of the mole fractions; flattened matrix such that \( dX_k/dx_n = \tt{ grad\_X[n*ldx+k]} \) is the gradient of species k in dimension n. Length is ldx * ndim . |
ldf | Leading dimension of the fluxes array (usually equal to the number of species) | |
[out] | fluxes | The diffusive mass fluxes; flattened matrix such that \( j_{kn} = \tt{ fluxes[n*ldf+k]} \) is the flux of species k in dimension n. Length is ldf * ndim . |
Reimplemented from Transport.
Definition at line 126 of file MixTransport.cpp.
|
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 17 of file MixTransport.cpp.
|
protected |
Update the temperature dependent parts of the species thermal conductivities.
These are evaluated from the polynomial fits of the temperature and are assumed to be independent of pressure
Definition at line 178 of file MixTransport.cpp.
|
protected |
vector of species thermal conductivities [W/m/K]
These are used in Wilke's rule to calculate the viscosity of the solution. length = m_nsp.
Definition at line 160 of file MixTransport.h.
|
protected |
Internal storage for the calculated mixture thermal conductivity [W/m/K].
Definition at line 163 of file MixTransport.h.
|
protected |
Update boolean for the species thermal conductivities.
Definition at line 166 of file MixTransport.h.
|
protected |
Update boolean for the mixture rule for the mixture thermal conductivity.
Definition at line 169 of file MixTransport.h.