Cantera

Table Of Contents

Previous topic

Zero-Dimensional Reactor Networks

Next topic

Physical Constants

This Page

Warning

This documentation is for an old version of Cantera. You can find docs for newer versions here.

One-dimensional Reacting Flows

Composite Domains

class cantera.FreeFlame(self, gas, grid=None)

Bases: cantera._cantera.FlameBase

A freely-propagating flat flame.

A domain of type FreeFlow named ‘flame’ will be created to represent the flame. The three domains comprising the stack are stored as self.inlet, self.flame, and self.outlet.

set_initial_guess(self)

Set the initial guess for the solution. The adiabatic flame temperature and equilibrium composition are computed for the inlet gas composition. The temperature profile rises linearly over 20% of the domain width to Tad, then is flat. The mass fraction profiles are set similarly.

class cantera.BurnerFlame(self, gas, grid=None)

Bases: cantera._cantera.FlameBase

A burner-stabilized flat flame.

Parameters:
  • gasSolution (using the IdealGas thermodynamic model) used to evaluate all gas properties and reaction rates.
  • grid – Array of initial grid points

A domain of class AxisymmetricStagnationFlow named flame will be created to represent the flame. The three domains comprising the stack are stored as self.burner, self.flame, and self.outlet.

set_initial_guess(self)

Set the initial guess for the solution. The adiabatic flame temperature and equilibrium composition are computed for the burner gas composition. The temperature profile rises linearly in the first 20% of the flame to Tad, then is flat. The mass fraction profiles are set similarly.

class cantera.CounterflowDiffusionFlame(self, gas, grid=None)

Bases: cantera._cantera.FlameBase

A counterflow diffusion flame

Parameters:
  • gasSolution (using the IdealGas thermodynamic model) used to evaluate all gas properties and reaction rates.
  • grid – Array of initial grid points

A domain of class AxisymmetricStagnationFlow named flame will be created to represent the flame. The three domains comprising the stack are stored as self.fuel_inlet, self.flame, and self.oxidizer_inlet.

set_initial_guess(self, fuel, oxidizer='O2', stoich=None)

Set the initial guess for the solution. The fuel species must be specified:

>>> f.set_initial_guess(fuel='CH4')

The oxidizer and corresponding stoichiometry must be specified if it is not ‘O2’. The initial guess is generated by assuming infinitely- fast chemistry.

Flow Domains

class cantera.FreeFlow

Bases: cantera._cantera._FlowBase

P

Pressure [Pa]

bounds(self, component)

Return the (lower, upper) bounds for a solution component.

>>> d.bounds('T')
(200.0, 5000.0)
component_index(self, str name)

Index of the component with name ‘name’

component_name(self, int n)

Name of the nth component.

component_names

List of the names of all components of this domain.

description

A description of this domain

energy_enabled

Determines whether or not to solve the energy equation.

grid

The grid for this domain

index

Index of this domain in a stack. Returns -1 if this domain is not part of a stack.

n_components

Number of solution components at each grid point.

n_points

Number of grid points belonging to this domain.

name

The name / id of this domain

set_bounds(self, *, default=None, Y=None, **kwargs)

Set the lower and upper bounds on the solution.

The argument list should consist of keyword/value pairs, with component names as keywords and (lower_bound, upper_bound) tuples as the values. The keyword default may be used to specify default bounds for all unspecified components. The keyword Y can be used to stand for all species mass fractions in flow domains.

>>> d.set_bounds(default=(0, 1), Y=(-1.0e-5, 2.0))
set_fixed_temp_profile(self, pos, T)

Set the fixed temperature profile. This profile is used whenever the energy equation is disabled.

Parameters:
  • pos – arrray of relative positions from 0 to 1
  • temp – array of temperature values
>>> d.set_fixed_temp_profile(array([0.0, 0.5, 1.0]),
...                          array([500.0, 1500.0, 2000.0])
set_steady_tolerances(self, *, default=None, Y=None, **kwargs)

Set the error tolerances for the steady-state problem.

The argument list should consist of keyword/value pairs, with component names as keywords and (rtol, atol) tuples as the values. The keyword default may be used to specify default bounds for all unspecified components. The keyword Y can be used to stand for all species mass fractions in flow domains.

set_transient_tolerances(self, *, default=None, Y=None, **kwargs)

Set the error tolerances for the steady-state problem.

The argument list should consist of keyword/value pairs, with component names as keywords and (rtol, atol) tuples as the values. The keyword default may be used to specify default bounds for all unspecified components. The keyword Y can be used to stand for all species mass fractions in flow domains.

set_transport(self, _SolutionBase phase)
soret_enabled

Determines whether or not to include diffusive mass fluxes due to the Soret effect. Enabling this option works only when using the multicomponent transport model.

tolerances(self, component)

Return the (relative, absolute) error tolerances for a solution component.

>>> rtol, atol = d.tolerances('u')
class cantera.AxisymmetricStagnationFlow

Bases: cantera._cantera._FlowBase

An axisymmetric flow domain.

In an axisymmetric flow domain, the equations solved are the similarity equations for the flow in a finite-height gap of infinite radial extent. The solution variables are:

u
axial velocity
V
radial velocity divided by radius
T
temperature
lambda
(1/r)(dP/dr)
Y_k
species mass fractions

It may be shown that if the boundary conditions on these variables are independent of radius, then a similarity solution to the exact governing equations exists in which these variables are all independent of radius. This solution holds only in in low-Mach-number limit, in which case (dP/dz) = 0, and lambda is a constant. (Lambda is treated as a spatially- varying solution variable for numerical reasons, but in the final solution it is always independent of z.) As implemented here, the governing equations assume an ideal gas mixture. Arbitrary chemistry is allowed, as well as arbitrary variation of the transport properties.

P

Pressure [Pa]

bounds(self, component)

Return the (lower, upper) bounds for a solution component.

>>> d.bounds('T')
(200.0, 5000.0)
component_index(self, str name)

Index of the component with name ‘name’

component_name(self, int n)

Name of the nth component.

component_names

List of the names of all components of this domain.

description

A description of this domain

energy_enabled

Determines whether or not to solve the energy equation.

grid

The grid for this domain

index

Index of this domain in a stack. Returns -1 if this domain is not part of a stack.

n_components

Number of solution components at each grid point.

n_points

Number of grid points belonging to this domain.

name

The name / id of this domain

set_bounds(self, *, default=None, Y=None, **kwargs)

Set the lower and upper bounds on the solution.

The argument list should consist of keyword/value pairs, with component names as keywords and (lower_bound, upper_bound) tuples as the values. The keyword default may be used to specify default bounds for all unspecified components. The keyword Y can be used to stand for all species mass fractions in flow domains.

>>> d.set_bounds(default=(0, 1), Y=(-1.0e-5, 2.0))
set_fixed_temp_profile(self, pos, T)

Set the fixed temperature profile. This profile is used whenever the energy equation is disabled.

Parameters:
  • pos – arrray of relative positions from 0 to 1
  • temp – array of temperature values
>>> d.set_fixed_temp_profile(array([0.0, 0.5, 1.0]),
...                          array([500.0, 1500.0, 2000.0])
set_steady_tolerances(self, *, default=None, Y=None, **kwargs)

Set the error tolerances for the steady-state problem.

The argument list should consist of keyword/value pairs, with component names as keywords and (rtol, atol) tuples as the values. The keyword default may be used to specify default bounds for all unspecified components. The keyword Y can be used to stand for all species mass fractions in flow domains.

set_transient_tolerances(self, *, default=None, Y=None, **kwargs)

Set the error tolerances for the steady-state problem.

The argument list should consist of keyword/value pairs, with component names as keywords and (rtol, atol) tuples as the values. The keyword default may be used to specify default bounds for all unspecified components. The keyword Y can be used to stand for all species mass fractions in flow domains.

set_transport(self, _SolutionBase phase)
soret_enabled

Determines whether or not to include diffusive mass fluxes due to the Soret effect. Enabling this option works only when using the multicomponent transport model.

tolerances(self, component)

Return the (relative, absolute) error tolerances for a solution component.

>>> rtol, atol = d.tolerances('u')

Boundaries

class cantera.Inlet1D

Bases: cantera._cantera.Boundary1D

A one-dimensional inlet. Note that an inlet can only be a terminal domain - it must be either the leftmost or rightmost domain in a stack.

spread_rate
class cantera.Outlet1D

Bases: cantera._cantera.Boundary1D

A one-dimensional outlet. An outlet imposes a zero-gradient boundary condition on the flow.

class cantera.OutletReservoir1D

Bases: cantera._cantera.Boundary1D

A one-dimensional outlet into a reservoir.

class cantera.SymmetryPlane1D

Bases: cantera._cantera.Boundary1D

A symmetry plane.

class cantera.Surface1D

Bases: cantera._cantera.Boundary1D

A solid surface.

class cantera.ReactingSurface1D

Bases: cantera._cantera.Boundary1D

A reacting solid surface.

coverage_enabled

Controls whether or not to solve the surface coverage equations.

set_kinetics(self, Kinetics kin)

Set the kinetics manager (surface reaction mechanism object).

Base Classes

class cantera.Domain1D

Bases: object

Domain1D(name=None, *args, **kwargs)

bounds(self, component)

Return the (lower, upper) bounds for a solution component.

>>> d.bounds('T')
(200.0, 5000.0)
component_index(self, str name)

Index of the component with name ‘name’

component_name(self, int n)

Name of the nth component.

component_names

List of the names of all components of this domain.

description

A description of this domain

grid

The grid for this domain

index

Index of this domain in a stack. Returns -1 if this domain is not part of a stack.

n_components

Number of solution components at each grid point.

n_points

Number of grid points belonging to this domain.

name

The name / id of this domain

set_bounds(self, *, default=None, Y=None, **kwargs)

Set the lower and upper bounds on the solution.

The argument list should consist of keyword/value pairs, with component names as keywords and (lower_bound, upper_bound) tuples as the values. The keyword default may be used to specify default bounds for all unspecified components. The keyword Y can be used to stand for all species mass fractions in flow domains.

>>> d.set_bounds(default=(0, 1), Y=(-1.0e-5, 2.0))
set_steady_tolerances(self, *, default=None, Y=None, **kwargs)

Set the error tolerances for the steady-state problem.

The argument list should consist of keyword/value pairs, with component names as keywords and (rtol, atol) tuples as the values. The keyword default may be used to specify default bounds for all unspecified components. The keyword Y can be used to stand for all species mass fractions in flow domains.

set_transient_tolerances(self, *, default=None, Y=None, **kwargs)

Set the error tolerances for the steady-state problem.

The argument list should consist of keyword/value pairs, with component names as keywords and (rtol, atol) tuples as the values. The keyword default may be used to specify default bounds for all unspecified components. The keyword Y can be used to stand for all species mass fractions in flow domains.

tolerances(self, component)

Return the (relative, absolute) error tolerances for a solution component.

>>> rtol, atol = d.tolerances('u')
class cantera.Boundary1D

Bases: cantera._cantera.Domain1D

Boundary1D(_SolutionBase phase, *args, **kwargs)

Base class for boundary domains.

Parameters:phase – The (gas) phase corresponding to the adjacent flow domain
T

The temperature [K] at this boundary.

X

Species mole fractions at this boundary. May be set as either a string or as an array.

Y

Species mass fractions at this boundary. May be set as either a string or as an array.

mdot

The mass flow rate per unit area [kg/m^2]

class cantera.Sim1D

Bases: object

Sim1D(domains, *args, **kwargs)

Class Sim1D is a container for one-dimensional domains. It also holds the multi-domain solution vector, and controls the process of finding the solution.

Domains are ordered left-to-right, with domain number 0 at the left.

domain_index(self, dom)

Get the index of a domain, specified either by name or as a Domain1D object.

domains
profile(self, domain, component)

Spatial profile of one component in one domain.

Parameters:
  • domain – Domain1D object, name, or index
  • component – component name or index
>>> T = s.profile(flow, 'T')
refine(self, loglevel=1)

Refine the grid, adding points where solution is not adequately resolved.

restore(self, filename='soln.xml', name='solution', loglevel=2)

Set the solution vector to a previously-saved solution.

Parameters:
  • filename – solution file
  • name – solution name within the file
  • loglevel – Amount of logging information to display while restoring, from 0 (disabled) to 2 (most verbose).
>>> s.restore(filename='save.xml', name='energy_off')
save(self, filename='soln.xml', name='solution', description='none', loglevel=1)

Save the solution in XML format. :param filename:

solution file
Parameters:
  • name – solution name within the file
  • description – custom description text
>>> s.save(filename='save.xml', name='energy_off',
...        description='solution with energy eqn. disabled')
set_fixed_temperature(self, T)

Set the temperature used to fix the spatial location of a freely propagating flame.

set_flat_profile(self, domain, component, value)

Set a flat profile for one component in one domain.

Parameters:
  • domain – Domain1D object, name, or index
  • component – component name or index
  • v – value
>>> s.set_flat_profile(d, 'u', -3.0)
set_grid_min(self, dz, domain=None)

Set the minimum grid spacing on domain. If domain is None, then set the grid spacing for all domains.

set_initial_guess(self)

Set the initial guess for the solution. Derived classes extend this function to set approximations for the temperature and composition profiles.

set_interrupt(self, f)

Set an interrupt function to be called each time that OneDim::eval is called. The signature of f is float f(float). The default interrupt function is used to trap KeyboardInterrupt exceptions so that ctrl-c can be used to break out of the C++ solver loop.

set_max_jac_age(self, ss_age, ts_age)

Set the maximum number of times the Jacobian will be used before it must be re-evaluated.

Parameters:
  • ss_age – age criterion during steady-state mode
  • ts_age – age criterion during time-stepping mode
set_max_time_step(self, tsmax)

Set the maximum time step.

set_min_time_step(self, tsmin)

Set the minimum time step.

set_profile(self, domain, component, positions, values)

Set an initial estimate for a profile of one component in one domain.

Parameters:
  • domain – Domain1D object, name, or index
  • component – component name or index
  • positions – sequence of relative positions, from 0 on the left to 1 on the right
  • values – sequence of values at the relative positions specified in positions
>>> s.set_profile(d, 'T', [0.0, 0.2, 1.0], [400.0, 800.0, 1500.0])
set_refine_criteria(self, domain, ratio=10.0, slope=0.8, curve=0.8, prune=0.05)

Set the criteria used to refine one domain.

Parameters:
  • domain – domain object, index, or name
  • ratio – additional points will be added if the ratio of the spacing on either side of a grid point exceeds this value
  • slope – maximum difference in value between two adjacent points, scaled by the maximum difference in the profile (0.0 < slope < 1.0). Adds points in regions of high slope.
  • curve – maximum difference in slope between two adjacent intervals, scaled by the maximum difference in the profile (0.0 < curve < 1.0). Adds points in regions of high curvature.
  • prune – if the slope or curve criteria are satisfied to the level of ‘prune’, the grid point is assumed not to be needed and is removed. Set prune significantly smaller than ‘slope’ and ‘curve’. Set to zero to disable pruning the grid.
>>> s.set_refine_criteria(d, ratio=5.0, slope=0.2, curve=0.3, prune=0.03)
set_time_step(self, stepsize, n_steps)

Set the sequence of time steps to try when Newton fails.

Parameters:
  • stepsize – initial time step size [s]
  • n_steps – sequence of integer step numbers
>>> s.set_time_step(1.0e-5, [1, 2, 5, 10])
set_time_step_factor(self, tfactor)

Set the factor by which the time step will be increased after a successful step, or decreased after an unsuccessful one.

set_value(self, domain, component, point, value)

Set the value of one component in one domain at one point to ‘value’.

Parameters:
  • domain – Domain1D object, name, or index
  • component – component name or index
  • point – grid point number within domain starting with 0 on the left
  • value – numerical value
>>> s.set(d, 3, 5, 6.7)
>>> s.set(1, 0, 5, 6.7)
>>> s.set('flow', 'T', 5, 500)
show_solution(self)

print the current solution.

show_stats(self, print_time=True)

Show the statistics for the last solution.

If invoked with no arguments or with a non-zero argument, the timing statistics will be printed. Otherwise, the timing will not be printed.

solve(self, loglevel=1, refine_grid=True)

Solve the problem.

Parameters:
  • loglevel – integer flag controlling the amount of diagnostic output. Zero suppresses all output, and 5 produces very verbose output.
  • refine_grid – if True, enable grid refinement.
value(self, domain, component, point)

Solution value at one point

Parameters:
  • domain – Domain1D object, name, or index
  • component – component name or index
  • point – grid point number within domain starting with 0 on the left
>>> t = s.value('flow', 'T', 6)
work_value(self, domain, component, point)

Internal work array value at one point. After calling eval, this array contains the values of the residual function.

Parameters:
  • domain – Domain1D object, name, or index
  • component – component name or index
  • point – grid point number in the domain, starting with zero at the left
>>> t = s.value(flow, 'T', 6)
class cantera.FlameBase(self, domains, gas, grid=None)

Bases: cantera._cantera.Sim1D

Base class for flames with a single flow domain

Parameters:
  • gas – object to use to evaluate all gas properties and reaction rates
  • grid – array of initial grid points
L

FlameBase.L(self)

Array containing the radial pressure gradient (1/r)(dP/dr) [N/m^4] at each point. Note: This value is named ‘lambda’ in the C++ code.

P

FlameBase.P(self)

T

FlameBase.T(self) Array containing the temperature [K] at each grid point.

V

FlameBase.V(self)

Array containing the tangential velocity gradient [1/s] at each point.

X

Get/Set the species mole fractions. Can be set as either an array 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
array([0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5])

Returns an array of size n_species x n_points.

Y

Get/Set the species mass fractions. Can be set as either an array 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
array([0.1, 0, 0, 0.4, 0, 0, 0, 0, 0.5])

Returns an array of size n_species x n_points.

chemical_potentials

Array of species chemical potentials [J/kmol]. Returns an array of size n_species x n_points.

concentrations

Get/Set the species concentrations [kmol/m^3]. Returns an array of size n_species x n_points.

cp

Heat capacity at constant pressure [J/kg/K or J/kmol/K] depending on basis. Returns an array of length n_points.

cp_mass

Specific heat capacity at constant pressure [J/kg/K]. Returns an array of length n_points.

cp_mole

Molar heat capacity at constant pressure [J/kmol/K]. Returns an array of length n_points.

creation_rates

Creation rates for each species. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases. Returns an array of size n_species x n_points.

cv

Heat capacity at constant volume [J/kg/K or J/kmol/K] depending on basis. Returns an array of length n_points.

cv_mass

Specific heat capacity at constant volume [J/kg/K]. Returns an array of length n_points.

cv_mole

Molar heat capacity at constant volume [J/kmol/K]. Returns an array of length n_points.

delta_enthalpy

Change in enthalpy for each reaction [J/kmol]. Returns an array of size n_reactions x n_points.

delta_entropy

Change in entropy for each reaction [J/kmol/K]. Returns an array of size n_reactions x n_points.

delta_gibbs

Change in Gibbs free energy for each reaction [J/kmol]. Returns an array of size n_reactions x n_points.

delta_standard_enthalpy

Change in standard-state enthalpy (independent of composition) for each reaction [J/kmol]. Returns an array of size n_reactions x n_points.

delta_standard_entropy

Change in standard-state entropy (independent of composition) for each reaction [J/kmol/K]. Returns an array of size n_reactions x n_points.

delta_standard_gibbs

Change in standard-state Gibbs free energy (independent of composition) for each reaction [J/kmol]. Returns an array of size n_reactions x n_points.

density

Density [kg/m^3 or kmol/m^3] depending on basis. Returns an array of length n_points.

density_mass

(Mass) density [kg/m^3]. Returns an array of length n_points.

density_mole

Molar density [kmol/m^3]. Returns an array of length n_points.

destruction_rates

Destruction rates for each species. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases. Returns an array of size n_species x n_points.

electrochemical_potentials

Array of species electrochemical potentials [J/kmol]. Returns an array of size n_species x n_points.

energy_enabled

FlameBase.energy_enabled(self)

enthalpy_mass

Specific enthalpy [J/kg]. Returns an array of length n_points.

enthalpy_mole

Molar enthalpy [J/kmol]. Returns an array of length n_points.

entropy_mass

Specific entropy [J/kg]. Returns an array of length n_points.

entropy_mole

Molar entropy [J/kmol/K]. Returns an array of length n_points.

equilibrium_constants

Equilibrium constants in concentration units for all reactions. Returns an array of size n_reactions x n_points.

forward_rate_constants

Forward rate constants for all reactions. Units are a combination of kmol, m^3 and s, that depend on the rate expression for the reaction. Returns an array of size n_reactions x n_points.

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. Returns an array of size n_reactions x n_points.

g

Gibbs free energy [J/kg or J/kmol] depending on basis. Returns an array of length n_points.

gibbs_mass

Specific Gibbs free energy [J/kg]. Returns an array of length n_points.

gibbs_mole

Molar Gibbs free energy [J/kmol]. Returns an array of length n_points.

grid

FlameBase.grid(self) Array of grid point positions along the flame.

h

Enthalpy [J/kg or J/kmol] depending on basis. Returns an array of length n_points.

int_energy

Internal energy in [J/kg or J/kmol]. Returns an array of length n_points.

int_energy_mass

Specific internal energy [J/kg]. Returns an array of length n_points.

int_energy_mole

Molar internal energy [J/kmol]. Returns an array of length n_points.

isothermal_compressibility

Isothermal compressibility [1/Pa]. Returns an array of length n_points.

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. Returns an array of size n_species x n_points.

mix_diff_coeffs_mass

Mixture-averaged diffusion coefficients [m^2/s] relating the diffusive mass fluxes to gradients in the species mass fractions. Returns an array of size n_species x n_points.

mix_diff_coeffs_mole

Mixture-averaged diffusion coefficients [m^2/s] relating the molar diffusive fluxes to gradients in the species mole fractions. Returns an array of size n_species x n_points.

net_production_rates

Net production rates for each species. [kmol/m^3/s] for bulk phases or [kmol/m^2/s] for surface phases. Returns an array of size n_species x n_points.

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. Returns an array of size n_reactions x n_points.

partial_molar_cp

Array of species partial molar specific heat capacities at constant pressure [J/kmol/K]. Returns an array of size n_species x n_points.

partial_molar_enthalpies

Array of species partial molar enthalpies [J/kmol]. Returns an array of size n_species x n_points.

partial_molar_entropies

Array of species partial molar entropies [J/kmol/K]. Returns an array of size n_species x n_points.

partial_molar_int_energies

Array of species partial molar internal energies [J/kmol]. Returns an array of size n_species x n_points.

partial_molar_volumes

Array of species partial molar volumes [m^3/kmol]. Returns an array of size n_species x n_points.

reverse_rate_constants

Reverse rate constants for all reactions. Units are a combination of kmol, m^3 and s, that depend on the rate expression for the reaction. Returns an array of size n_reactions x n_points.

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. Returns an array of size n_reactions x n_points.

s

Entropy [J/kg/K or J/kmol/K] depending on basis. Returns an array of length n_points.

set_gas_state(self, point)
set_profile(self, component, locations, values)
set_refine_criteria(self, ratio=10.0, slope=0.8, curve=0.8, prune=0.0)
solution(self, component, point=None)
soret_enabled

FlameBase.soret_enabled(self)

standard_cp_R

Array of nondimensional species standard-state specific heat capacities at constant pressure at the current temperature and pressure. Returns an array of size n_species x n_points.

standard_enthalpies_RT

Array of nondimensional species standard-state enthalpies at the current temperature and pressure. Returns an array of size n_species x n_points.

standard_entropies_R

Array of nondimensional species standard-state entropies at the current temperature and pressure. Returns an array of size n_species x n_points.

standard_gibbs_RT

Array of nondimensional species standard-state Gibbs free energies at the current temperature and pressure. Returns an array of size n_species x n_points.

standard_int_energies_RT

Array of nondimensional species standard-state internal energies at the current temperature and pressure. Returns an array of size n_species x n_points.

thermal_conductivity

Thermal conductivity. [W/m/K]. Returns an array of length n_points.

thermal_diff_coeffs

Return a one-dimensional array of the species thermal diffusion coefficients [kg/m/s]. Returns an array of size n_species x n_points.

thermal_expansion_coeff

Thermal expansion coefficient [1/K]. Returns an array of length n_points.

transport_model

FlameBase.transport_model(self)

u

FlameBase.u(self)

Array containing the velocity [m/s] normal to the flame at each point.

viscosity

Viscosity [Pa-s]. Returns an array of length n_points.

volume

Specific volume [m^3/kg or m^3/kmol] depending on basis. Returns an array of length n_points.

volume_mass

Specific volume [m^3/kg]. Returns an array of length n_points.

volume_mole

Molar volume [m^3/kmol]. Returns an array of length n_points.

write_csv(self, filename, species='X', quiet=True)

Write the velocity, temperature, density, and species profiles to a CSV file.

Parameters:
  • filename – Output file name
  • species – Attribute to use obtaining species profiles, e.g. X for mole fractions or Y for mass fractions.