Objects Representing Phases¶
Composite Phase Objects¶
These classes are composite representations of a substance which has thermodynamic, chemical kinetic, and (optionally) transport properties.
-
class
cantera.
Solution
(infile='', phaseid='', source=None, thermo=None, species=(), kinetics=None, reactions=(), **kwargs)¶ Bases:
cantera._cantera.ThermoPhase
,cantera._cantera.Kinetics
,cantera._cantera.Transport
A class for chemically-reacting solutions. Instances can be created to represent any type of solution – a mixture of gases, a liquid solution, or a solid solution, for example.
Class
Solution
derives from classesThermoPhase
,Kinetics
, andTransport
. It defines no methods of its own, and is provided so that a single object can be used to compute thermodynamic, kinetic, and transport properties of a solution.To skip initialization of the Transport object, pass the keyword argument
transport_model=None
to theSolution
constructor.The most common way to instantiate
Solution
objects is by using a phase definition, species and reactions defined in an input file:gas = ct.Solution('gri30.yaml')
If an input file defines multiple phases, the corresponding key in the phases map (in YAML), name (in CTI), or id (in XML) can be used to specify the desired phase via the
name
keyword argument of the constructor:gas = ct.Solution('diamond.yaml', name='gas') diamond = ct.Solution('diamond.yaml', name='diamond')
The name of the
Solution
object defaults to the phase identifier specified in the input file. Upon initialization of a ‘Solution’ object, a custom name can assigned via:gas.name = 'my_custom_name'
Solution
objects can also be constructed usingSpecies
andReaction
objects which can themselves either be imported from input files or defined directly in Python:spec = ct.Species.listFromFile('gri30.yaml') rxns = ct.Reaction.listFromFile('gri30.yaml') gas = ct.Solution(thermo='IdealGas', kinetics='GasKinetics', species=spec, reactions=rxns, name='my_custom_name')
where the
thermo
andkinetics
keyword arguments are strings specifying the thermodynamic and kinetics model, respectively, andspecies
andreactions
keyword arguments are lists ofSpecies
andReaction
objects, respectively.Types of underlying models that form the composite
Solution
object are queried using thethermo_model
,kinetics_model
andtransport_model
attributes; further, thecomposite
attribute is a shorthand returning a tuple containing the types of the three constitutive models.For non-trivial uses cases of this functionality, see the examples extract_submechanism.py and mechanism_reduction.py.
In addition,
Solution
objects can be constructed by passing the text of the CTI or XML phase definition in directly, using thesource
keyword argument:cti_def = ''' ideal_gas(name='gas', elements='O H Ar', species='gri30: all', reactions='gri30: all', options=['skip_undeclared_elements', 'skip_undeclared_species', 'skip_undeclared_third_bodies'], initial_state=state(temperature=300, pressure=101325))''' gas = ct.Solution(source=cti_def)
-
class
cantera.
Interface
(infile='', phaseid='', phases=(), thermo=None, species=(), kinetics=None, reactions=())¶ Bases:
cantera._cantera.InterfacePhase
,cantera._cantera.InterfaceKinetics
Two-dimensional interfaces.
Instances of class
Interface
represent reacting 2D interfaces between bulk 3D phases. ClassInterface
defines no methods of its own. All of its methods derive from eitherInterfacePhase
orInterfaceKinetics
.To construct an
Interface
object, adjacent bulk phases which participate in reactions need to be created and then passed in as a list in theadjacent
argument to the constructor:gas = ct.Solution('diamond.yaml', name='gas') diamond = ct.Solution('diamond.yaml', name='diamond') diamond_surf = ct.Interface('diamond.yaml', name='diamond_100', adjacent=[gas, diamond])
-
class
cantera.
DustyGas
(infile, phaseid='')¶ Bases:
cantera._cantera.ThermoPhase
,cantera._cantera.Kinetics
,cantera._cantera.DustyGasTransport
A composite class which models a gas in a stationary, solid, porous medium.
The only transport properties computed are the multicomponent diffusion coefficients. The model does not compute viscosity or thermal conductivity.
Pure Fluid Phases¶
The following convenience functions can be used to create PureFluid
objects
with the indicated equation of state:
-
cantera.
CarbonDioxide
()¶ Create a
PureFluid
object using the equation of state for carbon dioxide.The object returned by this method implements an accurate equation of state for carbon dioxide that can be used in the liquid, vapor, saturated liquid/vapor, and supercritical regions of the phase diagram. The equation of state is taken from
W. C. Reynolds, Thermodynamic Properties in SI: graphs, tables, and computational equations for forty substances. Stanford: Stanford University, 1979. Print.
For more details, see classes :ct:PureFluid and tpx::CarbonDioxide in the Cantera C++ source code documentation.
-
cantera.
Heptane
()¶ Create a
PureFluid
object using the equation of state for heptane.The object returned by this method implements an accurate equation of state for n-heptane that can be used in the liquid, vapor, saturated liquid/vapor, and supercritical regions of the phase diagram. The equation of state is taken from
W. C. Reynolds, Thermodynamic Properties in SI: graphs, tables, and computational equations for forty substances. Stanford: Stanford University, 1979. Print.
For more details, see classes :ct:PureFluid and tpx::Heptane in the Cantera C++ source code documentation.
-
cantera.
Hfc134a
()¶ Create a
PureFluid
object using the equation of state for HFC-134a.The object returned by this method implements an accurate equation of state for refrigerant HFC134a (R134a) that can be used in the liquid, vapor, saturated liquid/vapor, and supercritical regions of the phase diagram. Implements the equation of state given in:
R. Tillner-Roth, H. D. Baehr. An International Standard Formulation for The Thermodynamic Properties of 1,1,1,2-Tetrafluoroethane (HFC-134a) for Temperatures From 170 K to 455 K and Pressures up to 70 MPa. J. Phys. Chem. Ref. Data, Vol. 23, No. 5, 1994. pp. 657–729. http://dx.doi.org/10.1063/1.555958
For more details, see classes :ct:PureFluid and tpx::HFC134a in the Cantera C++ source code documentation.
-
cantera.
Hydrogen
()¶ Create a
PureFluid
object using the equation of state for hydrogen.The object returned by this method implements an accurate equation of state for hydrogen that can be used in the liquid, vapor, saturated liquid/vapor, and supercritical regions of the phase diagram. The equation of state is taken from
W. C. Reynolds, Thermodynamic Properties in SI: graphs, tables, and computational equations for forty substances Stanford: Stanford University, 1979. Print.
For more details, see classes :ct:PureFluid and tpx::hydrogen in the Cantera C++ source code documentation.
-
cantera.
Methane
()¶ Create a
PureFluid
object using the equation of state for methane.The object returned by this method implements an accurate equation of state for methane that can be used in the liquid, vapor, saturated liquid/vapor, and supercritical regions of the phase diagram. The equation of state is taken from
W. C. Reynolds, Thermodynamic Properties in SI: graphs, tables, and computational equations for forty substances Stanford: Stanford University, 1979. Print.
For more details, see classes :ct:PureFluid and tpx::methane in the Cantera C++ source code documentation.
-
cantera.
Nitrogen
()¶ Create a
PureFluid
object using the equation of state for nitrogen.The object returned by this method implements an accurate equation of state for nitrogen that can be used in the liquid, vapor, saturated liquid/vapor, and supercritical regions of the phase diagram. The equation of state is taken from
W. C. Reynolds, Thermodynamic Properties in SI: graphs, tables, and computational equations for forty substances Stanford: Stanford University, 1979. Print.
For more details, see classes :ct:PureFluid and tpx::nitrogen in the Cantera C++ source code documentation.
-
cantera.
Oxygen
()¶ Create a
PureFluid
object using the equation of state for oxygen.The object returned by this method implements an accurate equation of state for oxygen that can be used in the liquid, vapor, saturated liquid/vapor, and supercritical regions of the phase diagram. The equation of state is taken from
W. C. Reynolds, Thermodynamic Properties in SI: graphs, tables, and computational equations for forty substances Stanford: Stanford University, 1979. Print.
For more details, see classes :ct:PureFluid and tpx::oxygen in the Cantera C++ source code documentation.
-
cantera.
Water
(backend='Reynolds')¶ Create a
PureFluid
object using the equation of state for water and theWaterTransport
class for viscosity and thermal conductivity.The object returned by this method implements an accurate equation of state for water, where implementations are selected using the backend switch.
For the
Reynolds
backend, the equation of state is taken fromW. C. Reynolds, Thermodynamic Properties in SI: graphs, tables, and computational equations for forty substances. Stanford: Stanford University, 1979. Print.
which can be used in the liquid, vapor, saturated liquid/vapor, and supercritical regions of the phase diagram.
The
IAPWS95
backend implements an IAPWS (International Association for the Properties of Water and Steam) formulation for thermodynamic properties taken fromW. Wagner, A. Pruss, The IAPWS Formulation 1995 for the Thermodynamic Properties of Ordinary Water Substance for General and Scientific Use, J. Phys. Chem. Ref. Dat, 31, 387, 2002.
which currently only implements liquid and supercritical regions.
In both cases, formulas for transport are taken from
J. V. Sengers, J. T. R. Watson, Improved International Formulations for the Viscosity and Thermal Conductivity of Water Substance, J. Phys. Chem. Ref. Data, 15, 1291, 1986.
For more details, see classes :ct:PureFluid, tpx::water, :ct:WaterSSTP and :ct:WaterTransport in the Cantera C++ source code documentation.
Representing Quantities of Phases¶
-
class
cantera.
Quantity
(phase, mass=None, moles=None, constant='UV')¶ Bases:
object
A class representing a specific quantity of a
Solution
. In addition to the properties which can be computed for classSolution
, classQuantity
provides several additional capabilities. AQuantity
object is created from aSolution
with either the mass or number of moles specified:>>> gas = ct.Solution('gri30.xml') >>> gas.TPX = 300, 5e5, 'O2:1.0, N2:3.76' >>> q1 = ct.Quantity(gas, mass=5) # 5 kg of air
The state of a
Quantity
can be changed in the same way as aSolution
:>>> q1.TP = 500, 101325
Quantities have properties which provide access to extensive properties:
>>> q1.volume 7.1105094 >>> q1.enthalpy 1032237.84
The size of a
Quantity
can be changed by setting the mass or number of moles:>>> q1.moles = 3 >>> q1.mass 86.552196 >>> q1.volume 123.086
or by multiplication:
>>> q1 *= 2 >>> q1.moles 6.0
Finally, Quantities can be added, providing an easy way of calculating the state resulting from mixing two substances:
>>> q1.mass = 5 >>> q2 = ct.Quantity(gas) >>> q2.TPX = 300, 101325, 'CH4:1.0' >>> q2.mass = 1 >>> q3 = q1 + q2 # combine at constant UV >>> q3.T 432.31234 >>> q3.P 97974.9871 >>> q3.mole_fraction_dict() {'CH4': 0.26452900448117395, 'N2': 0.5809602821745349, 'O2': 0.1545107133442912}
If a different property pair should be held constant when combining, this can be specified as follows:
>>> q1.constant = q2.constant = 'HP' >>> q3 = q1 + q2 # combine at constant HP >>> q3.T 436.03320 >>> q3.P 101325.0
-
property
DP
¶ Get/Set density [kg/m^3] and pressure [Pa].
-
property
DPX
¶ Get/Set density [kg/m^3], pressure [Pa], and mole fractions.
-
property
DPY
¶ Get/Set density [kg/m^3], pressure [Pa], and mass fractions.
-
property
HP
¶ Get/Set enthalpy [J/kg or J/kmol] and pressure [Pa].
-
property
HPX
¶ Get/Set enthalpy [J/kg or J/kmol], pressure [Pa] and mole fractions.
-
property
HPY
¶ Get/Set enthalpy [J/kg or J/kmol], pressure [Pa] and mass fractions.
-
property
ID
¶ The identifier of the object. The default value corresponds to the CTI/XML/YAML input file phase entry.
Deprecated since version 2.5: To be deprecated with version 2.5, and removed thereafter. Usage merged with
name
.
-
property
P
¶ Pressure [Pa].
-
property
P_sat
¶ Saturation pressure [Pa] at the current temperature.
-
property
SP
¶ Get/Set entropy [J/kg/K or J/kmol/K] and pressure [Pa].
-
property
SPX
¶ Get/Set entropy [J/kg/K or J/kmol/K], pressure [Pa], and mole fractions.
-
property
SPY
¶ Get/Set entropy [J/kg/K or J/kmol/K], pressure [Pa], and mass fractions.
-
property
SV
¶ Get/Set entropy [J/kg/K or J/kmol/K] and specific volume [m^3/kg or m^3/kmol].
-
property
SVX
¶ Get/Set entropy [J/kg/K or J/kmol/K], specific volume [m^3/kg or m^3/kmol], and mole fractions.
-
property
SVY
¶ Get/Set entropy [J/kg/K or J/kmol/K], specific volume [m^3/kg or m^3/kmol], and mass fractions.
-
property
T
¶ Temperature [K].
-
property
TD
¶ Get/Set temperature [K] and density [kg/m^3 or kmol/m^3].
-
property
TDX
¶ Get/Set temperature [K], density [kg/m^3 or kmol/m^3], and mole fractions.
-
property
TDY
¶ Get/Set temperature [K] and density [kg/m^3 or kmol/m^3], and mass fractions.
-
property
TP
¶ Get/Set temperature [K] and pressure [Pa].
-
property
TPX
¶ Get/Set temperature [K], pressure [Pa], and mole fractions.
-
property
TPY
¶ Get/Set temperature [K], pressure [Pa], and mass fractions.
-
property
T_sat
¶ Saturation temperature [K] at the current pressure.
-
property
UV
¶ Get/Set internal energy [J/kg or J/kmol] and specific volume [m^3/kg or m^3/kmol].
-
property
UVX
¶ Get/Set internal energy [J/kg or J/kmol], specific volume [m^3/kg or m^3/kmol], and mole fractions.
-
property
UVY
¶ Get/Set internal energy [J/kg or J/kmol], specific volume [m^3/kg or m^3/kmol], and mass fractions.
-
property
X
¶ Get/Set the species mole fractions. Can be set as an array, as a dictionary, or as a string. Always returns an array:
>>> phase.X = [0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5] >>> phase.X = {'H2':0.1, 'O2':0.4, 'AR':0.5} >>> phase.X = 'H2:0.1, O2:0.4, AR:0.5' >>> phase.X array([0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5])
-
property
Y
¶ Get/Set the species mass fractions. Can be set as an array, as a dictionary, or as a string. Always returns an array:
>>> phase.Y = [0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5] >>> phase.Y = {'H2':0.1, 'O2':0.4, 'AR':0.5} >>> phase.Y = 'H2:0.1, O2:0.4, AR:0.5' >>> phase.Y array([0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5])
-
property
activities
¶ Array of nondimensional activities. Returns either molar or molal activities depending on the convention of the thermodynamic model.
-
property
activity_coefficients
¶ Array of nondimensional, molar activity coefficients.
-
property
add_reaction
¶ Add a new reaction to this phase.
-
property
add_species
¶ Add a new species to this phase. Missing elements will be added automatically.
-
property
add_species_alias
¶ Add the alternate species name alias for an original species name.
-
property
atomic_weight
¶ Atomic weight [kg/kmol] of element m
-
property
atomic_weights
¶ Array of atomic weight [kg/kmol] for each element in the mixture.
-
property
basis
¶ Determines whether intensive thermodynamic properties are treated on a
mass
(per kg) ormolar
(per kmol) basis. This affects the values returned by the propertiesh
,u
,s
,g
,v
,density
,cv
, andcp
, as well as the values used with the state-setting properties such asHPX
andUV
.
-
property
binary_diff_coeffs
¶ Binary diffusion coefficients [m^2/s].
-
property
case_sensitive_species_names
¶ Enforce case-sensitivity for look up of species names
-
property
charges
¶ Array of species charges [elem. charge].
-
property
chemical_potentials
¶ Array of species chemical potentials [J/kmol].
-
property
composite
¶ Returns tuple of thermo/kinetics/transport models associated with this SolutionBase object.
-
property
concentrations
¶ Get/Set the species concentrations [kmol/m^3].
-
property
cp_mass
¶ Specific heat capacity at constant pressure [J/kg/K].
-
property
cp_mole
¶ Molar heat capacity at constant pressure [J/kmol/K].
-
property
creation_rates
¶ Creation rates for each species. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
critical_pressure
¶ Critical pressure [Pa].
-
property
critical_temperature
¶ Critical temperature [K].
-
property
cv_mass
¶ Specific heat capacity at constant volume [J/kg/K].
-
property
cv_mole
¶ Molar heat capacity at constant volume [J/kmol/K].
-
property
delta_enthalpy
¶ Change in enthalpy for each reaction [J/kmol].
-
property
delta_entropy
¶ Change in entropy for each reaction [J/kmol/K].
-
property
delta_gibbs
¶ Change in Gibbs free energy for each reaction [J/kmol].
-
property
delta_standard_enthalpy
¶ Change in standard-state enthalpy (independent of composition) for each reaction [J/kmol].
-
property
delta_standard_entropy
¶ Change in standard-state entropy (independent of composition) for each reaction [J/kmol/K].
-
property
delta_standard_gibbs
¶ Change in standard-state Gibbs free energy (independent of composition) for each reaction [J/kmol].
-
property
density_mass
¶ (Mass) density [kg/m^3].
-
property
density_mole
¶ Molar density [kmol/m^3].
-
property
destruction_rates
¶ Destruction rates for each species. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
electric_potential
¶ Get/Set the electric potential [V] for this phase.
-
property
electrical_conductivity
¶ Electrical conductivity. [S/m].
-
property
electrochemical_potentials
¶ Array of species electrochemical potentials [J/kmol].
-
property
element_index
¶ The index of element element, which may be specified as a string or an integer. In the latter case, the index is checked for validity and returned. If no such element is present, an exception is thrown.
-
property
element_name
¶ Name of the element with index m.
-
property
element_names
¶ A list of all the element names.
-
property
elemental_mass_fraction
¶ Get the elemental mass fraction \(Z_{\mathrm{mass},m}\) of element \(m\) as defined by:
\[Z_{\mathrm{mass},m} = \sum_k \frac{a_{m,k} M_m}{M_k} Y_k \]with \(a_{m,k}\) being the number of atoms of element \(m\) in species \(k\), \(M_m\) the atomic weight of element \(m\), \(M_k\) the molecular weight of species \(k\), and \(Y_k\) the mass fraction of species \(k\):
>>> phase.elemental_mass_fraction('H') 1.0
- Parameters
m – Base element, may be specified by name or by index.
-
property
elemental_mole_fraction
¶ Get the elemental mole fraction \(Z_{\mathrm{mole},m}\) of element \(m\) (the number of atoms of element m divided by the total number of atoms) as defined by:
\[Z_{\mathrm{mole},m} = \frac{\sum_k a_{m,k} X_k} {\sum_k \sum_j a_{j,k} X_k} \]with \(a_{m,k}\) being the number of atoms of element \(m\) in species \(k\), \(\sum_j\) being a sum over all elements, and \(X_k\) being the mole fraction of species \(k\):
>>> phase.elemental_mole_fraction('H') 1.0
- Parameters
m – Base element, may be specified by name or by index.
-
property
enthalpy_mass
¶ Specific enthalpy [J/kg].
-
property
enthalpy_mole
¶ Molar enthalpy [J/kmol].
-
property
entropy_mass
¶ Specific entropy [J/kg/K].
-
property
entropy_mole
¶ Molar entropy [J/kmol/K].
-
equilibrate
(XY=None, *args, **kwargs)¶ Set the state to equilibrium. By default, the property pair
self.constant
is held constant. SeeThermoPhase.equilibrate
.
-
property
equilibrium_constants
¶ Equilibrium constants in concentration units for all reactions.
-
property
equivalence_ratio
¶ Get the equivalence ratio of the current mixture, which is a conserved quantity. Considers the oxidation of C to CO2, H to H2O and S to SO2. Other elements are assumed not to participate in oxidation (that is, N ends up as N2). If fuel and oxidizer are not specified, the equivalence ratio is computed from the available oxygen and the required oxygen for complete oxidation. The
basis
determines the composition of fuel and oxidizer:basis='mole'
(default) means mole fractions,basis='mass'
means mass fractions:>>> gas.set_equivalence_ratio(0.5, 'CH3:0.5, CH3OH:.5, N2:0.125', 'O2:0.21, N2:0.79, NO:0.01') >>> gas.equivalence_ratio('CH3:0.5, CH3OH:.5, N2:0.125', 'O2:0.21, N2:0.79, NO:0.01') 0.5
- Parameters
fuel – Fuel species name or mole/mass fractions as string, array, or dict.
oxidizer – Oxidizer species name or mole/mass fractions as a string, array, or dict.
basis – Determines if
fuel
andoxidizer
are given in mole fractions (basis='mole'
) or mass fractions (basis='mass'
)
-
property
find_isomers
¶ Find species/isomers matching a composition specified by comp.
-
property
forward_rate_constants
¶ Forward rate constants for all reactions. The computed values include all temperature-dependent, pressure-dependent, and third body contributions. Units are a combination of kmol, m^3 and s, that depend on the rate expression for the reaction.
-
property
forward_rates_of_progress
¶ Forward rates of progress for the reactions. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
get_equivalence_ratio
¶ Get the composition of a fuel/oxidizer mixture. This gives the equivalence ratio of an unburned mixture. This is not a quantity that is conserved after oxidation. Considers the oxidation of C to CO2, H to H2O and S to SO2. Other elements are assumed not to participate in oxidation (that is, N ends up as N2).
- Parameters
oxidizers – List of oxidizer species names as strings. Default: with
oxidizers=[]
, every species that contains O but does not contain H, C, or S is considered to be an oxidizer.ignore –
List of species names as strings to ignore.
>>> gas.set_equivalence_ratio(0.5, 'CH3:0.5, CH3OH:.5, N2:0.125', 'O2:0.21, N2:0.79, NO:0.01') >>> gas.get_equivalence_ratio() 0.5 >>> gas.get_equivalence_ratio(['O2']) # Only consider O2 as the oxidizer instead of O2 and NO 0.488095238095 >>> gas.X = 'CH4:1, O2:2, NO:0.1' >>> gas.get_equivalence_ratio(ignore=['NO']) 1.0
Deprecated since version 2.5: To be deprecated with version 2.5, and removed thereafter. Replaced by function
equivalence_ratio
.
-
property
gibbs_mass
¶ Specific Gibbs free energy [J/kg].
-
property
gibbs_mole
¶ Molar Gibbs free energy [J/kmol].
-
property
has_phase_transition
¶ Returns true if the phase represents a substance with phase transitions
-
property
heat_production_rates
¶ Get the volumetric heat production rates [W/m^3] on a per-reaction basis. The sum over all reactions results in the total volumetric heat release rate. Example: C. K. Law: Combustion Physics (2006), Fig. 7.8.6
>>> gas.heat_production_rates[1] # heat production rate of the 2nd reaction
-
property
heat_release_rate
¶ Get the total volumetric heat release rate [W/m^3].
-
property
int_energy_mass
¶ Specific internal energy [J/kg].
-
property
int_energy_mole
¶ Molar internal energy [J/kmol].
-
property
is_compressible
¶ Returns true if the density of the phase is an independent variable defining the thermodynamic state of a substance
-
property
is_pure
¶ Returns true if the phase represents a pure (fixed composition) substance
-
property
is_reversible
¶ True if reaction
i_reaction
is reversible.
-
property
isothermal_compressibility
¶ Isothermal compressibility [1/Pa].
-
property
kinetics_model
¶ Return type of kinetics.
-
property
kinetics_species_index
¶ The index of species species of phase phase within arrays returned by methods of class
Kinetics
. If species is a string, the phase argument is unused.
-
property
kinetics_species_name
¶ Name of the species with index k in the arrays returned by methods of class
Kinetics
.
-
property
kinetics_species_names
¶ A list of all species names, corresponding to the arrays returned by methods of class
Kinetics
.
-
property
mass_fraction_dict
¶
-
property
max_temp
¶ Maximum temperature for which the thermodynamic data for the phase are valid.
-
property
mean_molecular_weight
¶ The mean molecular weight (molar mass) [kg/kmol].
-
property
min_temp
¶ Minimum temperature for which the thermodynamic data for the phase are valid.
-
property
mix_diff_coeffs
¶ Mixture-averaged diffusion coefficients [m^2/s] relating the mass-averaged diffusive fluxes (with respect to the mass averaged velocity) to gradients in the species mole fractions.
-
property
mix_diff_coeffs_mass
¶ Mixture-averaged diffusion coefficients [m^2/s] relating the diffusive mass fluxes to gradients in the species mass fractions.
-
property
mix_diff_coeffs_mole
¶ Mixture-averaged diffusion coefficients [m^2/s] relating the molar diffusive fluxes to gradients in the species mole fractions.
-
property
mixture_fraction
¶ Get the mixture fraction of the current mixture (kg fuel / (kg oxidizer + kg fuel)) This is a quantity that is conserved after oxidation. Considers the oxidation of C to CO2, H to H2O and S to SO2. Other elements are assumed not to participate in oxidation (that is, N ends up as N2). The
basis
determines the composition of fuel and oxidizer:basis="mole"
(default) means mole fractions,basis="mass"
means mass fractions. The mixture fraction can be computed from a single element (for example, carbon withelement="C"
) or from all elements, which is the Bilger mixture fraction (element="Bilger"
). The Bilger mixture fraction is computed by default:>>> gas.set_mixture_fraction(0.5, 'CH3:0.5, CH3OH:.5, N2:0.125', 'O2:0.21, N2:0.79, NO:0.01') >>> gas.mixture_fraction('CH3:0.5, CH3OH:.5, N2:0.125', 'O2:0.21, N2:0.79, NO:0.01') 0.5
- Parameters
fuel – Fuel species name or mole/mass fractions as string, array, or dict.
oxidizer – Oxidizer species name or mole/mass fractions as a string, array, or dict.
basis – Determines if
fuel
andoxidizer
are given in mole fractions (basis='mole'
) or mass fractions (basis='mass'
)element – Computes the mixture fraction from the specified elemental mass fraction (given by element name or element index) or as the Bilger mixture fraction (default)
-
property
mobilities
¶ Electrical mobilities of charged species [m^2/s-V]
-
property
modify_reaction
¶ Modify the
Reaction
with indexirxn
to have the same rate parameters asrxn
.rxn
must have the same reactants and products and be of the same type (i.e.ElementaryReaction
,FalloffReaction
,PlogReaction
, etc.) as the existing reaction. This method does not modify the third-body efficiencies, reaction orders, or reversibility of the reaction.
-
property
modify_species
¶
-
property
mole_fraction_dict
¶
-
property
molecular_weights
¶ Array of species molecular weights (molar masses) [kg/kmol].
-
property
multi_diff_coeffs
¶ Multicomponent diffusion coefficients [m^2/s].
-
property
multiplier
¶ A scaling factor applied to the rate coefficient for reaction i_reaction. Can be used to carry out sensitivity analysis or to selectively disable a particular reaction. See
set_multiplier
.
-
property
n_atoms
¶ Number of atoms of element element in species species. The element and species may be specified by name or by index.
>>> phase.n_atoms('CH4','H') 4
-
property
n_elements
¶ Number of elements.
-
property
n_phases
¶ Number of phases in the reaction mechanism.
-
property
n_reactions
¶ Number of reactions in the reaction mechanism.
-
property
n_selected_species
¶ Number of species selected for output (by slicing of Solution object)
-
property
n_species
¶ Number of species.
-
property
n_total_species
¶ Total number of species in all phases participating in the kinetics mechanism.
-
property
name
¶ The name assigned to this object. The default value corresponds to the CTI/XML/YAML input file phase entry.
-
property
net_production_rates
¶ Net production rates for each species. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
net_rates_of_progress
¶ Net rates of progress for the reactions. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
partial_molar_cp
¶ Array of species partial molar specific heat capacities at constant pressure [J/kmol/K].
-
property
partial_molar_enthalpies
¶ Array of species partial molar enthalpies [J/kmol].
-
property
partial_molar_entropies
¶ Array of species partial molar entropies [J/kmol/K].
-
property
partial_molar_int_energies
¶ Array of species partial molar internal energies [J/kmol].
-
property
partial_molar_volumes
¶ Array of species partial molar volumes [m^3/kmol].
-
property
phase
¶ Get the underlying
Solution
object, with the state set to match the wrappingQuantity
object.
-
property
phase_of_matter
¶ Get the thermodynamic phase (gas, liquid, etc.) at the current conditions.
-
property
product_stoich_coeff
¶ The stoichiometric coefficient of species k_spec as a product in reaction i_reaction.
-
property
product_stoich_coeffs
¶ The array of product stoichiometric coefficients. Element [k,i] of this array is the product stoichiometric coefficient of species k in reaction i.
-
property
products
¶ The products portion of the reaction equation
-
property
reactant_stoich_coeff
¶ The stoichiometric coefficient of species k_spec as a reactant in reaction i_reaction.
-
property
reactant_stoich_coeffs
¶ The array of reactant stoichiometric coefficients. Element [k,i] of this array is the reactant stoichiometric coefficient of species k in reaction i.
-
property
reactants
¶ The reactants portion of the reaction equation
-
property
reaction
¶ Return a
Reaction
object representing the reaction with indexi_reaction
. Changes to this object do not affect theKinetics
orSolution
object until themodify_reaction
function is called.
-
property
reaction_equation
¶ The equation for the specified reaction. See also
reaction_equations
.
-
property
reaction_equations
¶ Returns a list containing the reaction equation for all reactions in the mechanism (if indices is unspecified) or the equations for each reaction in the sequence indices. For example:
>>> gas.reaction_equations() ['2 O + M <=> O2 + M', 'O + H + M <=> OH + M', 'O + H2 <=> H + OH', ...] >>> gas.reaction_equations([2,3]) ['O + H + M <=> OH + M', 'O + H2 <=> H + OH']
See also
reaction_equation
.
-
property
reaction_phase_index
¶ The index of the phase where the reactions occur.
-
property
reaction_type
¶ Type of reaction i_reaction.
-
property
reactions
¶ Return a list of all
Reaction
objects. Changes to these objects do not affect theKinetics
orSolution
object until themodify_reaction
function is called.
-
property
reference_pressure
¶ Reference state pressure [Pa].
-
property
report
¶ Generate a report describing the thermodynamic state of this phase. To print the report to the terminal, simply call the phase object. The following two statements are equivalent:
>>> phase() >>> print(phase.report())
-
property
reverse_rate_constants
¶ Reverse rate constants for all reactions. The computed values include all temperature-dependent, pressure-dependent, and third body contributions. Units are a combination of kmol, m^3 and s, that depend on the rate expression for the reaction.
-
property
reverse_rates_of_progress
¶ Reverse rates of progress for the reactions. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
selected_species
¶
-
property
set_equivalence_ratio
¶ Set the composition to a mixture of
fuel
andoxidizer
at the specified equivalence ratiophi
, holding temperature and pressure constant. Considers the oxidation of C to CO2, H to H2O and S to SO2. Other elements are assumed not to participate in oxidation (that is, N ends up as N2). Thebasis
determines the fuel and oxidizer compositions:basis='mole'
means mole fractions (default),basis='mass'
means mass fractions:>>> gas.set_equivalence_ratio(0.5, 'CH4', 'O2:1.0, N2:3.76', basis='mole') >>> gas.mass_fraction_dict() {'CH4': 0.02837633052851681, 'N2': 0.7452356312613029, 'O2': 0.22638803821018036} >>> gas.set_equivalence_ratio(1.2, {'NH3':0.8, 'CO':0.2}, 'O2:1.0', basis='mole') >>> gas.mass_fraction_dict() {'CO': 0.14784006249290754, 'NH3': 0.35956645545401045, 'O2': 0.49259348205308207}
- Parameters
phi – Equivalence ratio
fuel – Fuel species name or mole/mass fractions as string, array, or dict.
oxidizer – Oxidizer species name or mole/mass fractions as a string, array, or dict.
basis – Determines if
fuel
andoxidizer
are given in mole fractions (basis='mole'
) or mass fractions (basis='mass'
)
-
property
set_mixture_fraction
¶ Set the composition to a mixture of
fuel
andoxidizer
at the specified mixture fraction mixture_fraction (kg fuel / kg mixture), holding temperature and pressure constant. Considers the oxidation of C to CO2, H to H2O and S to SO2. Other elements are assumed not to participate in oxidation (that is, N ends up as N2). Thebasis
determines the composition of fuel and oxidizer:basis='mole'
(default) means mole fractions,basis='mass'
means mass fractions:>>> gas.set_mixture_fraction(0.5, 'CH4', 'O2:1.0, N2:3.76') >>> gas.mass_fraction_dict() {'CH4': 0.5, 'N2': 0.38350014242997776, 'O2': 0.11649985757002226} >>> gas.set_mixture_fraction(0.5, {'NH3':0.8, 'CO':0.2}, 'O2:1.0') >>> gas.mass_fraction_dict() {'CO': 0.145682068778996, 'NH3': 0.354317931221004, 'O2': 0.5}
- Parameters
mixture_fraction – Mixture fraction (kg fuel / kg mixture)
fuel – Fuel species name or mole/mass fractions as string, array, or dict.
oxidizer – Oxidizer species name or mole/mass fractions as a string, array, or dict.
basis – determines if
fuel
andoxidizer
are given in mole fractions (basis='mole'
) or mass fractions (basis='mass'
)
-
property
set_multiplier
¶ Set the multiplier for for reaction i_reaction to value. If i_reaction is not specified, then the multiplier for all reactions is set to value. See
multiplier
.
-
property
set_unnormalized_mass_fractions
¶ Set the mass fractions without normalizing to force sum(Y) == 1.0. Useful primarily when calculating derivatives with respect to Y[k] by finite difference.
-
property
set_unnormalized_mole_fractions
¶ Set the mole fractions without normalizing to force sum(X) == 1.0. Useful primarily when calculating derivatives with respect to X[k] by finite difference.
-
property
source
¶ The source of this object (such as a file name).
-
property
species
¶ Return the
Species
object for species k, where k is either the species index or the species name. If k is not specified, a list of all species objects is returned. Changes to this object do not affect theThermoPhase
orSolution
object until themodify_species
function is called.
-
property
species_index
¶ The index of species species, which may be specified as a string or an integer. In the latter case, the index is checked for validity and returned. If no such species is present, an exception is thrown.
-
property
species_name
¶ Name of the species with index k.
-
property
species_names
¶ A list of all the species names.
-
property
species_viscosities
¶ Pure species viscosities [Pa-s]
-
property
standard_cp_R
¶ Array of nondimensional species standard-state specific heat capacities at constant pressure at the current temperature and pressure.
-
property
standard_enthalpies_RT
¶ Array of nondimensional species standard-state enthalpies at the current temperature and pressure.
-
property
standard_entropies_R
¶ Array of nondimensional species standard-state entropies at the current temperature and pressure.
-
property
standard_gibbs_RT
¶ Array of nondimensional species standard-state Gibbs free energies at the current temperature and pressure.
-
property
standard_int_energies_RT
¶ Array of nondimensional species standard-state internal energies at the current temperature and pressure.
-
property
state_size
¶ Return size of vector defining internal state of the phase.
-
property
stoich_air_fuel_ratio
¶ Get the stoichiometric air to fuel ratio (kg oxidizer / kg fuel). Considers the oxidation of C to CO2, H to H2O and S to SO2. Other elements are assumed not to participate in oxidation (that is, N ends up as N2). The
basis
determines the composition of fuel and oxidizer:basis='mole'
(default) means mole fractions,basis='mass'
means mass fractions:>>> gas.set_mixture_fraction(0.5, 'CH3:0.5, CH3OH:.5, N2:0.125', 'O2:0.21, N2:0.79, NO:0.01') >>> gas.stoich_air_fuel_ratio('CH3:0.5, CH3OH:.5, N2:0.125', 'O2:0.21, N2:0.79, NO:0.01') 8.148040722239438
- Parameters
fuel – Fuel species name or mole/mass fractions as string, array, or dict.
oxidizer – Oxidizer species name or mole/mass fractions as a string, array, or dict.
basis – Determines if
fuel
andoxidizer
are given in mole fractions (basis='mole'
) or mass fractions (basis='mass'
)
-
property
thermal_conductivity
¶ Thermal conductivity. [W/m/K].
-
property
thermal_diff_coeffs
¶ Return a one-dimensional array of the species thermal diffusion coefficients [kg/m/s].
-
property
thermal_expansion_coeff
¶ Thermal expansion coefficient [1/K].
-
property
thermo_model
¶ Return thermodynamic model describing phase.
-
property
transport_model
¶ Get/Set the transport model associated with this transport model.
Setting a new transport model deletes the underlying C++ Transport object and replaces it with a new one implementing the specified model.
-
property
u
¶ Internal energy in [J/kg or J/kmol].
-
property
viscosity
¶ Viscosity [Pa-s].
-
property
volume_mass
¶ Specific volume [m^3/kg].
-
property
volume_mole
¶ Molar volume [m^3/kmol].
-
property
Representing Multiple States¶
-
class
cantera.
SolutionArray
(phase, shape=(0), states=None, extra=None, meta=None)¶ Bases:
object
A class providing a convenient interface for representing many thermodynamic states using the same
Solution
object and computing properties for that array of states.SolutionArray
can represent both 1D and multi-dimensional arrays of states, with shapes described in the same way as Numpy arrays. All of the states can be set in a single call:>>> gas = ct.Solution('gri30.yaml') >>> states = ct.SolutionArray(gas, (6, 10)) >>> T = np.linspace(300, 1000, 10) # row vector >>> P = ct.one_atm * np.linspace(0.1, 5.0, 6)[:,np.newaxis] # column vector >>> X = 'CH4:1.0, O2:1.0, N2:3.76' >>> states.TPX = T, P, X
Similar to Numpy arrays, input with fewer non-singleton dimensions than the
SolutionArray
is ‘broadcast’ to generate input of the appropriate shape. In the above example, the single value for the mole fraction input is applied to each input, while each row has a constant temperature and each column has a constant pressure.Computed properties are returned as Numpy arrays with the same shape as the array of states, with additional dimensions appended as necessary for non- scalar output (e.g. per-species or per-reaction properties):
>>> h = states.enthalpy_mass >>> h[i,j] # -> enthalpy at P[i] and T[j] >>> sk = states.partial_molar_entropies >>> sk[i,:,k] # -> entropy of species k at P[i] and each temperature >>> ropnet = states.net_rates_of_progress >>> ropnet[i,j,n] # -> net reaction rate for reaction n at P[i] and T[j]
In the case of 1D arrays, additional states can be appended one at a time:
>>> states = ct.SolutionArray(gas) # creates an empty SolutionArray >>> for phi in np.linspace(0.5, 2.0, 20): ... states.append(T=300, P=ct.one_atm, ... X={'CH4': phi, 'O2': 2, 'N2': 2*3.76}) >>> # 'states' now contains 20 elements >>> states.equilibrate('HP') >>> states.T # -> adiabatic flame temperature at various equivalence ratios
SolutionArray
objects can also be ‘sliced’ like Numpy arrays, which can be used both for accessing and setting properties:>>> states = ct.SolutionArray(gas, (6, 10)) >>> states[0].TP = 400, None # set the temperature of the first row to 400 K >>> cp = states[:,1].cp_mass # heat capacity of the second column
If many slices or elements of a property are going to be accessed (i.e. within a loop), it is generally more efficient to compute the property array once and access this directly, rather than repeatedly slicing the
SolutionArray
object, e.g.:>>> mu = states.viscosity >>> for i,j in np.ndindex(mu.shape): ... # do something with mu[i,j]
Information about a subset of species may also be accessed, using parentheses to specify the species:
>>> states('CH4').Y # -> mass fraction of CH4 in each state >>> states('H2','O2').partial_molar_cp # -> cp for H2 and O2
Properties and functions which are not dependent on the thermodynamic state act as pass-throughs to the underlying
Solution
object, and are not converted to arrays:>>> states.element_names ['O', 'H', 'C', 'N', 'Ar'] >>> s.reaction_equation(10) 'CH4 + O <=> CH3 + OH'
Data represented by a
SolutionArray
can be extracted and saved to a CSV file using thewrite_csv
method:>>> states.write_csv('somefile.csv', cols=('T', 'P', 'X', 'net_rates_of_progress'))
As long as stored columns specify a valid thermodynamic state, the contents of a
SolutionArray
can be restored using theread_csv
method:>>> states = ct.SolutionArray(gas) >>> states.read_csv('somefile.csv')
As an alternative to comma separated export and import, data extracted from
SolutionArray
objects can also be saved to and restored from a HDF container file using thewrite_hdf
:>>> states.write_hdf('somefile.h5', cols=('T', 'P', 'X'), group='some_key')
and
read_hdf
methods:>>> states = ct.SolutionArray(gas) >>> states.read_hdf('somefile.h5', key='some_key')
For HDF export and import, the (optional) keyword argument group allows for saving and accessing of multiple solutions in a single container file. Note that
write_hdf
andread_hdf
require a working installation of h5py. The packageh5py
can be installed using pip or conda.- Parameters
phase – The
Solution
object used to compute the thermodynamic, kinetic, and transport propertiesshape – A tuple or integer indicating the dimensions of the
SolutionArray
. If the shape is 1D, the array may be extended using theappend
method. Otherwise, the shape is fixed.states – The initial array of states. Used internally to provide slicing support.
-
property
DP
¶ Get/Set density [kg/m^3] and pressure [Pa].
-
property
DPQ
¶ Get the density [kg/m^3], pressure [Pa], and vapor fraction.
-
property
DPX
¶ Get/Set density [kg/m^3], pressure [Pa], and mole fractions.
-
property
DPY
¶ Get/Set density [kg/m^3], pressure [Pa], and mass fractions.
-
property
HP
¶ Get/Set enthalpy [J/kg or J/kmol] and pressure [Pa].
-
property
HPQ
¶ Get the enthalpy [J/kg or J/kmol], pressure [Pa] and vapor fraction.
-
property
HPX
¶ Get/Set enthalpy [J/kg or J/kmol], pressure [Pa] and mole fractions.
-
property
HPY
¶ Get/Set enthalpy [J/kg or J/kmol], pressure [Pa] and mass fractions.
-
property
ID
¶ The identifier of the object. The default value corresponds to the CTI/XML/YAML input file phase entry.
Deprecated since version 2.5: To be deprecated with version 2.5, and removed thereafter. Usage merged with
name
.
-
property
P
¶ Pressure [Pa].
-
property
PQ
¶ Get/Set the pressure [Pa] and vapor fraction of a two-phase state.
-
property
PV
¶ Get/Set the pressure [Pa] and specific volume [m^3/kg] of a PureFluid.
-
property
PX
¶ Get/Set the pressure [Pa] and vapor fraction of a two-phase state.
Deprecated since version 2.5: To be deprecated with version 2.5, and removed thereafter. Renamed to
PQ
.
-
property
P_sat
¶ Saturation pressure [Pa] at the current temperature.
-
property
Q
¶ Get/Set vapor fraction (quality). Can be set only when in the two-phase region.
-
property
SH
¶ Get/Set the specific entropy [J/kg/K] and the specific enthalpy [J/kg] of a PureFluid.
-
property
SP
¶ Get/Set entropy [J/kg/K or J/kmol/K] and pressure [Pa].
-
property
SPQ
¶ Get the entropy [J/kg/K or J/kmol/K], pressure [Pa], and vapor fraction.
-
property
SPX
¶ Get/Set entropy [J/kg/K or J/kmol/K], pressure [Pa], and mole fractions.
-
property
SPY
¶ Get/Set entropy [J/kg/K or J/kmol/K], pressure [Pa], and mass fractions.
-
property
ST
¶ Get/Set the entropy [J/kg/K] and temperature [K] of a PureFluid.
-
property
SV
¶ Get/Set entropy [J/kg/K or J/kmol/K] and specific volume [m^3/kg or m^3/kmol].
-
property
SVQ
¶ Get the entropy [J/kg/K or J/kmol/K], specific volume [m^3/kg or m^3/kmol], and vapor fraction.
-
property
SVX
¶ Get/Set entropy [J/kg/K or J/kmol/K], specific volume [m^3/kg or m^3/kmol], and mole fractions.
-
property
SVY
¶ Get/Set entropy [J/kg/K or J/kmol/K], specific volume [m^3/kg or m^3/kmol], and mass fractions.
-
property
T
¶ Temperature [K].
-
property
TD
¶ Get/Set temperature [K] and density [kg/m^3 or kmol/m^3].
-
property
TDQ
¶ Get the temperature [K], density [kg/m^3 or kmol/m^3], and vapor fraction.
-
property
TDX
¶ Get/Set temperature [K], density [kg/m^3 or kmol/m^3], and mole fractions.
-
property
TDY
¶ Get/Set temperature [K] and density [kg/m^3 or kmol/m^3], and mass fractions.
-
property
TH
¶ Get/Set the temperature [K] and the specific enthalpy [J/kg] of a PureFluid.
-
property
TP
¶ Get/Set temperature [K] and pressure [Pa].
-
property
TPQ
¶ Get/Set the temperature [K], pressure [Pa], and vapor fraction of a PureFluid.
An Exception is raised if the thermodynamic state is not consistent.
-
property
TPX
¶ Get/Set temperature [K], pressure [Pa], and mole fractions.
-
property
TPY
¶ Get/Set temperature [K], pressure [Pa], and mass fractions.
-
property
TQ
¶ Get/Set the temperature [K] and vapor fraction of a two-phase state.
-
property
TV
¶ Get/Set the temperature [K] and specific volume [m^3/kg] of a PureFluid.
-
property
TX
¶ Get/Set the temperature [K] and vapor fraction of a two-phase state.
Deprecated since version 2.5: To be deprecated with version 2.5, and removed thereafter. Renamed to
TQ
.
-
property
T_sat
¶ Saturation temperature [K] at the current pressure.
-
property
UP
¶ Get/Set the specific internal energy [J/kg] and the pressure [Pa] of a PureFluid.
-
property
UV
¶ Get/Set internal energy [J/kg or J/kmol] and specific volume [m^3/kg or m^3/kmol].
-
property
UVQ
¶ Get the internal energy [J/kg or J/kmol], specific volume [m^3/kg or m^3/kmol], and vapor fraction.
-
property
UVX
¶ Get/Set internal energy [J/kg or J/kmol], specific volume [m^3/kg or m^3/kmol], and mole fractions.
-
property
UVY
¶ Get/Set internal energy [J/kg or J/kmol], specific volume [m^3/kg or m^3/kmol], and mass fractions.
-
property
VH
¶ Get/Set the specific volume [m^3/kg] and the specific enthalpy [J/kg] of a PureFluid.
-
property
X
¶ Get/Set the species mole fractions. Can be set as an array, as a dictionary, or as a string. Always returns an array:
>>> phase.X = [0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5] >>> phase.X = {'H2':0.1, 'O2':0.4, 'AR':0.5} >>> phase.X = 'H2:0.1, O2:0.4, AR:0.5' >>> phase.X array([0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5])
-
property
Y
¶ Get/Set the species mass fractions. Can be set as an array, as a dictionary, or as a string. Always returns an array:
>>> phase.Y = [0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5] >>> phase.Y = {'H2':0.1, 'O2':0.4, 'AR':0.5} >>> phase.Y = 'H2:0.1, O2:0.4, AR:0.5' >>> phase.Y array([0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5])
-
append
(state=None, **kwargs)¶ Append an element to the array with the specified state. Elements can only be appended in cases where the array of states is one-dimensional.
The state may be specified in one of three ways:
as the array of [temperature, density, mass fractions] which is returned by
Solution.state
:mystates.append(gas.state)
as a tuple of three elements that corresponds to any of the full-state setters of
Solution
, e.g.TPY
orHPX
:mystates.append(TPX=(300, 101325, 'O2:1.0, N2:3.76'))
as separate keywords for each of the elements corresponding to one of the full-state setters:
mystates.append(T=300, P=101325, X={'O2':1.0, 'N2':3.76})
-
property
atomic_weight
¶ Atomic weight [kg/kmol] of element m
-
property
atomic_weights
¶ Array of atomic weight [kg/kmol] for each element in the mixture.
-
property
basis
¶ Determines whether intensive thermodynamic properties are treated on a
mass
(per kg) ormolar
(per kmol) basis. This affects the values returned by the propertiesh
,u
,s
,g
,v
,density
,cv
, andcp
, as well as the values used with the state-setting properties such asHPX
andUV
.
-
property
binary_diff_coeffs
¶ Binary diffusion coefficients [m^2/s].
-
property
charges
¶ Array of species charges [elem. charge].
-
property
chemical_potentials
¶ Array of species chemical potentials [J/kmol].
-
collect_data
(cols=None, tabular=False, threshold=0, species=None)¶ Returns the data specified by cols in an ordered dictionary, where keys correspond to SolutionArray attributes to be exported.
- Parameters
cols – A list of any properties of the solution that are scalars or which have a value for each species or reaction. If species names are specified, then either the mass or mole fraction of that species will be taken, depending on the value of species. cols may also include any arrays which were specified as ‘extra’ variables when defining the
SolutionArray
object. The special value ‘extra’ can be used to include all ‘extra’ variables.tabular – Split 2D data into separate 1D columns for each species / reaction
threshold – Relative tolerance for including a particular column if tabular output is enabled. The tolerance is applied by comparing the maximum absolute value for a particular column to the maximum absolute value in all columns for the same variable (e.g. mass fraction).
species – Specifies whether to use mass (‘Y’) or mole (‘X’) fractions for individual species specified in ‘cols’
-
property
concentrations
¶ Get/Set the species concentrations [kmol/m^3].
-
property
coverages
¶ Get/Set the fraction of sites covered by each species.
-
property
cp_mass
¶ Specific heat capacity at constant pressure [J/kg/K].
-
property
cp_mole
¶ Molar heat capacity at constant pressure [J/kmol/K].
-
property
creation_rates
¶ Creation rates for each species. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
critical_pressure
¶ Critical pressure [Pa].
-
property
critical_temperature
¶ Critical temperature [K].
-
property
cv_mass
¶ Specific heat capacity at constant volume [J/kg/K].
-
property
cv_mole
¶ Molar heat capacity at constant volume [J/kmol/K].
-
property
delta_enthalpy
¶ Change in enthalpy for each reaction [J/kmol].
-
property
delta_entropy
¶ Change in entropy for each reaction [J/kmol/K].
-
property
delta_gibbs
¶ Change in Gibbs free energy for each reaction [J/kmol].
-
property
delta_standard_enthalpy
¶ Change in standard-state enthalpy (independent of composition) for each reaction [J/kmol].
-
property
delta_standard_entropy
¶ Change in standard-state entropy (independent of composition) for each reaction [J/kmol/K].
-
property
delta_standard_gibbs
¶ Change in standard-state Gibbs free energy (independent of composition) for each reaction [J/kmol].
-
property
density_mass
¶ (Mass) density [kg/m^3].
-
property
density_mole
¶ Molar density [kmol/m^3].
-
property
destruction_rates
¶ Destruction rates for each species. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
electric_potential
¶ Get/Set the electric potential [V] for this phase.
-
property
electrical_conductivity
¶ Electrical conductivity. [S/m].
-
property
electrochemical_potentials
¶ Array of species electrochemical potentials [J/kmol].
-
property
element_index
¶ The index of element element, which may be specified as a string or an integer. In the latter case, the index is checked for validity and returned. If no such element is present, an exception is thrown.
-
property
element_name
¶ Name of the element with index m.
-
property
element_names
¶ A list of all the element names.
-
elemental_mass_fraction
(*args, **kwargs)¶
-
elemental_mole_fraction
(*args, **kwargs)¶
-
property
enthalpy_mass
¶ Specific enthalpy [J/kg].
-
property
enthalpy_mole
¶ Molar enthalpy [J/kmol].
-
property
entropy_mass
¶ Specific entropy [J/kg/K].
-
property
entropy_mole
¶ Molar entropy [J/kmol/K].
-
equilibrate
(*args, **kwargs)¶
-
property
equilibrium_constants
¶ Equilibrium constants in concentration units for all reactions.
-
property
forward_rate_constants
¶ Forward rate constants for all reactions. The computed values include all temperature-dependent, pressure-dependent, and third body contributions. Units are a combination of kmol, m^3 and s, that depend on the rate expression for the reaction.
-
property
forward_rates_of_progress
¶ Forward rates of progress for the reactions. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
from_pandas
(df)¶ Restores
SolutionArray
data from a pandas DataFrame df.This method is intendend for loading of data that were previously exported by
to_pandas
. The method requires a working pandas installation. The package ‘pandas’ can be installed using pip or conda.
-
property
gibbs_mass
¶ Specific Gibbs free energy [J/kg].
-
property
gibbs_mole
¶ Molar Gibbs free energy [J/kmol].
-
property
heat_production_rates
¶ Get the volumetric heat production rates [W/m^3] on a per-reaction basis. The sum over all reactions results in the total volumetric heat release rate. Example: C. K. Law: Combustion Physics (2006), Fig. 7.8.6
>>> gas.heat_production_rates[1] # heat production rate of the 2nd reaction
-
property
heat_release_rate
¶ Get the total volumetric heat release rate [W/m^3].
-
property
int_energy_mass
¶ Specific internal energy [J/kg].
-
property
int_energy_mole
¶ Molar internal energy [J/kmol].
-
property
is_reversible
¶ True if reaction
i_reaction
is reversible.
-
property
isothermal_compressibility
¶ Isothermal compressibility [1/Pa].
-
property
kinetics_species_index
¶ The index of species species of phase phase within arrays returned by methods of class
Kinetics
. If species is a string, the phase argument is unused.
-
property
max_temp
¶ Maximum temperature for which the thermodynamic data for the phase are valid.
-
property
mean_molecular_weight
¶ The mean molecular weight (molar mass) [kg/kmol].
-
property
meta
¶ Dictionary holding information describing the
SolutionArray
. Metadata should be provided for the creation ofSolutionArray
objects.
-
property
min_temp
¶ Minimum temperature for which the thermodynamic data for the phase are valid.
-
property
mix_diff_coeffs
¶ Mixture-averaged diffusion coefficients [m^2/s] relating the mass-averaged diffusive fluxes (with respect to the mass averaged velocity) to gradients in the species mole fractions.
-
property
mix_diff_coeffs_mass
¶ Mixture-averaged diffusion coefficients [m^2/s] relating the diffusive mass fluxes to gradients in the species mass fractions.
-
property
mix_diff_coeffs_mole
¶ Mixture-averaged diffusion coefficients [m^2/s] relating the molar diffusive fluxes to gradients in the species mole fractions.
-
property
modify_reaction
¶ Modify the
Reaction
with indexirxn
to have the same rate parameters asrxn
.rxn
must have the same reactants and products and be of the same type (i.e.ElementaryReaction
,FalloffReaction
,PlogReaction
, etc.) as the existing reaction. This method does not modify the third-body efficiencies, reaction orders, or reversibility of the reaction.
-
property
molecular_weights
¶ Array of species molecular weights (molar masses) [kg/kmol].
-
property
multi_diff_coeffs
¶ Multicomponent diffusion coefficients [m^2/s].
-
property
multiplier
¶ A scaling factor applied to the rate coefficient for reaction i_reaction. Can be used to carry out sensitivity analysis or to selectively disable a particular reaction. See
set_multiplier
.
-
property
n_atoms
¶ Number of atoms of element element in species species. The element and species may be specified by name or by index.
>>> phase.n_atoms('CH4','H') 4
-
property
n_elements
¶ Number of elements.
-
property
n_phases
¶ Number of phases in the reaction mechanism.
-
property
n_reactions
¶ Number of reactions in the reaction mechanism.
-
property
n_species
¶ Number of species.
-
property
n_total_species
¶ Total number of species in all phases participating in the kinetics mechanism.
-
property
name
¶ The name assigned to this object. The default value corresponds to the CTI/XML/YAML input file phase entry.
-
property
net_production_rates
¶ Net production rates for each species. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
net_rates_of_progress
¶ Net rates of progress for the reactions. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
property
partial_molar_cp
¶ Array of species partial molar specific heat capacities at constant pressure [J/kmol/K].
-
property
partial_molar_enthalpies
¶ Array of species partial molar enthalpies [J/kmol].
-
property
partial_molar_entropies
¶ Array of species partial molar entropies [J/kmol/K].
-
property
partial_molar_int_energies
¶ Array of species partial molar internal energies [J/kmol].
-
property
partial_molar_volumes
¶ Array of species partial molar volumes [m^3/kmol].
-
property
phase_of_matter
¶ Get the thermodynamic phase (gas, liquid, etc.) at the current conditions.
-
property
product_stoich_coeff
¶ The stoichiometric coefficient of species k_spec as a product in reaction i_reaction.
-
property
product_stoich_coeffs
¶ The array of product stoichiometric coefficients. Element [k,i] of this array is the product stoichiometric coefficient of species k in reaction i.
-
property
products
¶ The products portion of the reaction equation
-
property
reactant_stoich_coeff
¶ The stoichiometric coefficient of species k_spec as a reactant in reaction i_reaction.
-
property
reactant_stoich_coeffs
¶ The array of reactant stoichiometric coefficients. Element [k,i] of this array is the reactant stoichiometric coefficient of species k in reaction i.
-
property
reactants
¶ The reactants portion of the reaction equation
-
property
reaction
¶ Return a
Reaction
object representing the reaction with indexi_reaction
. Changes to this object do not affect theKinetics
orSolution
object until themodify_reaction
function is called.
-
property
reaction_equation
¶ The equation for the specified reaction. See also
reaction_equations
.
-
property
reaction_equations
¶ Returns a list containing the reaction equation for all reactions in the mechanism (if indices is unspecified) or the equations for each reaction in the sequence indices. For example:
>>> gas.reaction_equations() ['2 O + M <=> O2 + M', 'O + H + M <=> OH + M', 'O + H2 <=> H + OH', ...] >>> gas.reaction_equations([2,3]) ['O + H + M <=> OH + M', 'O + H2 <=> H + OH']
See also
reaction_equation
.
-
property
reaction_phase_index
¶ The index of the phase where the reactions occur.
-
property
reaction_type
¶ Type of reaction i_reaction.
-
property
reactions
¶ Return a list of all
Reaction
objects. Changes to these objects do not affect theKinetics
orSolution
object until themodify_reaction
function is called.
-
read_csv
(filename)¶ Read a CSV file named filename and restore data to the
SolutionArray
usingrestore_data
. This method allows for recreation of data previously exported bywrite_csv
.
-
read_hdf
(filename, group=None, subgroup=None, force=False)¶ Read a dataset from a HDF container file and restore data to the
SolutionArray
object. This method allows for recreation of data previously exported bywrite_hdf
.- Parameters
filename – name of the HDF container file; typical file extensions are
.hdf
,.hdf5
or.h5
.group – Identifier for the group in the container file. A group may contain a
SolutionArray
object or additional subgroups.subgroup – Optional name identifier for a subgroup representing a
SolutionArray
object to be read. If ‘None’, no subgroup is assumed to exist.force – If False, matching
SolutionArray
source identifiers are enforced (e.g. input file used for the creation of the underlyingSolution
object), with an error being raised if the current source does not match the original source. If True, the error is suppressed.
- Returns
User-defined attributes provided to describe the group holding the
SolutionArray
information.
The method imports data using
restore_data
and requires a working installation of h5py (h5py
can be installed using pip or conda).
-
property
reference_pressure
¶ Reference state pressure [Pa].
-
restore_data
(data)¶ Restores a
SolutionArray
based on data specified in an ordered dictionary. Thus, this method allows to restore data exported bycollect_data
.- Parameters
data – Dictionary holding data to be restored, where keys refer to thermodynamic states (e.g.
T
,density
) or extra entries, and values contain corresponding data.
The receiving
SolutionArray
either has to be empty or should have matching dimensions. Essential state properties and extra entries are detected automatically whereas stored information of calculated properties is omitted. If the receivingSolutionArray
has extra entries already specified, only those will be imported; if labels does not contain those entries, an error is raised.
-
property
reverse_rate_constants
¶ Reverse rate constants for all reactions. The computed values include all temperature-dependent, pressure-dependent, and third body contributions. Units are a combination of kmol, m^3 and s, that depend on the rate expression for the reaction.
-
property
reverse_rates_of_progress
¶ Reverse rates of progress for the reactions. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases.
-
set_equivalence_ratio
(phi, *args, **kwargs)¶ See
ThermoPhase.set_equivalence_ratio
Note that phi either needs to be a scalar value or dimensions have to be matched to the SolutionArray.
-
property
set_multiplier
¶ Set the multiplier for for reaction i_reaction to value. If i_reaction is not specified, then the multiplier for all reactions is set to value. See
multiplier
.
-
property
site_density
¶ Get/Set the site density. [kmol/m^2] for surface phases; [kmol/m] for edge phases.
-
sort
(col, reverse=False)¶ Sort SolutionArray by column col.
- Parameters
col – Column that is used to sort the SolutionArray.
reverse – If True, the sorted list is reversed (descending order).
-
property
source
¶ The source of this object (such as a file name).
-
property
species
¶ Return the
Species
object for species k, where k is either the species index or the species name. If k is not specified, a list of all species objects is returned. Changes to this object do not affect theThermoPhase
orSolution
object until themodify_species
function is called.
-
property
species_index
¶ The index of species species, which may be specified as a string or an integer. In the latter case, the index is checked for validity and returned. If no such species is present, an exception is thrown.
-
property
species_name
¶ Name of the species with index k.
-
property
species_names
¶ A list of all the species names.
-
property
standard_cp_R
¶ Array of nondimensional species standard-state specific heat capacities at constant pressure at the current temperature and pressure.
-
property
standard_enthalpies_RT
¶ Array of nondimensional species standard-state enthalpies at the current temperature and pressure.
-
property
standard_entropies_R
¶ Array of nondimensional species standard-state entropies at the current temperature and pressure.
-
property
standard_gibbs_RT
¶ Array of nondimensional species standard-state Gibbs free energies at the current temperature and pressure.
-
property
standard_int_energies_RT
¶ Array of nondimensional species standard-state internal energies at the current temperature and pressure.
-
property
thermal_conductivity
¶ Thermal conductivity. [W/m/K].
-
property
thermal_diff_coeffs
¶ Return a one-dimensional array of the species thermal diffusion coefficients [kg/m/s].
-
property
thermal_expansion_coeff
¶ Thermal expansion coefficient [1/K].
-
to_pandas
(cols=None, *args, **kwargs)¶ Returns the data specified by cols in a single pandas DataFrame.
Additional arguments are passed on to
collect_data
. This method works only with 1DSolutionArray
objects and requires a working pandas installation. Use pip or conda to installpandas
to enable this method.
-
property
transport_model
¶ Get/Set the transport model associated with this transport model.
Setting a new transport model deletes the underlying C++ Transport object and replaces it with a new one implementing the specified model.
-
property
u
¶ Internal energy in [J/kg or J/kmol].
-
property
viscosity
¶ Viscosity [Pa-s].
-
property
volume_mass
¶ Specific volume [m^3/kg].
-
property
volume_mole
¶ Molar volume [m^3/kmol].
-
write_csv
(filename, cols=None, *args, **kwargs)¶ Write a CSV file named filename containing the data specified by cols. The first row of the CSV file will contain column labels.
Additional arguments are passed on to
collect_data
. This method works only with 1DSolutionArray
objects.
-
write_hdf
(filename, *args, cols=None, group=None, subgroup=None, attrs={}, mode='a', append=False, compression=None, compression_opts=None, **kwargs)¶ Write data specified by cols to a HDF container file named filename. Note that it is possible to write multiple data entries to a single HDF container file, where group is used to differentiate data.
An example for the default HDF file structure is::
/ Group /group0 Group /group0/some_attr Attribute ... /group0/T Dataset ... /group0/phase Group /group0/phase/name Attribute /group0/phase/source Attribute
where
group0
is the default name for the top level HDF entry. In addition to datasets, information stored inSolutionArray.meta
is saved in form of HDF attributes. An additional intermediate layer may be created using the subgroup argument.- Parameters
filename – Name of the HDF container file; typical file extensions are
.hdf
,.hdf5
or.h5
.cols – A list of any properties of the solution being exported.
group – Identifier for the group in the container file. If no subgroup is specified, a group represents a
SolutionArray
. If ‘None’, group names default to ‘groupN’, with N being the number of pre-existing groups within the HDF container file.subgroup – Name identifier for an optional subgroup, with subgroups representing individual
SolutionArray
objects. If ‘None’, no subgroup is created.attrs – Dictionary of user-defined attributes added at the group level (typically used in conjunction with a subgroup argument).
mode – Mode h5py uses to open the output file {‘a’ to read/write if file exists, create otherwise (default); ‘w’ to create file, truncate if exists; ‘r+’ to read/write, file must exist}.
append – If False, the content of a pre-existing group is deleted before writing the
SolutionArray
in the first position. If True, the currentSolutionArray
objects is appended to the group.compression – Pre-defined h5py compression filters {None, ‘gzip’, ‘lzf’, ‘szip’} used for data compression.
compression_opts – Options for the h5py compression filter; for ‘gzip’, this corresponds to the compression level {None, 0-9}.
- Returns
Group identifier used for storing HDF data.
Arguments compression, and compression_opts are mapped to parameters for
h5py.create_dataset
; in both cases, the choices ofNone
results in default values set by h5py.Additional arguments (i.e. args and kwargs) are passed on to
collect_data
; seecollect_data
for further information. This method requires a working installation of h5py (h5py
can be installed using pip or conda).
Utility Functions¶
-
cantera.
add_directory
(directory)¶ Add a directory to search for Cantera data files.
-
cantera.
get_data_directories
()¶ Get a list of the directories Cantera searches for data files.