Cantera  2.1.2
Public Member Functions | Protected Member Functions | Private Types | Private Member Functions | Private Attributes | Friends | List of all members
LiquidTransport Class Reference

Class LiquidTransport implements models for transport properties for liquid phases. More...

#include <LiquidTransport.h>

Inheritance diagram for LiquidTransport:
[legend]
Collaboration diagram for LiquidTransport:
[legend]

Public Member Functions

 LiquidTransport (thermo_t *thermo=0, int ndim=1)
 Default constructor. More...
 
 LiquidTransport (const LiquidTransport &right)
 
LiquidTransportoperator= (const LiquidTransport &right)
 
virtual TransportduplMyselfAsTransport () const
 Duplication routine for objects which inherit from Transport. More...
 
virtual bool initLiquid (LiquidTransportParams &tr)
 Initialize the transport object. More...
 
virtual int model () const
 Transport model. More...
 
virtual doublereal viscosity ()
 Returns the viscosity of the solution. More...
 
virtual void getSpeciesViscosities (doublereal *const visc)
 Returns the pure species viscosities for all species. More...
 
virtual doublereal ionConductivity ()
 Returns the ionic conductivity of the solution. More...
 
virtual void getSpeciesIonConductivity (doublereal *const ionCond)
 Returns the pure species ionic conductivities for all species. More...
 
virtual void mobilityRatio (doublereal *mobRat)
 Returns the pointer to the mobility ratios of the binary combinations of the transported species for the solution Has size of the number of binary interactions = nsp*nsp. More...
 
virtual void getSpeciesMobilityRatio (doublereal **mobRat)
 Returns a double pointer to the mobility ratios of the transported species in each pure species phase. More...
 
virtual void selfDiffusion (doublereal *const selfDiff)
 Returns the self diffusion coefficients of the species in the phase. More...
 
virtual void getSpeciesSelfDiffusion (doublereal **selfDiff)
 Returns the self diffusion coefficients in the pure species phases. More...
 
virtual void getSpeciesHydrodynamicRadius (doublereal *const radius)
 Returns the hydrodynamic radius for all species. More...
 
virtual void getBinaryDiffCoeffs (const size_t ld, doublereal *const d)
 Returns the binary diffusion coefficients. More...
 
virtual void getMixDiffCoeffs (doublereal *const d)
 Get the Mixture diffusion coefficients. More...
 
virtual void getThermalDiffCoeffs (doublereal *const dt)
 Return the thermal diffusion coefficients. More...
 
virtual doublereal thermalConductivity ()
 Return the thermal conductivity of the solution. More...
 
virtual void getMobilities (doublereal *const mobil_e)
 Get the Electrical mobilities (m^2/V/s). More...
 
virtual void getFluidMobilities (doublereal *const mobil_f)
 Get the fluid mobilities (s kmol/kg). More...
 
virtual void set_Grad_V (const doublereal *const grad_V)
 Specify the value of the gradient of the voltage. More...
 
virtual void set_Grad_T (const doublereal *const grad_T)
 Specify the value of the gradient of the temperature. More...
 
virtual void set_Grad_X (const doublereal *const grad_X)
 Specify the value of the gradient of the MoleFractions. More...
 
virtual doublereal getElectricConduct ()
 Compute the mixture electrical conductivity from the Stefan-Maxwell equation. More...
 
virtual void getElectricCurrent (int ndim, const doublereal *grad_T, int ldx, const doublereal *grad_X, int ldf, const doublereal *grad_V, doublereal *current)
 Compute the electric current density in A/m^2. More...
 
virtual void getSpeciesVdiff (size_t ndim, const doublereal *grad_T, int ldx, const doublereal *grad_X, int ldf, doublereal *Vdiff)
 Get the species diffusive velocities wrt to the averaged velocity, given the gradients in mole fraction and temperature. More...
 
virtual void getSpeciesVdiffES (size_t ndim, const doublereal *grad_T, int ldx, const doublereal *grad_X, int ldf, const doublereal *grad_Phi, doublereal *Vdiff)
 Get the species diffusive velocities wrt to the averaged velocity, given the gradients in mole fraction, temperature and electrostatic potential. More...
 
virtual void getSpeciesFluxes (size_t ndim, const doublereal *const grad_T, size_t ldx, const doublereal *const grad_X, size_t ldf, doublereal *const fluxes)
 Return the species diffusive mass fluxes wrt to the averaged velocity in [kmol/m^2/s]. More...
 
virtual void getSpeciesFluxesES (size_t ndim, const doublereal *grad_T, size_t ldx, const doublereal *grad_X, size_t ldf, const doublereal *grad_Phi, doublereal *fluxes)
 Return the species diffusive mass fluxes wrt to the averaged velocity in [kmol/m^2/s]. More...
 
virtual void getSpeciesVdiffExt (size_t ldf, doublereal *Vdiff)
 Return the species diffusive velocities relative to the averaged velocity. More...
 
virtual void getSpeciesFluxesExt (size_t ldf, doublereal *fluxes)
 Return the species diffusive fluxes relative to the averaged velocity. More...
 
- Public Member Functions inherited from Transport
 Transport (thermo_t *thermo=0, size_t ndim=1)
 Constructor. More...
 
 Transport (const Transport &right)
 
Transportoperator= (const Transport &right)
 
thermo_tthermo ()
 
bool ready ()
 
void setNDim (const int ndim)
 Set the number of dimensions to be expected in flux expressions. More...
 
size_t nDim () const
 Return the number of dimensions in flux expressions. More...
 
void checkSpeciesIndex (size_t k) const
 Check that the specified species index is in range Throws an exception if k is greater than nSpecies() More...
 
void checkSpeciesArraySize (size_t kk) const
 Check that an array size is at least nSpecies() Throws an exception if kk is less than nSpecies(). More...
 
virtual void getMolarFluxes (const doublereal *const state1, const doublereal *const state2, const doublereal delta, doublereal *const cfluxes)
 Get the molar fluxes [kmol/m^2/s], given the thermodynamic state at two nearby points. More...
 
virtual void getMassFluxes (const doublereal *state1, const doublereal *state2, doublereal delta, doublereal *mfluxes)
 Get the mass fluxes [kg/m^2/s], given the thermodynamic state at two nearby points. More...
 
virtual void getMultiDiffCoeffs (const size_t ld, doublereal *const d)
 Return the Multicomponent diffusion coefficients. Units: [m^2/s]. More...
 
virtual void getMixDiffCoeffsMole (doublereal *const d)
 Returns a vector of mixture averaged diffusion coefficients. More...
 
virtual void getMixDiffCoeffsMass (doublereal *const d)
 Returns a vector of mixture averaged diffusion coefficients. More...
 
virtual void setParameters (const int type, const int k, const doublereal *const p)
 Set model parameters for derived classes. More...
 
void setVelocityBasis (VelocityBasis ivb)
 Sets the velocity basis. More...
 
VelocityBasis getVelocityBasis () const
 Gets the velocity basis. More...
 
virtual doublereal bulkViscosity ()
 The bulk viscosity in Pa-s. More...
 
virtual doublereal electricalConductivity ()
 
virtual bool initSolid (SolidTransportData &tr)
 Called by TransportFactory to set parameters. More...
 
virtual void setThermo (thermo_t &thermo)
 Specifies the ThermoPhase object. More...
 

Protected Member Functions

virtual bool update_T ()
 Returns true if temperature has changed, in which case flags are set to recompute transport properties. More...
 
virtual bool update_C ()
 Returns true if mixture composition has changed, in which case flags are set to recompute transport properties. More...
 
virtual void update_Grad_lnAC ()
 Updates the internal value of the gradient of the logarithm of the activity, which is used in the gradient of the chemical potential. More...
 
void stefan_maxwell_solve ()
 Solve the stefan_maxell equations for the diffusive fluxes. More...
 
void updateViscosity_T ()
 Updates the array of pure species viscosities internally. More...
 
void updateIonConductivity_T ()
 Update the temperature-dependent ionic conductivity terms for each species internally. More...
 
void updateMobilityRatio_T ()
 Updates the array of pure species mobility ratios internally. More...
 
void updateSelfDiffusion_T ()
 Updates the array of pure species self diffusion coeffs internally. More...
 
void updateHydrodynamicRadius_T ()
 Update the temperature-dependent hydrodynamic radius terms for each species internally. More...
 
void updateCond_T ()
 Update the temperature-dependent parts of the mixture-averaged thermal conductivity internally. More...
 
void updateViscosities_C ()
 Update the concentration parts of the viscosities. More...
 
void updateIonConductivity_C ()
 Update the concentration parts of the ionic conductivity. More...
 
void updateMobilityRatio_C ()
 Update the concentration parts of the mobility ratio. More...
 
void updateSelfDiffusion_C ()
 Update the concentration parts of the self diffusion. More...
 
void updateHydrodynamicRadius_C ()
 Update the concentration dependence of the hydrodynamic radius. More...
 
void updateDiff_T ()
 Update the binary Stefan-Maxwell diffusion coefficients wrt T using calls to the appropriate LTPspecies subclass. More...
 
- Protected Member Functions inherited from Transport
virtual bool initGas (GasTransportParams &tr)
 Called by TransportFactory to set parameters. More...
 
void finalize ()
 Enable the transport object for use. More...
 

Private Types

typedef std::vector< LTPspecies * > LTPvector
 Type def for LTPvector equating it with a vector of pointers to LTPspecies. More...
 

Private Member Functions

doublereal err (const std::string &msg) const
 Throw an exception if this method is invoked. More...
 

Private Attributes

size_t m_nsp2
 Number of species squared. More...
 
vector_fp m_mw
 Local copy of the molecular weights of the species. More...
 
std::vector< LTPspecies * > m_viscTempDep_Ns
 Viscosity for each species expressed as an appropriate subclass of LTPspecies. More...
 
LiquidTranInteractionm_viscMixModel
 Viscosity of the mixture expressed as a subclass of LiquidTranInteraction. More...
 
std::vector< LTPspecies * > m_ionCondTempDep_Ns
 Ionic conductivity for each species expressed as an appropriate subclass of LTPspecies. More...
 
LiquidTranInteractionm_ionCondMixModel
 Ionic Conductivity of the mixture expressed as a subclass of LiquidTranInteraction. More...
 
std::vector< LTPvectorm_mobRatTempDep_Ns
 Mobility ratio for the binary cominations of each species in each pure phase expressed as an appropriate subclass of LTPspecies. More...
 
std::vector
< LiquidTranInteraction * > 
m_mobRatMixModel
 Mobility ratio for each binary combination of mobile species in the mixture expressed as a subclass of LiquidTranInteraction. More...
 
std::vector< LTPvectorm_selfDiffTempDep_Ns
 Self Diffusion for each species in each pure species phase expressed as an appropriate subclass of LTPspecies. More...
 
std::vector
< LiquidTranInteraction * > 
m_selfDiffMixModel
 Self Diffusion for each species in the mixture expressed as a subclass of LiquidTranInteraction. More...
 
std::vector< LTPspecies * > m_lambdaTempDep_Ns
 Thermal conductivity for each species expressed as an appropriate subclass of LTPspecies. More...
 
LiquidTranInteractionm_lambdaMixModel
 Thermal conductivity of the mixture expressed as a subclass of LiquidTranInteraction. More...
 
std::vector< LTPspecies * > m_diffTempDep_Ns
 (NOT USED IN LiquidTransport.) Diffusion coefficient model for each species expressed as an appropriate subclass of LTPspecies More...
 
LiquidTranInteractionm_diffMixModel
 Species diffusivity of the mixture expressed as a subclass of LiquidTranInteraction. More...
 
DenseMatrix m_diff_Dij
 Setfan-Maxwell diffusion coefficients. More...
 
std::vector< LTPspecies * > m_radiusTempDep_Ns
 Hydrodynamic radius for each species expressed as an appropriate subclass of LTPspecies. More...
 
LiquidTranInteractionm_radiusMixModel
 (Not used in LiquidTransport) Hydrodynamic radius of the mixture expressed as a subclass of LiquidTranInteraction More...
 
vector_fp m_hydrodynamic_radius
 Species hydrodynamic radius. More...
 
vector_fp m_Grad_X
 Internal value of the gradient of the mole fraction vector. More...
 
vector_fp m_Grad_lnAC
 Gradient of the logarithm of the activity. More...
 
vector_fp m_Grad_T
 Internal value of the gradient of the Temperature vector. More...
 
vector_fp m_Grad_P
 Internal value of the gradient of the Pressure vector. More...
 
vector_fp m_Grad_V
 Internal value of the gradient of the Electric Voltage. More...
 
vector_fp m_Grad_mu
 Gradient of the electrochemical potential. More...
 
DenseMatrix m_bdiff
 Array of Binary Diffusivities. More...
 
vector_fp m_viscSpecies
 Internal value of the species viscosities. More...
 
vector_fp m_ionCondSpecies
 Internal value of the species ionic conductivities. More...
 
DenseMatrix m_mobRatSpecies
 Internal value of the species mobility ratios. More...
 
DenseMatrix m_selfDiffSpecies
 Internal value of the species self diffusion coefficients. More...
 
vector_fp m_lambdaSpecies
 Internal value of the species individual thermal conductivities. More...
 
int m_iStateMF
 State of the mole fraction vector. More...
 
vector_fp m_massfracs
 Local copy of the mass fractions of the species in the phase. More...
 
vector_fp m_massfracs_tran
 Local copy of the mass fractions of the species in the phase. More...
 
vector_fp m_molefracs
 Local copy of the mole fractions of the species in the phase. More...
 
vector_fp m_molefracs_tran
 Non-zero mole fraction vector used in transport property calculations. More...
 
vector_fp m_concentrations
 Local copy of the concentrations of the species in the phase. More...
 
doublereal concTot_
 Local copy of the total concentration. More...
 
doublereal concTot_tran_
 Local copy of the total concentration. More...
 
doublereal meanMolecularWeight_
 Mean molecular mass. More...
 
doublereal dens_
 Density. More...
 
vector_fp m_chargeSpecies
 Local copy of the charge of each species. More...
 
vector_fp m_volume_spec
 Specific volume for each species. Local copy from thermo object. More...
 
vector_fp m_actCoeff
 Vector of activity coefficients. More...
 
DenseMatrix m_B
 RHS to the stefan-maxwell equation. More...
 
DenseMatrix m_A
 Matrix for the stefan maxwell equation. More...
 
doublereal m_temp
 Current Temperature -> locally stored. More...
 
doublereal m_press
 Current value of the pressure. More...
 
Array2D m_flux
 Solution of the Stefan Maxwell equation in terms of flux. More...
 
Array2D m_Vdiff
 Solution of the Stefan Maxwell equation. More...
 
doublereal m_lambda
 Saved value of the mixture thermal conductivity. More...
 
doublereal m_viscmix
 Saved value of the mixture viscosity. More...
 
doublereal m_ionCondmix
 Saved value of the mixture ionic conductivity. More...
 
vector_fp m_mobRatMix
 Saved values of the mixture mobility ratios. More...
 
vector_fp m_selfDiffMix
 Saved values of the mixture self diffusion coefficients. More...
 
vector_fp m_spwork
 work space More...
 
bool m_visc_mix_ok
 Boolean indicating that the top-level mixture viscosity is current. More...
 
bool m_visc_temp_ok
 Boolean indicating that weight factors wrt viscosity is current. More...
 
bool m_visc_conc_ok
 Flag to indicate that the pure species viscosities are current wrt the concentration. More...
 
bool m_ionCond_mix_ok
 Boolean indicating that the top-level mixture ionic conductivity is current. More...
 
bool m_ionCond_temp_ok
 Boolean indicating that weight factors wrt ionic conductivty is current. More...
 
bool m_ionCond_conc_ok
 Flag to indicate that the pure species ionic conductivities are current wrt the concentration. More...
 
bool m_cond_mix_ok
 Flag to indicate that the mixture conductivity is current. More...
 
bool m_mobRat_mix_ok
 Boolean indicating that the top-level mixture mobility ratio is current. More...
 
bool m_mobRat_temp_ok
 Boolean indicating that weight factors wrt mobility ratio is current. More...
 
bool m_mobRat_conc_ok
 Flag to indicate that the pure species mobility ratios are current wrt the concentration. More...
 
bool m_selfDiff_mix_ok
 Boolean indicating that the top-level mixture self diffusion is current. More...
 
bool m_selfDiff_temp_ok
 Boolean indicating that weight factors wrt self diffusion is current. More...
 
bool m_selfDiff_conc_ok
 Flag to indicate that the pure species self diffusion are current wrt the concentration. More...
 
bool m_radi_mix_ok
 Boolean indicating that mixture diffusion coeffs are current. More...
 
bool m_radi_temp_ok
 Boolean indicating that temperature dependence of hydrodynamic radius is current. More...
 
bool m_radi_conc_ok
 Flag to indicate that the hydrodynamic radius is current is current wrt the concentration. More...
 
bool m_diff_mix_ok
 Boolean indicating that mixture diffusion coeffs are current. More...
 
bool m_diff_temp_ok
 Boolean indicating that binary diffusion coeffs are current. More...
 
bool m_lambda_temp_ok
 Flag to indicate that the pure species conductivities are current wrt the temperature. More...
 
bool m_lambda_mix_ok
 Boolean indicating that mixture conductivity is current. More...
 
int m_mode
 Mode indicator for transport models – currently unused. More...
 
bool m_debug
 Debugging flags. More...
 

Friends

class TransportFactory
 

Additional Inherited Members

- Protected Attributes inherited from Transport
thermo_tm_thermo
 pointer to the object representing the phase More...
 
bool m_ready
 true if finalize has been called More...
 
size_t m_nsp
 Number of species. More...
 
size_t m_nDim
 Number of dimensions used in flux expressions. More...
 
int m_velocityBasis
 Velocity basis from which diffusion velocities are computed. More...
 

Detailed Description

Class LiquidTransport implements models for transport properties for liquid phases.

Liquid Transport is set up with some flexibility in this class. Transport properties like viscosity and thermal conductivity are allowed flexibility within the constraints of the LiquidTransportProperty and LiquidTransportInteractions classes. For species diffusion, the LiquidTransport class focuses on the Stefan-Maxwell equation to determine the diffusion velocities. Other options for liquid diffusion include solvent-dominated diffusion, and a class SolventTransport should be forthcoming.

The class LiquidTransport has several roles.

  1. It brings together the individual species transport properties, expressed as subclasses of LTPspecies (Liquid Transport Properties of Species) through LiquidTransportData, with models for the composition dependence of liquid transport properties expressed as subclasses of LiquidTranInteraction (mixing rules) through LiquidTransportParams. Calculating mixture properties generally consists of calling the getMixTansProp member of LiquidTranInteraction by passing a vector of LTPSpecies
  2. It calculates the bulk velocity \( \vec{v} \) and individual species diffusion velocities, \( \vec{V_i} \) using the Stefan-Maxwell equations. It is possible to set a flag to calculate relative to a mass-averaged bulk velocity, relative to a mole-averaged bulk velocity or relative to a single species velocity using the <velocityBasis basis="mass">, <velocityBasis basis="mass">, or <velocityBasis basis="Cl-"> keyword. Mass-averaged velocities are the default for which the diffusion velocities satisfy

    \[ \sum_{i} Y_i \vec{V_i} = 0 \]

    for mass fraction \( Y_i \). For mole-averaged velocities

    \[ \sum_{i} X_i \vec{V_i} = 0 \]

    for mole fraction \( X_i \). or

    \[ \vec{V_i} = 0 \]

    for reference species \( i \).
  3. It provides access to a number of derived quantities related to transport properties as described in the various methods below.

Within LiquidTransport, the state is presumed to be defined in terms of the species mole fraction, temperature and pressure. Charged species are expected and quantities like the electric current are computed based on a combined electrochemcial potential.

Definition at line 76 of file LiquidTransport.h.

Member Typedef Documentation

typedef std::vector<LTPspecies*> LTPvector
private

Type def for LTPvector equating it with a vector of pointers to LTPspecies.

Definition at line 842 of file LiquidTransport.h.

Constructor & Destructor Documentation

LiquidTransport ( thermo_t thermo = 0,
int  ndim = 1 
)

Default constructor.

This requires call to initLiquid(LiquidTransportParams& tr) after filling LiquidTransportParams to complete instantiation. The filling of LiquidTransportParams is currently carried out in the TransportFactory class, but might be moved at some point.

Parameters
thermoThermoPhase object holding species information.
ndimNumber of spatial dimensions.

Definition at line 20 of file LiquidTransport.cpp.

Referenced by LiquidTransport::duplMyselfAsTransport().

Member Function Documentation

Transport * duplMyselfAsTransport ( ) const
virtual

Duplication routine for objects which inherit from Transport.

This virtual routine can be used to duplicate objects derived from Transport even if the application only has a pointer to Transport to work with.

These routines are basically wrappers around the derived copy constructor.

Reimplemented from Transport.

Definition at line 191 of file LiquidTransport.cpp.

References LiquidTransport::LiquidTransport().

bool initLiquid ( LiquidTransportParams tr)
virtual

Initialize the transport object.

Here we change all of the internal dimensions to be sufficient. We get the object ready to do property evaluations. A lot of the input required to do property evaluations is contained in the LiquidTransportParams class that is filled in TransportFactory.

Parameters
trTransport parameters for all of the species in the phase.

Reimplemented from Transport.

Definition at line 226 of file LiquidTransport.cpp.

References Phase::charge(), LiquidTransportData::hydroRadius, LiquidTransportParams::ionConductivity, LiquidTransportData::ionConductivity, LiquidTransportParams::LTData, LiquidTransport::m_actCoeff, LiquidTransport::m_bdiff, LiquidTransport::m_chargeSpecies, LiquidTransport::m_concentrations, LiquidTransport::m_diff_mix_ok, LiquidTransport::m_diff_temp_ok, LiquidTransport::m_diffMixModel, LiquidTransport::m_diffTempDep_Ns, LiquidTransport::m_flux, LiquidTransport::m_Grad_lnAC, LiquidTransport::m_Grad_mu, LiquidTransport::m_Grad_T, LiquidTransport::m_Grad_V, LiquidTransport::m_Grad_X, LiquidTransport::m_hydrodynamic_radius, LiquidTransport::m_ionCond_conc_ok, LiquidTransport::m_ionCond_mix_ok, LiquidTransport::m_ionCond_temp_ok, LiquidTransport::m_ionCondMixModel, LiquidTransport::m_ionCondSpecies, LiquidTransport::m_ionCondTempDep_Ns, LiquidTransport::m_lambda_mix_ok, LiquidTransport::m_lambda_temp_ok, LiquidTransport::m_lambdaMixModel, LiquidTransport::m_lambdaSpecies, LiquidTransport::m_lambdaTempDep_Ns, LiquidTransport::m_massfracs, LiquidTransport::m_massfracs_tran, LiquidTransport::m_mobRat_conc_ok, LiquidTransport::m_mobRat_mix_ok, LiquidTransport::m_mobRat_temp_ok, LiquidTransport::m_mobRatMix, LiquidTransport::m_mobRatMixModel, LiquidTransport::m_mobRatSpecies, LiquidTransport::m_mobRatTempDep_Ns, LiquidTransport::m_mode, LiquidTransport::m_molefracs, LiquidTransport::m_molefracs_tran, LiquidTransport::m_mw, Transport::m_nDim, Transport::m_nsp, LiquidTransport::m_nsp2, LiquidTransport::m_radi_conc_ok, LiquidTransport::m_radi_temp_ok, LiquidTransport::m_radiusTempDep_Ns, LiquidTransport::m_selfDiff_conc_ok, LiquidTransport::m_selfDiff_mix_ok, LiquidTransport::m_selfDiff_temp_ok, LiquidTransport::m_selfDiffMix, LiquidTransport::m_selfDiffMixModel, LiquidTransport::m_selfDiffSpecies, LiquidTransport::m_selfDiffTempDep_Ns, LiquidTransport::m_spwork, Transport::m_thermo, LiquidTransport::m_Vdiff, Transport::m_velocityBasis, LiquidTransport::m_visc_conc_ok, LiquidTransport::m_visc_mix_ok, LiquidTransport::m_visc_temp_ok, LiquidTransport::m_viscMixModel, LiquidTransport::m_viscSpecies, LiquidTransport::m_viscTempDep_Ns, LiquidTransport::m_volume_spec, LiquidTransportParams::mobilityRatio, LiquidTransportData::mobilityRatio, TransportParams::mode_, Phase::molecularWeights(), Phase::nSpecies(), DenseMatrix::resize(), Array2D::resize(), LiquidTransportParams::selfDiffusion, LiquidTransportData::selfDiffusion, LiquidTransportParams::speciesDiffusivity, LiquidTransportData::speciesDiffusivity, Phase::speciesName(), LiquidTransportParams::thermalCond, LiquidTransportData::thermalCond, TransportParams::thermo, TransportParams::velocityBasis_, LiquidTransportParams::viscosity, and LiquidTransportData::viscosity.

virtual int model ( ) const
inlinevirtual

Transport model.

The transport model is the set of equations used to compute the transport properties. This method returns an integer flag that identifies the transport model implemented. The base class returns 0.

Reimplemented from Transport.

Definition at line 110 of file LiquidTransport.h.

Referenced by LiquidTransport::err().

doublereal viscosity ( )
virtual

Returns the viscosity of the solution.

The viscosity calculation is handled by subclasses of LiquidTranInteraction as specified in the input file. These in turn employ subclasses of LTPspecies to determine the individual species viscosities.

Reimplemented from Transport.

Definition at line 402 of file LiquidTransport.cpp.

References LiquidTranInteraction::getMixTransProp(), LiquidTransport::m_visc_mix_ok, LiquidTransport::m_viscmix, LiquidTransport::m_viscMixModel, LiquidTransport::m_viscTempDep_Ns, LiquidTransport::update_C(), and LiquidTransport::update_T().

void getSpeciesViscosities ( doublereal *const  visc)
virtual

Returns the pure species viscosities for all species.

The pure species viscosities are evaluated using the appropriate subclasses of LTPspecies as specified in the input file.

Parameters
viscarray of length "number of species" to hold returned viscosities.

Reimplemented from Transport.

Definition at line 417 of file LiquidTransport.cpp.

References LiquidTransport::m_visc_temp_ok, LiquidTransport::m_viscSpecies, LiquidTransport::update_T(), and LiquidTransport::updateViscosity_T().

doublereal ionConductivity ( )
virtual

Returns the ionic conductivity of the solution.

The ionic conductivity calculation is handled by subclasses of LiquidTranInteraction as specified in the input file. These in turn employ subclasses of LTPspecies to determine the individual species ionic conductivities.

Reimplemented from Transport.

Definition at line 426 of file LiquidTransport.cpp.

References LiquidTranInteraction::getMixTransProp(), LiquidTransport::m_ionCond_mix_ok, LiquidTransport::m_ionCondmix, LiquidTransport::m_ionCondMixModel, LiquidTransport::m_ionCondTempDep_Ns, LiquidTransport::update_C(), and LiquidTransport::update_T().

void getSpeciesIonConductivity ( doublereal *const  ionCond)
virtual

Returns the pure species ionic conductivities for all species.

The pure species ionic conductivities are evaluated using the appropriate subclasses of LTPspecies as specified in the input file.

Parameters
ionCondArray of length "number of species" to hold returned ionic conductivities.

Reimplemented from Transport.

Definition at line 452 of file LiquidTransport.cpp.

References LiquidTransport::m_ionCond_temp_ok, LiquidTransport::m_ionCondSpecies, LiquidTransport::update_T(), and LiquidTransport::updateIonConductivity_T().

void mobilityRatio ( doublereal *  mobRat)
virtual

Returns the pointer to the mobility ratios of the binary combinations of the transported species for the solution Has size of the number of binary interactions = nsp*nsp.

The mobility ratio calculation is handled by subclasses of LiquidTranInteraction as specified in the input file. These in turn employ subclasses of LTPspecies to determine the mobility ratios in the pure species.

Parameters
mobRatVector of mobility ratios

Reimplemented from Transport.

Definition at line 461 of file LiquidTransport.cpp.

References LiquidTransport::m_mobRat_mix_ok, LiquidTransport::m_mobRatMix, LiquidTransport::m_mobRatMixModel, LiquidTransport::m_mobRatTempDep_Ns, Transport::m_nsp, LiquidTransport::m_nsp2, LiquidTransport::update_C(), and LiquidTransport::update_T().

void getSpeciesMobilityRatio ( doublereal **  mobRat)
virtual

Returns a double pointer to the mobility ratios of the transported species in each pure species phase.

Has size of the number of binary interactions by the number of species (nsp*nsp X nsp). The pure species mobility ratios are evaluated using the appropriate subclasses of LTPspecies as specified in the input file.

Parameters
mobRatarray of length "number of species" to hold returned mobility ratios.

Reimplemented from Transport.

Definition at line 483 of file LiquidTransport.cpp.

References LiquidTransport::m_mobRat_temp_ok, LiquidTransport::m_mobRatSpecies, Transport::m_nsp, LiquidTransport::m_nsp2, LiquidTransport::update_T(), and LiquidTransport::updateMobilityRatio_T().

void selfDiffusion ( doublereal *const  selfDiff)
virtual

Returns the self diffusion coefficients of the species in the phase.

Has size of nsp(coeffs)

The self diffusion coefficient is the diffusion coefficient of a tracer species at the current temperature and composition of the species. Therefore, the dilute limit of transport is assumed for the tracer species. The effective formula may be calculated from the stefan-maxwell formulation by adding another row for the tracer species, assigning all D's to be equal to the respective species D's, and then taking the limit as the tracer species mole fraction goes to zero. The corresponding flux equation for the tracer species k in units of kmol m-2 s-1 is.

\[ J_k = - D^{sd}_k \frac{C_k}{R T} \nabla \mu_k \]

The derivative is taken at constant T and P.

The self diffusion calculation is handled by subclasses of LiquidTranInteraction as specified in the input file. These in turn employ subclasses of LTPspecies to determine the individual species self diffusion coeffs.

Parameters
selfDiffVector of self-diffusion coefficients Length = number of species in phase units = m**2 s-1

Reimplemented from Transport.

Definition at line 496 of file LiquidTransport.cpp.

References Transport::m_nsp, LiquidTransport::m_selfDiff_mix_ok, LiquidTransport::m_selfDiffMix, LiquidTransport::m_selfDiffMixModel, LiquidTransport::m_selfDiffTempDep_Ns, LiquidTransport::update_C(), and LiquidTransport::update_T().

void getSpeciesSelfDiffusion ( doublereal **  selfDiff)
virtual

Returns the self diffusion coefficients in the pure species phases.

Has size of nsp(coeffs) x nsp(phases)

The pure species molar volumes are evaluated using the appropriate subclasses of LTPspecies as specified in the input file.

Parameters
selfDiffarray of length "number of species" to hold returned self diffusion coeffs.

Reimplemented from Transport.

Definition at line 510 of file LiquidTransport.cpp.

References Transport::m_nsp, LiquidTransport::m_selfDiff_temp_ok, LiquidTransport::m_selfDiffSpecies, LiquidTransport::update_T(), and LiquidTransport::updateSelfDiffusion_T().

void getSpeciesHydrodynamicRadius ( doublereal *const  radius)
virtual

Returns the hydrodynamic radius for all species.

The species hydrodynamic radii are evaluated using the appropriate subclasses of LTPspecies as specified in the input file.

Parameters
radiusarray of length "number of species" to hold returned radii.

Definition at line 523 of file LiquidTransport.cpp.

References LiquidTransport::m_hydrodynamic_radius, LiquidTransport::m_radi_temp_ok, LiquidTransport::update_T(), and LiquidTransport::updateHydrodynamicRadius_T().

void getBinaryDiffCoeffs ( const size_t  ld,
doublereal *const  d 
)
virtual

Returns the binary diffusion coefficients.

The binary diffusion coefficients are specified in the input file through the LiquidTransportInteractions class. These are the binary interaction coefficients employed in the Stefan-Maxwell equation.

Parameters
ldnumber of species in system
dvector of binary diffusion coefficients units = m2 s-1. length = ld*ld = (number of species)^2

Reimplemented from Transport.

Definition at line 554 of file LiquidTransport.cpp.

References LiquidTransport::m_bdiff, LiquidTransport::m_diff_temp_ok, Transport::m_nsp, LiquidTransport::update_T(), and LiquidTransport::updateDiff_T().

void getMixDiffCoeffs ( doublereal *const  d)
virtual

Get the Mixture diffusion coefficients.

The mixture diffusion coefficients are not well defined in the context of LiquidTransport because the Stefan Maxwell equation is solved. Here the mixture diffusion coefficients are defined according to Ficks law:

\[ X_i \vec{V_i} = -D_i \nabla X_i. \]

Solving Ficks Law for \( D_i \) gives a mixture diffusion coefficient

\[ D_i = - X_i \vec{V_i} / ( \nabla X_i ). \]

If \( \nabla X_i = 0 \) this is undefined and the nonsensical value -1 is returned.

Note that this evaluation of \( \vec{V_i} \) requires a solve of the Stefan Maxwell equation making this determination of the mixture averaged diffusion coefficients a slow method for obtaining diffusion coefficients.

Also note that the Stefan Maxwell solve will be based upon the thermodynamic state (including gradients) most recently set. Gradients can be set specifically using set_Grad_V, set_Grad_X and set_Grad_T or through calls to getSpeciesFluxes, getSpeciesFluxesES, getSpeciesVdiff, getSpeciesVdiffES, etc.

Parameters
dvector of mixture diffusion coefficients units = m2 s-1. length = number of species

Reimplemented from Transport.

Definition at line 748 of file LiquidTransport.cpp.

References LiquidTransport::m_Grad_X, LiquidTransport::m_molefracs, Transport::m_nDim, Transport::m_nsp, LiquidTransport::m_Vdiff, and LiquidTransport::stefan_maxwell_solve().

Referenced by LiquidTransport::getFluidMobilities(), and LiquidTransport::getMobilities().

void getThermalDiffCoeffs ( doublereal *const  dt)
virtual

Return the thermal diffusion coefficients.

These are all zero for this simple implementation

Parameters
dtthermal diffusion coefficients

Reimplemented from Transport.

Definition at line 547 of file LiquidTransport.cpp.

References Transport::m_nsp.

doublereal thermalConductivity ( )
virtual

Return the thermal conductivity of the solution.

The thermal conductivity calculation is handled by subclasses of LiquidTranInteraction as specified in the input file. These in turn employ subclasses of LTPspecies to determine the individual species thermal conductivities.

Reimplemented from Transport.

Definition at line 533 of file LiquidTransport.cpp.

References LiquidTranInteraction::getMixTransProp(), LiquidTransport::m_cond_mix_ok, LiquidTransport::m_lambda, LiquidTransport::m_lambda_mix_ok, LiquidTransport::m_lambdaMixModel, LiquidTransport::m_lambdaTempDep_Ns, LiquidTransport::update_C(), and LiquidTransport::update_T().

void getMobilities ( doublereal *const  mobil_e)
virtual

Get the Electrical mobilities (m^2/V/s).

The electrical mobilities are not well defined in the context of LiquidTransport because the Stefan Maxwell equation is solved. Here the electrical mobilities are calculated from the mixture-averaged diffusion coefficients through a call to getMixDiffCoeffs() using the Einstein relation

\[ \mu^e_k = \frac{F D_k}{R T} \]

Note that this call to getMixDiffCoeffs() requires a solve of the Stefan Maxwell equation making this determination of the mixture averaged diffusion coefficients a slow method for obtaining diffusion coefficients.

Also note that the Stefan Maxwell solve will be based upon the thermodynamic state (including gradients) most recently set. Gradients can be set specifically using set_Grad_V, set_Grad_X and set_Grad_T or through calls to getSpeciesFluxes, getSpeciesFluxesES, getSpeciesVdiff, getSpeciesVdiffES, etc.

Parameters
mobil_eReturns the electrical mobilities of the species in array mobil_e. The array must be dimensioned at least as large as the number of species.

Reimplemented from Transport.

Definition at line 578 of file LiquidTransport.cpp.

References Cantera::Boltzmann, DATA_PTR, LiquidTransport::getMixDiffCoeffs(), Transport::m_nsp, LiquidTransport::m_spwork, and LiquidTransport::m_temp.

void getFluidMobilities ( doublereal *const  mobil_f)
virtual

Get the fluid mobilities (s kmol/kg).

The fluid mobilities are not well defined in the context of LiquidTransport because the Stefan Maxwell equation is solved. Here the fluid mobilities are calculated from the mixture-averaged diffusion coefficients through a call to getMixDiffCoeffs() using the Einstein relation

\[ \mu^f_k = \frac{D_k}{R T} \]

Note that this call to getMixDiffCoeffs() requires a solve of the Stefan Maxwell equation making this determination of the mixture averaged diffusion coefficients a slow method for obtaining diffusion coefficients.

Also note that the Stefan Maxwell solve will be based upon the thermodynamic state (including gradients) most recently set. Gradients can be set specifically using set_Grad_V, set_Grad_X and set_Grad_T or through calls to getSpeciesFluxes, getSpeciesFluxesES, getSpeciesVdiff, getSpeciesVdiffES, etc.

Parameters
mobil_fReturns the fluid mobilities of the species in array mobil_f. The array must be dimensioned at least as large as the number of species.

Reimplemented from Transport.

Definition at line 587 of file LiquidTransport.cpp.

References DATA_PTR, Cantera::GasConstant, LiquidTransport::getMixDiffCoeffs(), Transport::m_nsp, LiquidTransport::m_spwork, and LiquidTransport::m_temp.

void set_Grad_V ( const doublereal *const  grad_V)
virtual

Specify the value of the gradient of the voltage.

Parameters
grad_VGradient of the voltage (length num dimensions);

Definition at line 603 of file LiquidTransport.cpp.

References LiquidTransport::m_Grad_V, and Transport::m_nDim.

Referenced by LiquidTransport::getElectricConduct(), LiquidTransport::getElectricCurrent(), LiquidTransport::getSpeciesFluxesES(), and LiquidTransport::getSpeciesVdiffES().

void set_Grad_T ( const doublereal *const  grad_T)
virtual

Specify the value of the gradient of the temperature.

Parameters
grad_TGradient of the temperature (length num dimensions);

Definition at line 596 of file LiquidTransport.cpp.

References LiquidTransport::m_Grad_T, and Transport::m_nDim.

Referenced by LiquidTransport::getElectricConduct(), LiquidTransport::getElectricCurrent(), LiquidTransport::getSpeciesFluxes(), LiquidTransport::getSpeciesFluxesES(), LiquidTransport::getSpeciesVdiff(), and LiquidTransport::getSpeciesVdiffES().

void set_Grad_X ( const doublereal *const  grad_X)
virtual

Specify the value of the gradient of the MoleFractions.

Parameters
grad_XGradient of the mole fractions(length nsp * num dimensions);

Definition at line 610 of file LiquidTransport.cpp.

References LiquidTransport::m_Grad_X, Transport::m_nDim, and Transport::m_nsp.

Referenced by LiquidTransport::getElectricConduct(), LiquidTransport::getElectricCurrent(), LiquidTransport::getSpeciesFluxes(), LiquidTransport::getSpeciesFluxesES(), LiquidTransport::getSpeciesVdiff(), and LiquidTransport::getSpeciesVdiffES().

doublereal getElectricConduct ( )
virtual

Compute the mixture electrical conductivity from the Stefan-Maxwell equation.

To compute the mixture electrical conductance, the Stefan Maxwell equation is solved for zero species gradients and for unit potential gradient, \( \nabla V \). The species fluxes are converted to current by summing over the charge-weighted fluxes according to

\[ \vec{i} = \sum_{i} z_i F \rho \vec{V_i} / W_i \]

where \( z_i \) is the charge on species i, \( F \) is Faradays constant, \( \rho \) is the density, \( W_i \) is the molecular mass of species i. The conductance, \( \kappa \) is obtained from

\[ \kappa = \vec{i} / \nabla V. \]

Reimplemented from Transport.

Definition at line 618 of file LiquidTransport.cpp.

References LiquidTransport::getSpeciesFluxesExt(), LiquidTransport::m_chargeSpecies, LiquidTransport::m_mw, Transport::m_nDim, Transport::m_nsp, LiquidTransport::set_Grad_T(), LiquidTransport::set_Grad_V(), and LiquidTransport::set_Grad_X().

void getElectricCurrent ( int  ndim,
const doublereal *  grad_T,
int  ldx,
const doublereal *  grad_X,
int  ldf,
const doublereal *  grad_V,
doublereal *  current 
)
virtual

Compute the electric current density in A/m^2.

The electric current is computed first by computing the species diffusive fluxes using the Stefan Maxwell solution and then the current, \( \vec{i} \) by summing over the charge-weighted fluxes according to

\[ \vec{i} = \sum_{i} z_i F \rho \vec{V_i} / W_i \]

where \( z_i \) is the charge on species i, \( F \) is Faradays constant, \( \rho \) is the density, \( W_i \) is the molecular mass of species i.

Parameters
ndimThe number of spatial dimensions (1, 2, or 3).
grad_TThe temperature gradient (ignored in this model).
ldxLeading dimension of the grad_X array.
grad_XGradients of the mole fraction Flat vector with the m_nsp in the inner loop. length = ldx * ndim
ldfLeading dimension of the grad_V and current vectors.
grad_VThe electrostatic potential gradient.
currentThe electric current in A/m^2.

Reimplemented from Transport.

Definition at line 652 of file LiquidTransport.cpp.

References LiquidTransport::getSpeciesFluxesExt(), LiquidTransport::m_chargeSpecies, LiquidTransport::m_mw, Transport::m_nDim, Transport::m_nsp, LiquidTransport::set_Grad_T(), LiquidTransport::set_Grad_V(), and LiquidTransport::set_Grad_X().

void getSpeciesVdiff ( size_t  ndim,
const doublereal *  grad_T,
int  ldx,
const doublereal *  grad_X,
int  ldf,
doublereal *  Vdiff 
)
virtual

Get the species diffusive velocities wrt to the averaged velocity, given the gradients in mole fraction and temperature.

The average velocity can be computed on a mole-weighted or mass-weighted basis, or the diffusion velocities may be specified as relative to a specific species (i.e. a solvent) all according to the velocityBasis input parameter.

Units for the returned velocities are m s-1.

Parameters
ndimNumber of dimensions in the flux expressions
grad_TGradient of the temperature (length = ndim)
ldxLeading dimension of the grad_X array (usually equal to m_nsp but not always)
grad_XGradients of the mole fraction Flat vector with the m_nsp in the inner loop. length = ldx * ndim
ldfLeading dimension of the fluxes array (usually equal to m_nsp but not always)
VdiffOutput of the diffusive velocities. Flat vector with the m_nsp in the inner loop. length = ldx * ndim

Reimplemented from Transport.

Definition at line 678 of file LiquidTransport.cpp.

References LiquidTransport::getSpeciesVdiffExt(), LiquidTransport::set_Grad_T(), and LiquidTransport::set_Grad_X().

void getSpeciesVdiffES ( size_t  ndim,
const doublereal *  grad_T,
int  ldx,
const doublereal *  grad_X,
int  ldf,
const doublereal *  grad_Phi,
doublereal *  Vdiff 
)
virtual

Get the species diffusive velocities wrt to the averaged velocity, given the gradients in mole fraction, temperature and electrostatic potential.

The average velocity can be computed on a mole-weighted or mass-weighted basis, or the diffusion velocities may be specified as relative to a specific species (i.e. a solvent) all according to the velocityBasis input parameter.

Units for the returned velocities are m s-1.

Parameters
ndimNumber of dimensions in the flux expressions
grad_TGradient of the temperature (length = ndim)
ldxLeading dimension of the grad_X array (usually equal to m_nsp but not always)
grad_XGradients of the mole fraction Flat vector with the m_nsp in the inner loop. length = ldx * ndim
ldfLeading dimension of the fluxes array (usually equal to m_nsp but not always)
grad_PhiGradients of the electrostatic potential (length = ndim)
VdiffOutput of the species diffusion velocities Flat vector with the m_nsp in the inner loop. length = ldx * ndim

Reimplemented from Transport.

Definition at line 688 of file LiquidTransport.cpp.

References LiquidTransport::getSpeciesVdiffExt(), LiquidTransport::set_Grad_T(), LiquidTransport::set_Grad_V(), and LiquidTransport::set_Grad_X().

void getSpeciesFluxes ( size_t  ndim,
const doublereal *const  grad_T,
size_t  ldx,
const doublereal *const  grad_X,
size_t  ldf,
doublereal *const  fluxes 
)
virtual

Return the species diffusive mass fluxes wrt to the averaged velocity in [kmol/m^2/s].

The diffusive mass flux of species k [kmol/m^2/s] is computed using the Stefan-Maxwell equation

\[ X_i \nabla \mu_i = RT \sum_i \frac{X_i X_j}{D_{ij}} ( \vec{V}_j - \vec{V}_i ) \]

to determine the diffusion velocity and

\[ \vec{N}_i = C_T X_i \vec{V}_i \]

to determine the diffusion flux. Here \( C_T \) is the total concentration of the mixture [kmol/m^3], \( D_{ij} \) are the Stefa-Maxwell interaction parameters in [m^2/s], \( \vec{V}_{i} \) is the diffusion velocity of species i, \( \mu_i \) is the electrochemical potential of species i.

Note that for this method, there is no argument for the gradient of the electric potential (voltage). Electric potential gradients can be set with set_Grad_V() or method getSpeciesFluxesES() can be called.x

The diffusion velocity is relative to an average velocity that can be computed on a mole-weighted or mass-weighted basis, or the diffusion velocities may be specified as relative to a specific species (i.e. a solvent) all according to the velocityBasis input parameter.

Parameters
ndimThe number of spatial dimensions (1, 2, or 3).
grad_TThe temperature gradient (ignored in this model). (length = ndim)
ldxLeading dimension of the grad_X array. (usually equal to m_nsp but not always)
grad_XGradients of the mole fraction Flat vector with the m_nsp in the inner loop. length = ldx * ndim
ldfLeading dimension of the fluxes array (usually equal to m_nsp but not always)
fluxesOutput of the diffusive mass fluxes Flat vector with the m_nsp in the inner loop. length = ldx * ndim

Reimplemented from Transport.

Definition at line 702 of file LiquidTransport.cpp.

References LiquidTransport::getSpeciesFluxesExt(), LiquidTransport::set_Grad_T(), and LiquidTransport::set_Grad_X().

void getSpeciesFluxesES ( size_t  ndim,
const doublereal *  grad_T,
size_t  ldx,
const doublereal *  grad_X,
size_t  ldf,
const doublereal *  grad_Phi,
doublereal *  fluxes 
)
virtual

Return the species diffusive mass fluxes wrt to the averaged velocity in [kmol/m^2/s].

The diffusive mass flux of species k is computed using the Stefan-Maxwell equation

\[ X_i \nabla \mu_i = RT \sum_i \frac{X_i X_j}{D_{ij}} ( \vec{V}_j - \vec{V}_i ) \]

to determine the diffusion velocity and

\[ \vec{N}_i = C_T X_i \vec{V}_i \]

to determine the diffusion flux. Here \( C_T \) is the total concentration of the mixture [kmol/m^3], \( D_{ij} \) are the Stefa-Maxwell interaction parameters in [m^2/s], \( \vec{V}_{i} \) is the diffusion velocity of species i, \( \mu_i \) is the electrochemical potential of species i.

The diffusion velocity is relative to an average velocity that can be computed on a mole-weighted or mass-weighted basis, or the diffusion velocities may be specified as relative to a specific species (i.e. a solvent) all according to the velocityBasis input parameter.

Parameters
ndimThe number of spatial dimensions (1, 2, or 3).
grad_TThe temperature gradient (ignored in this model). (length = ndim)
ldxLeading dimension of the grad_X array. (usually equal to m_nsp but not always)
grad_XGradients of the mole fraction Flat vector with the m_nsp in the inner loop. length = ldx * ndim
ldfLeading dimension of the fluxes array (usually equal to m_nsp but not always)
grad_PhiGradients of the electrostatic potential length = ndim
fluxesOutput of the diffusive mass fluxes Flat vector with the m_nsp in the inner loop. length = ldx * ndim

Reimplemented from Transport.

Definition at line 712 of file LiquidTransport.cpp.

References LiquidTransport::getSpeciesFluxesExt(), LiquidTransport::set_Grad_T(), LiquidTransport::set_Grad_V(), and LiquidTransport::set_Grad_X().

void getSpeciesVdiffExt ( size_t  ldf,
doublereal *  Vdiff 
)
virtual

Return the species diffusive velocities relative to the averaged velocity.

This method acts similarly to getSpeciesVdiffES() but requires all gradients to be preset using methods set_Grad_X(), set_Grad_V(), set_Grad_T(). See the documentation of getSpeciesVdiffES() for details.

Parameters
ldfLeading dimension of the Vdiff array.
VdiffOutput of the diffusive velocities. Flat vector with the m_nsp in the inner loop. length = ldx * ndim

Definition at line 726 of file LiquidTransport.cpp.

References Transport::m_nDim, Transport::m_nsp, LiquidTransport::m_Vdiff, and LiquidTransport::stefan_maxwell_solve().

Referenced by LiquidTransport::getSpeciesVdiff(), and LiquidTransport::getSpeciesVdiffES().

void getSpeciesFluxesExt ( size_t  ldf,
doublereal *  fluxes 
)
virtual

Return the species diffusive fluxes relative to the averaged velocity.

This method acts similarly to getSpeciesFluxesES() but requires all gradients to be preset using methods set_Grad_X(), set_Grad_V(), set_Grad_T(). See the documentation of getSpeciesFluxesES() for details.

units = kg/m2/s

Parameters
ldfLeading dimension of the Vdiff array.
fluxesOutput of the diffusive fluxes. Flat vector with the m_nsp in the inner loop. length = ldx * ndim

Definition at line 737 of file LiquidTransport.cpp.

References LiquidTransport::m_flux, Transport::m_nDim, Transport::m_nsp, and LiquidTransport::stefan_maxwell_solve().

Referenced by LiquidTransport::getElectricConduct(), LiquidTransport::getElectricCurrent(), LiquidTransport::getSpeciesFluxes(), and LiquidTransport::getSpeciesFluxesES().

bool update_T ( void  )
protectedvirtual

Returns true if temperature has changed, in which case flags are set to recompute transport properties.

This is called whenever a transport property is requested. The first task is to check whether the temperature has changed since the last call to update_T(). If it hasn't then an immediate return is carried out.

Note this should be a lightweight function since it's part of all of the interfaces.

Returns
Returns true if the temperature has changed, and false otherwise

Definition at line 765 of file LiquidTransport.cpp.

References Cantera::fp2str(), LiquidTransport::m_diff_mix_ok, LiquidTransport::m_diff_temp_ok, LiquidTransport::m_ionCond_conc_ok, LiquidTransport::m_ionCond_mix_ok, LiquidTransport::m_ionCond_temp_ok, LiquidTransport::m_lambda_mix_ok, LiquidTransport::m_lambda_temp_ok, LiquidTransport::m_mobRat_conc_ok, LiquidTransport::m_mobRat_mix_ok, LiquidTransport::m_mobRat_temp_ok, LiquidTransport::m_radi_temp_ok, LiquidTransport::m_selfDiff_conc_ok, LiquidTransport::m_selfDiff_mix_ok, LiquidTransport::m_selfDiff_temp_ok, LiquidTransport::m_temp, Transport::m_thermo, LiquidTransport::m_visc_conc_ok, LiquidTransport::m_visc_mix_ok, LiquidTransport::m_visc_temp_ok, and Phase::temperature().

Referenced by LiquidTransport::getBinaryDiffCoeffs(), LiquidTransport::getSpeciesHydrodynamicRadius(), LiquidTransport::getSpeciesIonConductivity(), LiquidTransport::getSpeciesMobilityRatio(), LiquidTransport::getSpeciesSelfDiffusion(), LiquidTransport::getSpeciesViscosities(), LiquidTransport::ionConductivity(), LiquidTransport::mobilityRatio(), LiquidTransport::selfDiffusion(), LiquidTransport::stefan_maxwell_solve(), LiquidTransport::thermalConductivity(), and LiquidTransport::viscosity().

bool update_C ( )
protectedvirtual

Returns true if mixture composition has changed, in which case flags are set to recompute transport properties.

This is called for every interface call to check whether the concentrations have changed. Concentrations change whenever the pressure or the mole fraction has changed. If it has changed, the recalculations should be done.

Note this should be a lightweight function since it's part of all of the interfaces.

Returns
Returns true if the mixture composition has changed, and false otherwise.

Definition at line 809 of file LiquidTransport.cpp.

References LiquidTransport::concTot_, LiquidTransport::concTot_tran_, DATA_PTR, LiquidTransport::dens_, Phase::density(), Phase::getConcentrations(), Phase::getMassFractions(), Phase::getMoleFractions(), LiquidTransport::m_concentrations, LiquidTransport::m_diff_mix_ok, LiquidTransport::m_ionCond_conc_ok, LiquidTransport::m_ionCond_mix_ok, LiquidTransport::m_iStateMF, LiquidTransport::m_lambda_mix_ok, LiquidTransport::m_massfracs, LiquidTransport::m_massfracs_tran, LiquidTransport::m_mobRat_conc_ok, LiquidTransport::m_mobRat_mix_ok, LiquidTransport::m_molefracs, LiquidTransport::m_molefracs_tran, Transport::m_nsp, LiquidTransport::m_press, LiquidTransport::m_selfDiff_conc_ok, LiquidTransport::m_selfDiff_mix_ok, Transport::m_thermo, LiquidTransport::m_visc_conc_ok, LiquidTransport::m_visc_mix_ok, Phase::meanMolecularWeight(), LiquidTransport::meanMolecularWeight_, ThermoPhase::pressure(), Phase::stateMFNumber(), and Cantera::Tiny.

Referenced by LiquidTransport::ionConductivity(), LiquidTransport::mobilityRatio(), LiquidTransport::selfDiffusion(), LiquidTransport::stefan_maxwell_solve(), LiquidTransport::thermalConductivity(), and LiquidTransport::viscosity().

void update_Grad_lnAC ( )
protectedvirtual

Updates the internal value of the gradient of the logarithm of the activity, which is used in the gradient of the chemical potential.

Evaluate the gradients of the activity as they alter the diffusion coefficient.

The gradient of the chemical potential can be written in terms of gradient of the logarithm of the mole fraction times a correction associated with the gradient of the activity coefficient relative to that of the mole fraction. Specifically, the gradients of the logarithms of each are involved according to the formula

\[ \nabla \mu_k = RT \left[ \nabla ( \ln X_k ) + \nabla ( \ln \gamma_k ) \right] = RT \left[ \nabla ( \ln a_k ) \right] \]

The gradient in the activity coefficient requires the use of thermophase getdlnActCoeff that calculates its change based on a chane in the state (i.e. temperature and composition of each species) which was first implemented in MargulesVPSSTP.cpp (LiquidTransport.h doxygen)

Definition at line 951 of file LiquidTransport.cpp.

References ThermoPhase::getdlnActCoeffds(), LiquidTransport::m_Grad_lnAC, LiquidTransport::m_Grad_T, LiquidTransport::m_Grad_X, LiquidTransport::m_molefracs, Transport::m_nDim, Transport::m_nsp, and Transport::m_thermo.

Referenced by LiquidTransport::stefan_maxwell_solve().

void stefan_maxwell_solve ( )
protected

Solve the stefan_maxell equations for the diffusive fluxes.

The diffusive mass flux of species k is computed using the Stefan-Maxwell equation

\[ X_i \nabla \mu_i = RT \sum_i \frac{X_i X_j}{D_{ij}} ( \vec{V}_j - \vec{V}_i ) \]

to determine the diffusion velocity and

\[ \vec{N}_i = C_T X_i \vec{V}_i \]

to determine the diffusion flux. Here \( C_T \) is the total concentration of the mixture [kmol/m^3], \( D_{ij} \) are the Stefa-Maxwell interaction parameters in [m^2/s], \( \vec{V}_{i} \) is the diffusion velocity of species i, \( \mu_i \) is the electrochemical potential of species i.

The diffusion velocity is relative to an average velocity that can be computed on a mole-weighted or mass-weighted basis, or the diffusion velocities may be specified as relative to a specific species (i.e. a solvent) all according to the velocityBasis input parameter.

The gradient in the activity coefficient requires the use of thermophase getdlnActCoeff that calculates its change based on a change in the state i.e. temperature and composition of each species. First implemented in MargulesVPSSTP.cppmeter.

One of the Stefan Maxwell equations is replaced by the appropriate definition of the mass-averaged velocity, the mole-averaged velocity or the specification that velocities are relative to that of one species.

Definition at line 976 of file LiquidTransport.cpp.

References ThermoPhase::activityConvention(), Cantera::cAC_CONVENTION_MOLALITY, LiquidTransport::concTot_, DATA_PTR, Cantera::GasConstant, ThermoPhase::getActivityCoefficients(), LiquidTransport::m_A, LiquidTransport::m_actCoeff, LiquidTransport::m_B, LiquidTransport::m_bdiff, LiquidTransport::m_chargeSpecies, LiquidTransport::m_diff_temp_ok, LiquidTransport::m_flux, LiquidTransport::m_Grad_lnAC, LiquidTransport::m_Grad_mu, LiquidTransport::m_Grad_T, LiquidTransport::m_Grad_V, LiquidTransport::m_massfracs_tran, LiquidTransport::m_molefracs, LiquidTransport::m_molefracs_tran, Transport::m_nDim, Transport::m_nsp, Transport::m_thermo, LiquidTransport::m_Vdiff, Transport::m_velocityBasis, Phase::molarVolume(), Phase::molecularWeight(), Phase::molecularWeights(), DenseMatrix::resize(), Cantera::solve(), Phase::temperature(), LiquidTransport::update_C(), LiquidTransport::update_Grad_lnAC(), LiquidTransport::update_T(), LiquidTransport::updateDiff_T(), Cantera::VB_MASSAVG, and Cantera::VB_MOLEAVG.

Referenced by LiquidTransport::getMixDiffCoeffs(), LiquidTransport::getSpeciesFluxesExt(), and LiquidTransport::getSpeciesVdiffExt().

void updateViscosity_T ( )
protected

Updates the array of pure species viscosities internally.

The flag m_visc_ok is set to true.

Note that for viscosity, a positive activation energy corresponds to the typical case of a positive argument to the exponential so that the Arrhenius expression is

\[ \mu = A T^n \exp( + E / R T ) \]

Definition at line 882 of file LiquidTransport.cpp.

References Transport::m_nsp, LiquidTransport::m_visc_mix_ok, LiquidTransport::m_visc_temp_ok, LiquidTransport::m_viscSpecies, and LiquidTransport::m_viscTempDep_Ns.

Referenced by LiquidTransport::getSpeciesViscosities().

void updateIonConductivity_T ( )
protected

Update the temperature-dependent ionic conductivity terms for each species internally.

The flag m_ionCond_temp_ok is set to true.

Definition at line 896 of file LiquidTransport.cpp.

References LiquidTransport::m_ionCond_mix_ok, LiquidTransport::m_ionCond_temp_ok, LiquidTransport::m_ionCondSpecies, LiquidTransport::m_ionCondTempDep_Ns, and Transport::m_nsp.

Referenced by LiquidTransport::getSpeciesIonConductivity().

void updateMobilityRatio_T ( )
protected

Updates the array of pure species mobility ratios internally.

The flag m_mobRat_ok is set to true.

Definition at line 910 of file LiquidTransport.cpp.

References LiquidTransport::m_mobRat_mix_ok, LiquidTransport::m_mobRat_temp_ok, LiquidTransport::m_mobRatSpecies, LiquidTransport::m_mobRatTempDep_Ns, Transport::m_nsp, and LiquidTransport::m_nsp2.

Referenced by LiquidTransport::getSpeciesMobilityRatio().

void updateSelfDiffusion_T ( )
protected

Updates the array of pure species self diffusion coeffs internally.

The flag m_selfDiff_ok is set to true.

Definition at line 926 of file LiquidTransport.cpp.

References Transport::m_nsp, LiquidTransport::m_nsp2, LiquidTransport::m_selfDiff_mix_ok, LiquidTransport::m_selfDiff_temp_ok, LiquidTransport::m_selfDiffSpecies, and LiquidTransport::m_selfDiffTempDep_Ns.

Referenced by LiquidTransport::getSpeciesSelfDiffusion().

void updateHydrodynamicRadius_T ( )
protected

Update the temperature-dependent hydrodynamic radius terms for each species internally.

The flag m_radi_temp_ok is set to true.

Definition at line 942 of file LiquidTransport.cpp.

References LiquidTransport::m_hydrodynamic_radius, Transport::m_nsp, LiquidTransport::m_radi_mix_ok, LiquidTransport::m_radi_temp_ok, and LiquidTransport::m_radiusTempDep_Ns.

Referenced by LiquidTransport::getSpeciesHydrodynamicRadius().

void updateCond_T ( )
protected

Update the temperature-dependent parts of the mixture-averaged thermal conductivity internally.

Definition at line 861 of file LiquidTransport.cpp.

References LiquidTransport::m_lambda_mix_ok, LiquidTransport::m_lambda_temp_ok, LiquidTransport::m_lambdaSpecies, LiquidTransport::m_lambdaTempDep_Ns, and Transport::m_nsp.

void updateViscosities_C ( )
protected

Update the concentration parts of the viscosities.

Internal routine is run whenever the update_boolean m_visc_conc_ok is false. Currently there is no concentration dependence for the pure species viscosities.

Definition at line 877 of file LiquidTransport.cpp.

References LiquidTransport::m_visc_conc_ok.

void updateIonConductivity_C ( )
protected

Update the concentration parts of the ionic conductivity.

Internal routine is run whenever the update_boolean m_ionCond_conc_ok is false. Currently there is no concentration dependence for the pure species ionic conductivity.

Definition at line 891 of file LiquidTransport.cpp.

References LiquidTransport::m_ionCond_conc_ok.

void updateMobilityRatio_C ( )
protected

Update the concentration parts of the mobility ratio.

Internal routine is run whenever the update_boolean m_mobRat_conc_ok is false. Currently there is no concentration dependence for the pure species mobility ratio.

Definition at line 905 of file LiquidTransport.cpp.

References LiquidTransport::m_mobRat_conc_ok.

void updateSelfDiffusion_C ( )
protected

Update the concentration parts of the self diffusion.

Internal routine is run whenever the update_boolean m_selfDiff_conc_ok is false. Currently there is no concentration dependence for the pure species self diffusion.

Definition at line 921 of file LiquidTransport.cpp.

References LiquidTransport::m_selfDiff_conc_ok.

void updateHydrodynamicRadius_C ( )
protected

Update the concentration dependence of the hydrodynamic radius.

Internal routine is run whenever the update_boolean m_radi_conc_ok is false. Currently there is no concentration dependence for the hydrodynamic radius.

Definition at line 937 of file LiquidTransport.cpp.

References LiquidTransport::m_radi_conc_ok.

void updateDiff_T ( )
protected

Update the binary Stefan-Maxwell diffusion coefficients wrt T using calls to the appropriate LTPspecies subclass.

Definition at line 870 of file LiquidTransport.cpp.

References LiquidTransport::m_bdiff, LiquidTransport::m_diff_mix_ok, LiquidTransport::m_diff_temp_ok, and LiquidTransport::m_diffMixModel.

Referenced by LiquidTransport::getBinaryDiffCoeffs(), and LiquidTransport::stefan_maxwell_solve().

doublereal err ( const std::string &  msg) const
private

Throw an exception if this method is invoked.

This probably indicates something is not yet implemented.

Parameters
msgIndicates the member function which is not implemented

Definition at line 1221 of file LiquidTransport.cpp.

References Cantera::int2str(), and LiquidTransport::model().

Member Data Documentation

size_t m_nsp2
private
vector_fp m_mw
private

Local copy of the molecular weights of the species.

Length is equal to the number of species in the phase

Definition at line 803 of file LiquidTransport.h.

Referenced by LiquidTransport::getElectricConduct(), LiquidTransport::getElectricCurrent(), and LiquidTransport::initLiquid().

std::vector<LTPspecies*> m_viscTempDep_Ns
private

Viscosity for each species expressed as an appropriate subclass of LTPspecies.

These subclasses of LTPspecies evaluate the species-specific transport properties according to the parameters parsed in TransportFactory::getLiquidSpeciesTransportData().

Definition at line 812 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::updateViscosity_T(), and LiquidTransport::viscosity().

LiquidTranInteraction* m_viscMixModel
private

Viscosity of the mixture expressed as a subclass of LiquidTranInteraction.

These subclasses of LiquidTranInteraction evaluate the mixture transport properties according to the parameters parsed in TransportFactory::getLiquidInteractionsTransportData().

Definition at line 821 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::viscosity().

std::vector<LTPspecies*> m_ionCondTempDep_Ns
private

Ionic conductivity for each species expressed as an appropriate subclass of LTPspecies.

These subclasses of LTPspecies evaluate the species-specific transport properties according to the parameters parsed in TransportFactory::getLiquidSpeciesTransportData().

Definition at line 830 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::ionConductivity(), and LiquidTransport::updateIonConductivity_T().

LiquidTranInteraction* m_ionCondMixModel
private

Ionic Conductivity of the mixture expressed as a subclass of LiquidTranInteraction.

These subclasses of LiquidTranInteraction evaluate the mixture transport properties according to the parameters parsed in TransportFactory::getLiquidInteractionsTransportData().

Definition at line 839 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::ionConductivity().

std::vector<LTPvector> m_mobRatTempDep_Ns
private

Mobility ratio for the binary cominations of each species in each pure phase expressed as an appropriate subclass of LTPspecies.

These subclasses of LTPspecies evaluate the species-specific transport properties according to the parameters parsed in TransportFactory::getLiquidSpeciesTransportData().

Definition at line 851 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::mobilityRatio(), and LiquidTransport::updateMobilityRatio_T().

std::vector<LiquidTranInteraction*> m_mobRatMixModel
private

Mobility ratio for each binary combination of mobile species in the mixture expressed as a subclass of LiquidTranInteraction.

These subclasses of LiquidTranInteraction evaluate the mixture transport properties according to the parameters parsed in TransportFactory::getLiquidInteractionsTransportData().

Definition at line 860 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::mobilityRatio().

std::vector<LTPvector> m_selfDiffTempDep_Ns
private

Self Diffusion for each species in each pure species phase expressed as an appropriate subclass of LTPspecies.

These subclasses of LTPspecies evaluate the species-specific transport properties according to the parameters parsed in TransportFactory::getLiquidSpeciesTransportData().

Definition at line 869 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::selfDiffusion(), and LiquidTransport::updateSelfDiffusion_T().

std::vector<LiquidTranInteraction*> m_selfDiffMixModel
private

Self Diffusion for each species in the mixture expressed as a subclass of LiquidTranInteraction.

These subclasses of LiquidTranInteraction evaluate the mixture transport properties according to the parameters parsed in TransportFactory::getLiquidInteractionsTransportData().

Definition at line 878 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::selfDiffusion().

std::vector<LTPspecies*> m_lambdaTempDep_Ns
private

Thermal conductivity for each species expressed as an appropriate subclass of LTPspecies.

These subclasses of LTPspecies evaluate the species-specific transport properties according to the parameters parsed in TransportFactory::getLiquidSpeciesTransportData().

Definition at line 887 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::thermalConductivity(), and LiquidTransport::updateCond_T().

LiquidTranInteraction* m_lambdaMixModel
private

Thermal conductivity of the mixture expressed as a subclass of LiquidTranInteraction.

These subclasses of LiquidTranInteraction evaluate the mixture transport properties according to the parameters parsed in TransportFactory::getLiquidInteractionsTransportData().

Definition at line 896 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::thermalConductivity().

std::vector<LTPspecies*> m_diffTempDep_Ns
private

(NOT USED IN LiquidTransport.) Diffusion coefficient model for each species expressed as an appropriate subclass of LTPspecies

These subclasses of LTPspecies evaluate the species-specific transport properties according to the parameters parsed in TransportFactory::getLiquidSpeciesTransportData().

Since the LiquidTransport class uses the Stefan-Maxwell equation to describe species diffusivity, the species-specific diffusivity is irrelevant.

Definition at line 910 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid().

LiquidTranInteraction* m_diffMixModel
private

Species diffusivity of the mixture expressed as a subclass of LiquidTranInteraction.

This will return an array of Stefan-Maxwell interaction parameters for use in the Stefan-Maxwell solution.

These subclasses of LiquidTranInteraction evaluate the mixture transport properties according to the parameters parsed in TransportFactory::getLiquidInteractionsTransportData().

Definition at line 920 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::updateDiff_T().

DenseMatrix m_diff_Dij
private

Setfan-Maxwell diffusion coefficients.

Definition at line 923 of file LiquidTransport.h.

std::vector<LTPspecies*> m_radiusTempDep_Ns
private

Hydrodynamic radius for each species expressed as an appropriate subclass of LTPspecies.

These subclasses of LTPspecies evaluate the species-specific transport properties according to the parameters parsed in TransportFactory::getLiquidSpeciesTransportData(). length = nsp

Definition at line 932 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::updateHydrodynamicRadius_T().

LiquidTranInteraction* m_radiusMixModel
private

(Not used in LiquidTransport) Hydrodynamic radius of the mixture expressed as a subclass of LiquidTranInteraction

These subclasses of LiquidTranInteraction evaluate the mixture transport properties according to the parameters parsed in TransportFactory::getLiquidInteractionsTransportData().

Definition at line 942 of file LiquidTransport.h.

vector_fp m_hydrodynamic_radius
private
vector_fp m_Grad_X
private

Internal value of the gradient of the mole fraction vector.

Note, this is the only gradient value that can and perhaps should reflect the true state of the mole fractions in the application solution vector. In other words no cropping or massaging of the values to make sure they are above zero should occur. - developing ....

m_nsp is the number of species in the fluid k is the species index n is the dimensional index (x, y, or z). It has a length equal to m_nDim

m_Grad_X[n*m_nsp + k]

Definition at line 962 of file LiquidTransport.h.

Referenced by LiquidTransport::getMixDiffCoeffs(), LiquidTransport::initLiquid(), LiquidTransport::set_Grad_X(), and LiquidTransport::update_Grad_lnAC().

vector_fp m_Grad_lnAC
private

Gradient of the logarithm of the activity.

This quantity appears in the gradient of the chemical potential. It replaces the gradient of the mole fraction, and in this way serves to "modify" the diffusion coefficient.

\[ m\_Grad\_lnAC[k] = \nabla ( \ln X_k ) + \nabla ( \ln \gamma_k ) \]

k is the species index n is the dimensional index (x, y, or z). It has a length equal to m_nDim

m_Grad_X[n*m_nsp + k]

Definition at line 981 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::stefan_maxwell_solve(), and LiquidTransport::update_Grad_lnAC().

vector_fp m_Grad_T
private

Internal value of the gradient of the Temperature vector.

Generally, if a transport property needs this in its evaluation it will look to this place to get it.

No internal property is precalculated based on gradients. Gradients are assumed to be freshly updated before every property call.

Definition at line 991 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::set_Grad_T(), LiquidTransport::stefan_maxwell_solve(), and LiquidTransport::update_Grad_lnAC().

vector_fp m_Grad_P
private

Internal value of the gradient of the Pressure vector.

Generally, if a transport property needs this in its evaluation it will look to this place to get it.

No internal property is precalculated based on gradients. Gradients are assumed to be freshly updated before every property call.

Definition at line 1001 of file LiquidTransport.h.

vector_fp m_Grad_V
private

Internal value of the gradient of the Electric Voltage.

Generally, if a transport property needs this in its evaluation it will look to this place to get it.

No internal property is precalculated based on gradients. Gradients are assumed to be freshly updated before every property call.

Definition at line 1011 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::set_Grad_V(), and LiquidTransport::stefan_maxwell_solve().

vector_fp m_Grad_mu
private

Gradient of the electrochemical potential.

m_nsp is the number of species in the fluid. k is the species index. n is the dimensional index (x, y, or z)

\[ m\_Grad\_mu[n*m_nsp + k] \]

Definition at line 1022 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::stefan_maxwell_solve().

DenseMatrix m_bdiff
private

Array of Binary Diffusivities.

These are evaluated according to the subclass of LiquidTranInteraction stored in m_diffMixModel.

This has a size equal to nsp x nsp. It is a symmetric matrix. D_ii is the self diffusion coefficient. D_ii is not needed except for when there is one species in the mixture.

units m2/sec

Definition at line 1037 of file LiquidTransport.h.

Referenced by LiquidTransport::getBinaryDiffCoeffs(), LiquidTransport::initLiquid(), LiquidTransport::stefan_maxwell_solve(), and LiquidTransport::updateDiff_T().

vector_fp m_viscSpecies
private

Internal value of the species viscosities.

Viscosity of the species evaluated using subclass of LTPspecies held in m_viscTempDep_Ns.

Length = number of species

controlling update boolean -> m_visc_temp_ok

Definition at line 1048 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesViscosities(), LiquidTransport::initLiquid(), and LiquidTransport::updateViscosity_T().

vector_fp m_ionCondSpecies
private

Internal value of the species ionic conductivities.

Ionic conductivity of the species evaluated using subclass of LTPspecies held in m_ionCondTempDep_Ns.

Length = number of species

controlling update boolean -> m_ionCond_temp_ok

Definition at line 1059 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesIonConductivity(), LiquidTransport::initLiquid(), and LiquidTransport::updateIonConductivity_T().

DenseMatrix m_mobRatSpecies
private

Internal value of the species mobility ratios.

Mobility ratio of the species evaluated using subclass of LTPspecies held in m_mobRatTempDep_Ns.

Length = number of species

controlling update boolean -> m_mobRat_temp_ok

Definition at line 1070 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesMobilityRatio(), LiquidTransport::initLiquid(), and LiquidTransport::updateMobilityRatio_T().

DenseMatrix m_selfDiffSpecies
private

Internal value of the species self diffusion coefficients.

Self diffusion of the species evaluated using subclass of LTPspecies held in m_selfDiffTempDep_Ns.

Length = number of species

controlling update boolean -> m_selfDiff_temp_ok

Definition at line 1081 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesSelfDiffusion(), LiquidTransport::initLiquid(), and LiquidTransport::updateSelfDiffusion_T().

vector_fp m_lambdaSpecies
private

Internal value of the species individual thermal conductivities.

Thermal conductivities of the species evaluated using subclass of LTPspecies held in m_lambdaTempDep_Ns.

Length = number of species

Definition at line 1091 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::updateCond_T().

int m_iStateMF
private

State of the mole fraction vector.

Definition at line 1094 of file LiquidTransport.h.

Referenced by LiquidTransport::update_C().

vector_fp m_massfracs
private

Local copy of the mass fractions of the species in the phase.

The mass fraction vector comes from the ThermoPhase object.

length = m_nsp

Definition at line 1102 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::update_C().

vector_fp m_massfracs_tran
private

Local copy of the mass fractions of the species in the phase.

This version of the mass fraction vector is adjusted to a minimum lower bound of Tiny for use in transport calculations.

Definition at line 1109 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::stefan_maxwell_solve(), and LiquidTransport::update_C().

vector_fp m_molefracs
private

Local copy of the mole fractions of the species in the phase.

The mole fractions here are assumed to be bounded by 0.0 and 1.0 and they are assumed to add up to one exactly. This mole fraction vector comes from the ThermoPhase object. Derivative quantities from this are referred to as bounded.

Update info? length = m_nsp

Definition at line 1121 of file LiquidTransport.h.

Referenced by LiquidTransport::getMixDiffCoeffs(), LiquidTransport::initLiquid(), LiquidTransport::stefan_maxwell_solve(), LiquidTransport::update_C(), and LiquidTransport::update_Grad_lnAC().

vector_fp m_molefracs_tran
private

Non-zero mole fraction vector used in transport property calculations.

The mole fractions here are assumed to be bounded by Tiny and 1.0 and they may not be assumed to add up to one. This mole fraction vector is created from the ThermoPhase object. Derivative quantities of this use the _tran suffix.

Update info? length = m_nsp

Definition at line 1133 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::stefan_maxwell_solve(), and LiquidTransport::update_C().

vector_fp m_concentrations
private

Local copy of the concentrations of the species in the phase.

The concentrations are consistent with the m_molefracs vector which is bounded and sums to one.

Update info? length = m_nsp

Definition at line 1143 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::update_C().

doublereal concTot_
private

Local copy of the total concentration.

This is consistent with the m_concentrations[] and m_molefracs[] vector.

Definition at line 1150 of file LiquidTransport.h.

Referenced by LiquidTransport::stefan_maxwell_solve(), and LiquidTransport::update_C().

doublereal concTot_tran_
private

Local copy of the total concentration.

This is consistent with the x_molefracs_tran vector and with the concTot_ number;

Definition at line 1157 of file LiquidTransport.h.

Referenced by LiquidTransport::update_C().

doublereal meanMolecularWeight_
private

Mean molecular mass.

Definition at line 1160 of file LiquidTransport.h.

Referenced by LiquidTransport::update_C().

doublereal dens_
private

Density.

Definition at line 1163 of file LiquidTransport.h.

Referenced by LiquidTransport::update_C().

vector_fp m_chargeSpecies
private

Local copy of the charge of each species.

Contains the charge of each species (length m_nsp)

Definition at line 1169 of file LiquidTransport.h.

Referenced by LiquidTransport::getElectricConduct(), LiquidTransport::getElectricCurrent(), LiquidTransport::initLiquid(), and LiquidTransport::stefan_maxwell_solve().

vector_fp m_volume_spec
private

Specific volume for each species. Local copy from thermo object.

Definition at line 1172 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid().

vector_fp m_actCoeff
private

Vector of activity coefficients.

Definition at line 1175 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::stefan_maxwell_solve().

DenseMatrix m_B
private

RHS to the stefan-maxwell equation.

Definition at line 1178 of file LiquidTransport.h.

Referenced by LiquidTransport::stefan_maxwell_solve().

DenseMatrix m_A
private

Matrix for the stefan maxwell equation.

Definition at line 1181 of file LiquidTransport.h.

Referenced by LiquidTransport::stefan_maxwell_solve().

doublereal m_temp
private

Current Temperature -> locally stored.

This is used to test whether new temperature computations should be performed.

Definition at line 1188 of file LiquidTransport.h.

Referenced by LiquidTransport::getFluidMobilities(), LiquidTransport::getMobilities(), and LiquidTransport::update_T().

doublereal m_press
private

Current value of the pressure.

Definition at line 1191 of file LiquidTransport.h.

Referenced by LiquidTransport::update_C().

Array2D m_flux
private

Solution of the Stefan Maxwell equation in terms of flux.

This is the mass flux of species k in units of kg m-3 s-1.

Definition at line 1197 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesFluxesExt(), LiquidTransport::initLiquid(), and LiquidTransport::stefan_maxwell_solve().

Array2D m_Vdiff
private

Solution of the Stefan Maxwell equation.

This is the diffusion velocity of species k in units of m/s and relative to the mole-averaged velocity.

Definition at line 1204 of file LiquidTransport.h.

Referenced by LiquidTransport::getMixDiffCoeffs(), LiquidTransport::getSpeciesVdiffExt(), LiquidTransport::initLiquid(), and LiquidTransport::stefan_maxwell_solve().

doublereal m_lambda
private

Saved value of the mixture thermal conductivity.

Definition at line 1207 of file LiquidTransport.h.

Referenced by LiquidTransport::thermalConductivity().

doublereal m_viscmix
private

Saved value of the mixture viscosity.

Definition at line 1210 of file LiquidTransport.h.

Referenced by LiquidTransport::viscosity().

doublereal m_ionCondmix
private

Saved value of the mixture ionic conductivity.

Definition at line 1213 of file LiquidTransport.h.

Referenced by LiquidTransport::ionConductivity().

vector_fp m_mobRatMix
private

Saved values of the mixture mobility ratios.

Definition at line 1216 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::mobilityRatio().

vector_fp m_selfDiffMix
private

Saved values of the mixture self diffusion coefficients.

Definition at line 1219 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::selfDiffusion().

vector_fp m_spwork
private

work space

Length is equal to m_nsp

Definition at line 1225 of file LiquidTransport.h.

Referenced by LiquidTransport::getFluidMobilities(), LiquidTransport::getMobilities(), and LiquidTransport::initLiquid().

bool m_visc_mix_ok
private

Boolean indicating that the top-level mixture viscosity is current.

This is turned false for every change in T, P, or C.

Definition at line 1232 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::update_C(), LiquidTransport::update_T(), LiquidTransport::updateViscosity_T(), and LiquidTransport::viscosity().

bool m_visc_temp_ok
private

Boolean indicating that weight factors wrt viscosity is current.

Definition at line 1235 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesViscosities(), LiquidTransport::initLiquid(), LiquidTransport::update_T(), and LiquidTransport::updateViscosity_T().

bool m_visc_conc_ok
private

Flag to indicate that the pure species viscosities are current wrt the concentration.

Definition at line 1239 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::update_C(), LiquidTransport::update_T(), and LiquidTransport::updateViscosities_C().

bool m_ionCond_mix_ok
private

Boolean indicating that the top-level mixture ionic conductivity is current.

This is turned false for every change in T, P, or C.

Definition at line 1245 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::ionConductivity(), LiquidTransport::update_C(), LiquidTransport::update_T(), and LiquidTransport::updateIonConductivity_T().

bool m_ionCond_temp_ok
private

Boolean indicating that weight factors wrt ionic conductivty is current.

Definition at line 1248 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesIonConductivity(), LiquidTransport::initLiquid(), LiquidTransport::update_T(), and LiquidTransport::updateIonConductivity_T().

bool m_ionCond_conc_ok
private

Flag to indicate that the pure species ionic conductivities are current wrt the concentration.

Definition at line 1252 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::update_C(), LiquidTransport::update_T(), and LiquidTransport::updateIonConductivity_C().

bool m_cond_mix_ok
private

Flag to indicate that the mixture conductivity is current.

Definition at line 1255 of file LiquidTransport.h.

Referenced by LiquidTransport::thermalConductivity().

bool m_mobRat_mix_ok
private

Boolean indicating that the top-level mixture mobility ratio is current.

This is turned false for every change in T, P, or C.

Definition at line 1261 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::mobilityRatio(), LiquidTransport::update_C(), LiquidTransport::update_T(), and LiquidTransport::updateMobilityRatio_T().

bool m_mobRat_temp_ok
private

Boolean indicating that weight factors wrt mobility ratio is current.

Definition at line 1264 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesMobilityRatio(), LiquidTransport::initLiquid(), LiquidTransport::update_T(), and LiquidTransport::updateMobilityRatio_T().

bool m_mobRat_conc_ok
private

Flag to indicate that the pure species mobility ratios are current wrt the concentration.

Definition at line 1268 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::update_C(), LiquidTransport::update_T(), and LiquidTransport::updateMobilityRatio_C().

bool m_selfDiff_mix_ok
private

Boolean indicating that the top-level mixture self diffusion is current.

This is turned false for every change in T, P, or C.

Definition at line 1274 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::selfDiffusion(), LiquidTransport::update_C(), LiquidTransport::update_T(), and LiquidTransport::updateSelfDiffusion_T().

bool m_selfDiff_temp_ok
private

Boolean indicating that weight factors wrt self diffusion is current.

Definition at line 1277 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesSelfDiffusion(), LiquidTransport::initLiquid(), LiquidTransport::update_T(), and LiquidTransport::updateSelfDiffusion_T().

bool m_selfDiff_conc_ok
private

Flag to indicate that the pure species self diffusion are current wrt the concentration.

Definition at line 1281 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::update_C(), LiquidTransport::update_T(), and LiquidTransport::updateSelfDiffusion_C().

bool m_radi_mix_ok
private

Boolean indicating that mixture diffusion coeffs are current.

Definition at line 1284 of file LiquidTransport.h.

Referenced by LiquidTransport::updateHydrodynamicRadius_T().

bool m_radi_temp_ok
private

Boolean indicating that temperature dependence of hydrodynamic radius is current.

Definition at line 1288 of file LiquidTransport.h.

Referenced by LiquidTransport::getSpeciesHydrodynamicRadius(), LiquidTransport::initLiquid(), LiquidTransport::update_T(), and LiquidTransport::updateHydrodynamicRadius_T().

bool m_radi_conc_ok
private

Flag to indicate that the hydrodynamic radius is current is current wrt the concentration.

Definition at line 1292 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), and LiquidTransport::updateHydrodynamicRadius_C().

bool m_diff_mix_ok
private

Boolean indicating that mixture diffusion coeffs are current.

Definition at line 1295 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::update_C(), LiquidTransport::update_T(), and LiquidTransport::updateDiff_T().

bool m_diff_temp_ok
private
bool m_lambda_temp_ok
private

Flag to indicate that the pure species conductivities are current wrt the temperature.

Definition at line 1302 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::update_T(), and LiquidTransport::updateCond_T().

bool m_lambda_mix_ok
private

Boolean indicating that mixture conductivity is current.

Definition at line 1305 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid(), LiquidTransport::thermalConductivity(), LiquidTransport::update_C(), LiquidTransport::update_T(), and LiquidTransport::updateCond_T().

int m_mode
private

Mode indicator for transport models – currently unused.

Definition at line 1308 of file LiquidTransport.h.

Referenced by LiquidTransport::initLiquid().

bool m_debug
private

Debugging flags.

Turn on to get debugging information

Definition at line 1314 of file LiquidTransport.h.


The documentation for this class was generated from the following files: