Cantera  3.1.0a1
Transport Class Reference

Base class for transport property managers. More...

#include <Transport.h>

Inheritance diagram for Transport:
[legend]

Detailed Description

Base class for transport property managers.

All classes that compute transport properties for a single phase derive from this class. Class Transport is meant to be used as a base class only. It is possible to instantiate it, but its methods throw exceptions if called.

Relationship of the Transport class to the ThermoPhase Class

This section describes how calculations are carried out within the Transport class. The Transport class and derived classes of the the Transport class necessarily use the ThermoPhase class to obtain the list of species and the thermodynamic state of the phase.

No state information is stored within Transport classes. Queries to the underlying ThermoPhase object must be made to obtain the state of the system.

An exception to this however is the state information concerning the the gradients of variables. This information is not stored within the ThermoPhase objects. It may be collected within the Transport objects. In fact, the meaning of const operations within the Transport class refers to calculations which do not change the state of the system nor the state of the first order gradients of the system.

When a const operation is evoked within the Transport class, it is also implicitly assumed that the underlying state within the ThermoPhase object has not changed its values.

Todo:
Provide a general mechanism to store the gradients of state variables within the system.

Definition at line 71 of file Transport.h.

Public Member Functions

 Transport ()=default
 Constructor. More...
 
 Transport (const Transport &)=delete
 
Transportoperator= (const Transport &)=delete
 
virtual string transportModel () const
 Identifies the model represented by this Transport object. More...
 
ThermoPhasethermo ()
 Phase object. More...
 
void checkSpeciesIndex (size_t k) const
 Check that the specified species index is in range. More...
 
void checkSpeciesArraySize (size_t kk) const
 Check that an array size is at least nSpecies(). More...
 
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. More...
 
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. More...
 
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. More...
 
virtual void getThermalDiffCoeffs (double *const dt)
 Return a vector of Thermal diffusion coefficients [kg/m/sec]. More...
 
virtual void getBinaryDiffCoeffs (const size_t ld, double *const d)
 Returns the matrix of binary diffusion coefficients [m^2/s]. More...
 
virtual void getMultiDiffCoeffs (const size_t ld, double *const d)
 Return the Multicomponent diffusion coefficients. Units: [m^2/s]. More...
 
virtual void getMixDiffCoeffs (double *const d)
 Returns a vector of mixture averaged diffusion coefficients. More...
 
virtual void getMixDiffCoeffsMole (double *const d)
 Returns a vector of mixture averaged diffusion coefficients. More...
 
virtual void getMixDiffCoeffsMass (double *const d)
 Returns a vector of mixture averaged diffusion coefficients. More...
 
virtual void getViscosityPolynomial (size_t i, double *coeffs) const
 Return the polynomial fits to the viscosity of species i. More...
 
virtual void getConductivityPolynomial (size_t i, double *coeffs) const
 Return the temperature fits of the heat conductivity of species i. More...
 
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) More...
 
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) More...
 
virtual void setViscosityPolynomial (size_t i, double *coeffs)
 Modify the polynomial fits to the viscosity of species i. More...
 
virtual void setConductivityPolynomial (size_t i, double *coeffs)
 Modify the temperature fits of the heat conductivity of species i. More...
 
virtual void setBinDiffusivityPolynomial (size_t i, size_t j, double *coeffs)
 Modify the polynomial fits to the binary diffusivity of species pair (i, j) More...
 
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) More...
 
AnyMap parameters () const
 Return the parameters for a phase definition which are needed to reconstruct an identical object using the newTransport function. More...
 
Transport Properties
virtual double viscosity ()
 The viscosity in Pa-s. More...
 
virtual void getSpeciesViscosities (double *const visc)
 Returns the pure species viscosities. More...
 
virtual double bulkViscosity ()
 The bulk viscosity in Pa-s. More...
 
virtual double thermalConductivity ()
 Returns the mixture thermal conductivity in W/m/K. More...
 
virtual double electricalConductivity ()
 The electrical conductivity (Siemens/m). More...
 
virtual void getMobilities (double *const mobil_e)
 Get the Electrical mobilities (m^2/V/s). More...
 
Transport manager construction

These methods are used during construction.

virtual void init (ThermoPhase *thermo, int mode=0, int log_level=0)
 Initialize a transport manager. More...
 
virtual bool CKMode () const
 Boolean indicating the form of the transport properties polynomial fits. More...
 

Protected Attributes

ThermoPhasem_thermo
 pointer to the object representing the phase More...
 
size_t m_nsp = 0
 Number of species. More...
 

Constructor & Destructor Documentation

◆ Transport()

Transport ( )
default

Constructor.

New transport managers should be created using TransportFactory, not by calling the constructor directly.

See also
TransportFactory

Member Function Documentation

◆ transportModel()

virtual string transportModel ( ) const
inlinevirtual

Identifies the model represented by this Transport object.

Each derived class should override this method to return a meaningful identifier.

Since
New in Cantera 3.0. The name returned by this method corresponds to the canonical name used in the YAML input format.

Reimplemented in WaterTransport, UnityLewisTransport, MultiTransport, MixTransport, IonGasTransport, HighPressureGasTransport, and DustyGasTransport.

Definition at line 93 of file Transport.h.

◆ thermo()

ThermoPhase& thermo ( )
inline

Phase object.

Every transport manager is designed to compute properties for a specific phase of a mixture, which might be a liquid solution, a gas mixture, a surface, etc. This method returns a reference to the object representing the phase itself.

Definition at line 103 of file Transport.h.

◆ checkSpeciesIndex()

void checkSpeciesIndex ( size_t  k) const

Check that the specified species index is in range.

Throws an exception if k is greater than nSpecies()

Definition at line 16 of file Transport.cpp.

◆ checkSpeciesArraySize()

void checkSpeciesArraySize ( size_t  kk) const

Check that an array size is at least nSpecies().

Throws an exception if kk is less than nSpecies(). Used before calls which take an array pointer.

Definition at line 23 of file Transport.cpp.

◆ viscosity()

virtual double viscosity ( )
inlinevirtual

The viscosity in Pa-s.

Reimplemented in WaterTransport, IonGasTransport, HighPressureGasTransport, and GasTransport.

Definition at line 122 of file Transport.h.

◆ getSpeciesViscosities()

virtual void getSpeciesViscosities ( double *const  visc)
inlinevirtual

Returns the pure species viscosities.

The units are Pa-s and the length is the number of species

Parameters
viscVector of viscosities

Reimplemented in GasTransport.

Definition at line 133 of file Transport.h.

◆ bulkViscosity()

virtual double bulkViscosity ( )
inlinevirtual

The bulk viscosity in Pa-s.

The bulk viscosity is only non-zero in rare cases. Most transport managers either overload this method to return zero, or do not implement it, in which case an exception is thrown if called.

Reimplemented in WaterTransport.

Definition at line 144 of file Transport.h.

◆ thermalConductivity()

virtual double thermalConductivity ( )
inlinevirtual

Returns the mixture thermal conductivity in W/m/K.

Units are in W / m K or equivalently kg m / s3 K

Returns
thermal conductivity in W/m/K.

Reimplemented in WaterTransport, MultiTransport, MixTransport, IonGasTransport, and HighPressureGasTransport.

Definition at line 155 of file Transport.h.

◆ electricalConductivity()

virtual double electricalConductivity ( )
inlinevirtual

The electrical conductivity (Siemens/m).

Reimplemented in IonGasTransport.

Definition at line 161 of file Transport.h.

◆ getMobilities()

virtual void getMobilities ( double *const  mobil_e)
inlinevirtual

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.

Frequently, but not always, the mobility is calculated from the diffusion coefficient using the Einstein relation

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

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

Reimplemented in MixTransport, and IonGasTransport.

Definition at line 182 of file Transport.h.

◆ getSpeciesFluxes()

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 
)
inlinevirtual

Get the species diffusive mass fluxes wrt to the specified solution averaged velocity, given the gradients in mole fraction and temperature.

Units for the returned fluxes are kg m-2 s-1.

Usually the specified solution average velocity is the mass averaged velocity. This is changed in some subclasses, however.

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)
fluxesOutput of the diffusive mass fluxes. Flat vector with the m_nsp in the inner loop. length = ldx * ndim

Reimplemented in MultiTransport, and MixTransport.

Definition at line 208 of file Transport.h.

◆ getMolarFluxes()

virtual void getMolarFluxes ( const double *const  state1,
const double *const  state2,
const double  delta,
double *const  cfluxes 
)
inlinevirtual

Get the molar fluxes [kmol/m^2/s], given the thermodynamic state at two nearby points.

Parameters
[in]state1Array of temperature, density, and mass fractions for state 1.
[in]state2Array of temperature, density, and mass fractions for state 2.
[in]deltaDistance from state 1 to state 2 (m).
[out]cfluxesOutput array containing the diffusive molar fluxes of species from state1 to state2. This is a flat vector with m_nsp in the inner loop. length = ldx * ndim. Units are [kmol/m^2/s].

Reimplemented in MultiTransport, and DustyGasTransport.

Definition at line 228 of file Transport.h.

◆ getMassFluxes()

virtual void getMassFluxes ( const double *  state1,
const double *  state2,
double  delta,
double *  mfluxes 
)
inlinevirtual

Get the mass fluxes [kg/m^2/s], given the thermodynamic state at two nearby points.

Parameters
[in]state1Array of temperature, density, and mass fractions for state 1.
[in]state2Array of temperature, density, and mass fractions for state 2.
[in]deltaDistance from state 1 to state 2 (m).
[out]mfluxesOutput array containing the diffusive mass fluxes of species from state1 to state2. This is a flat vector with m_nsp in the inner loop. length = ldx * ndim. Units are [kg/m^2/s].

Reimplemented in MultiTransport.

Definition at line 248 of file Transport.h.

◆ getThermalDiffCoeffs()

virtual void getThermalDiffCoeffs ( double *const  dt)
inlinevirtual

Return a vector of Thermal diffusion coefficients [kg/m/sec].

The thermal diffusion coefficient \( D^T_k \) is defined so that the diffusive mass flux of species k induced by the local temperature gradient is given by the following formula:

\[ M_k J_k = -D^T_k \nabla \ln T. \]

The thermal diffusion coefficient can be either positive or negative.

Parameters
dtOn return, dt will contain the species thermal diffusion coefficients. Dimension dt at least as large as the number of species. Units are kg/m/s.

Reimplemented in MultiTransport, MixTransport, and HighPressureGasTransport.

Definition at line 271 of file Transport.h.

◆ getBinaryDiffCoeffs()

virtual void getBinaryDiffCoeffs ( const size_t  ld,
double *const  d 
)
inlinevirtual

Returns the matrix of binary diffusion coefficients [m^2/s].

Parameters
[in]ldInner stride for writing the two dimension diffusion coefficients into a one dimensional vector
[out]dDiffusion coefficient matrix (must be at least m_k * m_k in length.

Reimplemented in HighPressureGasTransport, and GasTransport.

Definition at line 283 of file Transport.h.

◆ getMultiDiffCoeffs()

virtual void getMultiDiffCoeffs ( const size_t  ld,
double *const  d 
)
inlinevirtual

Return the Multicomponent diffusion coefficients. Units: [m^2/s].

If the transport manager implements a multicomponent diffusion model, then this method returns the array of multicomponent diffusion coefficients. Otherwise it throws an exception.

Parameters
[in]ldThe dimension of the inner loop of d (usually equal to m_nsp)
[out]dflat vector of diffusion coefficients, fortran ordering. d[ld*j+i] is the D_ij diffusion coefficient (the diffusion coefficient for species i due to concentration gradients in species j). Units: m^2/s

Reimplemented in MultiTransport, HighPressureGasTransport, and DustyGasTransport.

Definition at line 300 of file Transport.h.

◆ getMixDiffCoeffs()

virtual void getMixDiffCoeffs ( double *const  d)
inlinevirtual

Returns a vector of mixture averaged diffusion coefficients.

Mixture-averaged diffusion coefficients [m^2/s]. If the transport manager implements a mixture-averaged diffusion model, then this method returns the array of mixture-averaged diffusion coefficients. Otherwise it throws an exception.

Parameters
dReturn vector of mixture averaged diffusion coefficients Units = m2/s. Length = n_sp

Reimplemented in UnityLewisTransport, IonGasTransport, and GasTransport.

Definition at line 315 of file Transport.h.

◆ getMixDiffCoeffsMole()

virtual void getMixDiffCoeffsMole ( double *const  d)
inlinevirtual

Returns a vector of mixture averaged diffusion coefficients.

Reimplemented in UnityLewisTransport, and GasTransport.

Definition at line 321 of file Transport.h.

◆ getMixDiffCoeffsMass()

virtual void getMixDiffCoeffsMass ( double *const  d)
inlinevirtual

Returns a vector of mixture averaged diffusion coefficients.

Reimplemented in UnityLewisTransport, and GasTransport.

Definition at line 327 of file Transport.h.

◆ getViscosityPolynomial()

virtual void getViscosityPolynomial ( size_t  i,
double *  coeffs 
) const
inlinevirtual

Return the polynomial fits to the viscosity of species i.

Reimplemented in GasTransport.

Definition at line 333 of file Transport.h.

◆ getConductivityPolynomial()

virtual void getConductivityPolynomial ( size_t  i,
double *  coeffs 
) const
inlinevirtual

Return the temperature fits of the heat conductivity of species i.

Reimplemented in GasTransport.

Definition at line 339 of file Transport.h.

◆ getBinDiffusivityPolynomial()

virtual void getBinDiffusivityPolynomial ( size_t  i,
size_t  j,
double *  coeffs 
) const
inlinevirtual

Return the polynomial fits to the binary diffusivity of species pair (i, j)

Reimplemented in GasTransport.

Definition at line 345 of file Transport.h.

◆ getCollisionIntegralPolynomial()

virtual void getCollisionIntegralPolynomial ( size_t  i,
size_t  j,
double *  astar_coeffs,
double *  bstar_coeffs,
double *  cstar_coeffs 
) const
inlinevirtual

Return the polynomial fits to the collision integral of species pair (i, j)

Reimplemented in GasTransport.

Definition at line 351 of file Transport.h.

◆ setViscosityPolynomial()

virtual void setViscosityPolynomial ( size_t  i,
double *  coeffs 
)
inlinevirtual

Modify the polynomial fits to the viscosity of species i.

Reimplemented in GasTransport.

Definition at line 360 of file Transport.h.

◆ setConductivityPolynomial()

virtual void setConductivityPolynomial ( size_t  i,
double *  coeffs 
)
inlinevirtual

Modify the temperature fits of the heat conductivity of species i.

Reimplemented in GasTransport.

Definition at line 366 of file Transport.h.

◆ setBinDiffusivityPolynomial()

virtual void setBinDiffusivityPolynomial ( size_t  i,
size_t  j,
double *  coeffs 
)
inlinevirtual

Modify the polynomial fits to the binary diffusivity of species pair (i, j)

Reimplemented in GasTransport.

Definition at line 372 of file Transport.h.

◆ setCollisionIntegralPolynomial()

virtual void setCollisionIntegralPolynomial ( size_t  i,
size_t  j,
double *  astar_coeffs,
double *  bstar_coeffs,
double *  cstar_coeffs,
bool  flag 
)
inlinevirtual

Modify the polynomial fits to the collision integral of species pair (i, j)

Reimplemented in GasTransport.

Definition at line 378 of file Transport.h.

◆ parameters()

AnyMap parameters ( ) const

Return the parameters for a phase definition which are needed to reconstruct an identical object using the newTransport function.

This excludes the individual species transport properties, which are handled separately.

Definition at line 30 of file Transport.cpp.

◆ init()

virtual void init ( ThermoPhase thermo,
int  mode = 0,
int  log_level = 0 
)
inlinevirtual

Initialize a transport manager.

This routine sets up a transport manager. It calculates the collision integrals and populates species-dependent data structures.

Parameters
thermoPointer to the ThermoPhase object
modeChemkin compatible mode or not. This alters the specification of the collision integrals. defaults to no.
log_levelDefaults to zero, no logging

Reimplemented in WaterTransport, MultiTransport, MixTransport, GasTransport, and IonGasTransport.

Definition at line 407 of file Transport.h.

◆ CKMode()

virtual bool CKMode ( ) const
inlinevirtual

Boolean indicating the form of the transport properties polynomial fits.

Returns true if the Chemkin form is used.

Reimplemented in GasTransport.

Definition at line 411 of file Transport.h.

Member Data Documentation

◆ m_thermo

ThermoPhase* m_thermo
protected

pointer to the object representing the phase

Definition at line 420 of file Transport.h.

◆ m_nsp

size_t m_nsp = 0
protected

Number of species.

Definition at line 423 of file Transport.h.


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