Thermochemistry and Transport#

Pyrometheus generates code to evaluate chemical source terms. These appear in the conservation equations of reacting flows. Here, we lay out the corresponding equations. We focus on a homogeneous adiabatic reactor for simplicify. Yet, the systems explained here can easily be adapted to other configurations (e.g., isochoric or inhomogeneous reactors).

Chemical Kinetics and Thermodynamics#

Our goal is to express, in as much detail, the chemical kinetics of reactive flows. We assume a homogeneous mixture of ideal gases evolves at constant pressure \(p\) and enthalpy \(h_{0}\). We characterize its chemical composition by the species mass fractions \(\boldsymbol{y} = \{ y_{i} \}_{i = 1}^{N}\). These evolve from an initial condition \(y_{i}(0) = y_{i}^{0}\) according to

\[\label{eq:species_conservation} \frac{dy_{i}}{dt} = S_{i} \equiv \frac{ W_{i}\dot{\omega}_{i} }{ \rho },\]

where \(S_{i}\) is the chemical source term of the \(i^{\mathrm{th}}\) species (in \(\mathrm{s}^{-1}\)), \(W_{i}\) its molecular weight (in \(\mathrm{kg/kmol})\) and \(\dot{\omega}_{i}\) its molar production rate (in \(\mathrm{kmol/m^{3}-s}\)). The mixture density \(\rho\) (in \(\mathrm{kg/m^{3}}\)) is obtained from

\[pW = \rho RT,\]

where

\[W = \sum_{i = 1}^{N}W_{i}y_{i}\]

is the mixture molecular weight, \(R\) the universal gas constant (in \(\mathrm{J/kmol\cdot K}\)), and \(T\) the temperature (in \(\mathrm{K}\)). We explain how to obtain the mixture temperature in detail in Section 1.2.

To evaluate [species_conservation], we need the net production rates \(\dot{\boldsymbol{\omega}} = \{ \dot{\omega}_{i} \}_{i = 1}^{N}\). These represent changes in composition due to chemical reactions

\[\label{eq:reactions} \sum_{\ell = 1}^{N}\nu_{i\ell}^{\prime}\mathcal{S}_{\ell} \rightleftharpoons \sum_{k = 1}^{N}\nu_{ik}^{\prime\prime}\mathcal{S}_{k},\qquad j = 1,\dots,M,\]

where \(\nu_{ij}^{\prime}\) and \(\nu_{ij}^{\prime\prime}\) are the forward and reverse stoichiometric coefficients of species \(\mathcal{S}_{i}\) in the \(j^{\mathrm{th}}\) reaction. Per [reactions], species \(\mathcal{S}_{i}\) can only be produced (or destroyed) by an amount \(\nu_{ij}^{\prime\prime}\) (or \(\nu_{ij}^{\prime}\)) in the \(j^{\mathrm{th}}\) reaction. Thus, \(\{ \dot{\omega}_{i} \}_{i = 1}^{N}\) are linear combinations of the reaction rates of progress \(R_{j}\),

\[\label{eq:production_rates} \dot{\omega}_{i} = \sum_{j = 1}^{M}\nu_{ij}R_{j},\qquad i = 1,\dots,N,\]

where \(\nu_{ij} = \nu_{ij}^{\prime\prime} - \nu_{ij}^{\prime}\) is the net stoichiometric coefficient of the \(i^{\mathrm{th}}\) species in the \(j^{\mathrm{j}}\) reaction. The rates of progress are given by the law of mass-action,

\[\label{eq:reaction_rates} R_{j} = k_{j}(T)\left[ \prod_{\ell = 1}^{N}\left(\frac{ \rho y_{\ell} }{ W_{\ell} }\right)^{\nu_{mj}^{\prime}} - \frac{1}{K_{j}(T)}\prod_{k = 1}^{N}\left(\frac{ \rho y_{k} }{ W_{k} }\right)^{\nu_{mj}^{\prime\prime}} \right],\qquad j = 1,\dots,M,\]

where \(k_{j}(T)\) is the rate coefficient of the \(j^{\mathrm{th}}\) reaction and \(K_{j}(T)\) its equilibrium constant. Depending on the reaction, the rate coefficient \(k_{j}(T)\) may take different forms (and even become a function of pressure). Its simplest form is the Arrhenius expression,

\[\label{eq:rate_coeff} k_{j}(T) = A_{j}T^{b_{j}}\exp\left({ -\frac{\theta_{a,j}}{T} }\right),\qquad j = 1,\dots,M\]

where \(A_{j}\) is the pre-exponential, \(b_{j}\) is the temperature exponent, and \(\theta_{a,j}\) is the activation temperature.

The equilibrium constant is evaluated through equilibrium thermodynamics

\[\label{eq:equil_constants} K_{j}(T) = \left( \frac{p_{0}}{RT} \right)^{\sum_{i = 0}^{\nu_{ij}}}\exp\left( -\sum_{i = 1}^{N}\frac{\nu_{ij}g_{i}(T)}{RT} \right),\qquad j = 1,\dots,M,\]

where \(p_{0} = 1\) \(\mathrm{atm}\) and

\[g_{i}(T) = h_{i}(T) - T\,s_{i}(T),\qquad i = 1,\dots,N\]

are the species Gibbs functions, with \(h_{i}\) and \(s_{i}\) the species enthalpies and entropies.

Species Thermodynamics#

Conservation of Energy#

To evaluate the rates of progress [reaction_rates], we need the temperature. Yet, we have defered any discussion on how to compute it from other state variables.

Transport Properties#

Pyrometheus-generated code provides routines to evaluate species and mixture transport properties. These follow most closely the Cantera implementation, which is based on polynomial fits to collision integrals. This approach is based on the kinetic theory of gases, for which a complete overview can be found in chapter 12 of [Kee_2003].

The viscosity of the \(n^{\mathrm{th}}\) species in the mixture is:

\[\mu_n = \sqrt{T} \left[\sum_{m = 0}^{4} a_{m, n}\, (\log\, T)^{m}\right]^2,\]

where the coefficients \(a_{m, n}\) are provided by Cantera. The viscosity of the mixture is then obtained via the mixture rule

\[\mu = \sum_{n = 1}^{N} \frac{X_n \mu_n}{\sum_{j = 1}^{N} X_j\Phi_{nj}}\]

where \(X_{n} = W Y_{n} / W_{(n)}\) is the mole fraction of species \(n\), and

\[\Phi_{nj} = \frac{1}{\sqrt{8}} \left( 1 + \frac{W_n}{W_j} \right)^{-\frac{1}{2}} \left( 1 + \left[ \frac{\mu_n}{\mu_j} \right]^{\frac{1}{2}} \left[ \frac{W_j}{W_n} \right]^{\frac{1}{4}} \right)^2.\]

The thermal conductivity of species \(n\) is

\[\lambda_n = \sqrt{T} \sum_{m = 0}^{0} b_{m, n}\, (\log\, T)^{m}.\]

The mixture viscosity is

\[\lambda = \frac{1}{2} \left( \sum_{n = 1}^{N} X_n \lambda_n + \frac{1}{\sum_{n = 1}^{N} \frac{X_n}{\lambda_n} } \right).\]

The binary mass diffusivities, in \(\frac{m^2}{s}\), for species $i$ and $j$ are

\[D_{i,j}(T) = \frac{T^{3/2}}{p} \sum_{m = 0}^{4}c_{i,j,m}\, (\log\, T)^m\]

The mixture-averaged diffusivity of species \(n\) is

\[\mathscr{D}_{n} = \frac{W - X_{(n)}W_{n}}{W}\left\lbrace \sum_{m \neq n}\frac{X_{m}}{D_{nm}} \right\rbrace^{-1}\]

This expression becomes singular for \(X_n = 1\) (for any \(n\), so \(\sum_{m \neq n} X_m/D_{nm} = 0\)). Thus, following Cantera, it only returns the mixture-averaged diffusivity if

\[\sum_{m \neq n} \frac{X_{m}}{D_nm} > 0,\]

and \(D_{nn}\) otherwise. The conditional is implemented using numpy.where() and, of course, it is difficult to satisfy in finite-precision calculations. It can lead to round-off errors in \(\mathscr{D}_{n}\) but, like Cantera, Pyrometheus does not attempt to correct this behavior to avoid the use of arbitrary thresholds.