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. [17].
The viscosity is computed using the Wilke mixture rule (kg /m /s)
\[ \mu = \sum_k \frac{\mu_k X_k}{\sum_j \Phi_{k,j} X_j}. \]
Here \( \mu_k \) is the viscosity of pure species k, and
\[ \Phi_{k,j} = \frac{\left[1 + \sqrt{\left(\frac{\mu_k}{\mu_j}\sqrt{\frac{M_j}{M_k}}\right)}\right]^2} {\sqrt{8}\sqrt{1 + M_k/M_j}} \]
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
\[ j_T = - \lambda \nabla T \]
The flux of energy has units of energy (kg m2 /s2) per second per area.
The units of lambda are W / m K which is equivalent to kg m / s^3 K.
Definition at line 54 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. | |
double | thermalConductivity () override |
Returns the mixture thermal conductivity (W/m /K) | |
void | getMobilities (double *const mobil) override |
Get the Electrical mobilities (m^2/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 wrt to the mass averaged velocity, given the gradients in mole fraction and temperature. | |
void | init (ThermoPhase *thermo, int mode=0, int log_level=0) override |
Initialize a transport manager. | |
Public Member Functions inherited from GasTransport | |
double | viscosity () override |
Viscosity of the mixture (kg /m /s) | |
void | getSpeciesViscosities (double *const visc) override |
Get the pure-species viscosities. | |
void | getBinaryDiffCoeffs (const size_t ld, double *const d) override |
Returns the matrix of binary diffusion coefficients. | |
void | getMixDiffCoeffs (double *const d) override |
Returns the Mixture-averaged diffusion coefficients [m^2/s]. | |
void | getMixDiffCoeffsMole (double *const d) override |
Returns the mixture-averaged diffusion coefficients [m^2/s]. | |
void | getMixDiffCoeffsMass (double *const d) override |
Returns the mixture-averaged diffusion coefficients [m^2/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, int log_level=0) override |
Initialize a transport manager. | |
bool | CKMode () const override |
Boolean indicating the form of the transport properties polynomial fits. | |
Public Member Functions inherited from Transport | |
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 nSpecies(). | |
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 wrt to the specified solution averaged velocity, given the gradients in mole fraction and temperature. | |
virtual void | getMolarFluxes (const double *const state1, const double *const state2, const double delta, double *const cfluxes) |
Get the molar fluxes [kmol/m^2/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^2/s], given the thermodynamic state at two nearby points. | |
virtual void | getThermalDiffCoeffs (double *const dt) |
Return a vector of Thermal diffusion coefficients [kg/m/sec]. | |
virtual void | getBinaryDiffCoeffs (const size_t ld, double *const d) |
Returns the matrix of binary diffusion coefficients [m^2/s]. | |
virtual void | getMultiDiffCoeffs (const size_t ld, double *const d) |
Return the Multicomponent diffusion coefficients. Units: [m^2/s]. | |
virtual void | getMixDiffCoeffs (double *const d) |
Returns a vector of mixture averaged diffusion coefficients. | |
virtual void | getMixDiffCoeffsMole (double *const d) |
Returns a vector of mixture averaged diffusion coefficients. | |
virtual void | getMixDiffCoeffsMass (double *const d) |
Returns a vector of mixture averaged diffusion coefficients. | |
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. | |
virtual double | bulkViscosity () |
The bulk viscosity in Pa-s. | |
virtual double | electricalConductivity () |
The electrical conductivity (Siemens/m). | |
Protected Member Functions | |
void | updateCond_T () |
Update the temperature dependent parts of the species thermal conductivities. | |
Protected Member Functions inherited from GasTransport | |
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 [29]. | |
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. | |
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. | |
Protected Attributes inherited from GasTransport | |
vector< double > | m_molefracs |
Vector of species mole fractions. | |
double | m_viscmix = 0.0 |
Internal storage for the viscosity of the mixture (kg /m /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 |
m_phi is a Viscosity Weighting Function. size = m_nsp * n_nsp | |
vector< double > | m_spwork |
work space length = m_kk | |
vector< double > | m_visc |
vector of species viscosities (kg /m /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 sqrt(kg /m /s). | |
vector< double > | m_polytempvec |
Powers of the ln temperature, up to fourth order. | |
double | m_temp = -1.0 |
Current value of the temperature at which the properties in this object are calculated (Kelvin). | |
double | m_kbt = 0.0 |
Current value of Boltzmann constant times the temperature (Joules) | |
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 nsp x 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 of each species in the phase. | |
vector< double > | m_eps |
Lennard-Jones well-depth of the species in the current phase. | |
vector< double > | m_sigma |
Lennard-Jones diameter of the species in the current phase. | |
DenseMatrix | m_reducedMass |
This is the reduced mass of the interaction between species i and j. | |
DenseMatrix | m_diam |
hard-sphere diameter for (i,j) collision | |
DenseMatrix | m_epsilon |
The effective well depth for (i,j) collisions. | |
DenseMatrix | m_dipole |
The effective dipole moment for (i,j) collisions. | |
DenseMatrix | m_delta |
Reduced dipole moment of the interaction between two species. | |
vector< double > | m_w_ac |
Pitzer acentric factor. | |
vector< double > | m_disp |
Dispersion coefficient. | |
vector< double > | m_quad_polar |
Quadrupole polarizability. | |
int | m_log_level = 0 |
Level of verbose printing during initialization. | |
Protected Attributes inherited from Transport | |
ThermoPhase * | m_thermo |
pointer to the object representing the phase | |
size_t | m_nsp = 0 |
Number of species. | |
|
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 60 of file MixTransport.h.
|
overridevirtual |
Return the thermal diffusion coefficients.
For this approximation, these are all zero.
dt | Vector of thermal diffusion coefficients. Units = kg/m/s |
Reimplemented from Transport.
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
\[ j_T = - \lambda \nabla T \]
The flux of energy has units of energy (kg m2 /s2) per second per area.
The units of lambda are W / m K which is equivalent to kg m / s^3 K.
Reimplemented from Transport.
Definition at line 32 of file MixTransport.cpp.
|
overridevirtual |
Get the Electrical mobilities (m^2/V/s).
This function returns the mobilities. In some formulations this is equal to the normal mobility multiplied by Faraday's constant.
Here, the mobility is calculated from the diffusion coefficient using the Einstein relation
\[ \mu^e_k = \frac{F D_k}{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 83 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 100 of file MixTransport.cpp.
|
overridevirtual |
Get the species diffusive mass fluxes wrt to the mass averaged velocity, given the gradients in mole fraction and temperature.
Units for the returned fluxes are kg m-2 s-1.
The diffusive mass flux of species k is computed from
\[ \vec{j}_k = -n M_k D_k \nabla X_k. \]
ndim | Number of dimensions in the flux expressions |
grad_T | Gradient of the temperature (length = ndim) |
ldx | Leading dimension of the grad_X array (usually equal to m_nsp but not always) |
grad_X | Gradients of the mole fraction. Flat vector with the m_nsp in the inner loop. length = ldx * ndim |
ldf | Leading dimension of the fluxes array (usually equal to m_nsp but not always) |
fluxes | Output of the diffusive mass fluxes. Flat vector with the m_nsp in the inner loop. length = ldx * ndim |
Reimplemented from Transport.
Definition at line 58 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. |
log_level | Defaults to zero, no logging |
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 114 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. units = W /m /K = kg m /s^3 /K. length = m_kk.
Definition at line 166 of file MixTransport.h.
|
protected |
Internal storage for the calculated mixture thermal conductivity.
Units = W /m /K
Definition at line 172 of file MixTransport.h.
|
protected |
Update boolean for the species thermal conductivities.
Definition at line 175 of file MixTransport.h.
|
protected |
Update boolean for the mixture rule for the mixture thermal conductivity.
Definition at line 178 of file MixTransport.h.