© Springer Nature Switzerland AG 2020
E. Lindner et al. (eds.)Mathematical Modelling in Real Life ProblemsMathematics in Industry33https://doi.org/10.1007/978-3-030-50388-8_3

3. Uncertainty Quantification of Chemical Kinetic Reaction Rate Coefficients

É. Valkó1   and T. Turányi2
(1)
Institute of Mathematics, ELTE Eötvös Loránd University, Budapest, Hungary
(2)
Institute of Chemistry, ELTE Eötvös Loránd University, Budapest, Hungary
 
Keywords
Chemical kineticsRate coefficientsTemperature dependent uncertainty limitsUncertainty analysis

3.1 Introduction

In chemical kinetics, the reaction rate coefficients (also called rate constants) characterize the speed of a chemical reaction. The temperature dependence of the rate coefficients can be defined by Arrhenius parameters A, n, and E. Chemical kinetics databases for many elementary gas-phase reactions provide the recommended values of the Arrhenius parameters, the temperature range of their validity and the uncertainty of rate coefficient k defined by uncertainty parameter f. However, in the kinetics databases the uncertainty parameter f is usually considered to be temperature independent. The goal of this project is to calculate temperature dependent uncertainty limits of rate coefficients of elementary reactions in such a way that these limits are consistent with the temperature dependence of the rate coefficient. The project consists of three tasks:
  1. 1.

    Visualization of the available kinetics data and calculation of the f original uncertainty parameters of the investigated rate coefficient in selected temperature points.

     
  2. 2.

    Calculation of the f extreme uncertainty parameters in selected temperature points.

     
  3. 3.

    Determination of a continuous, temperature dependent uncertainty function $$f_{\text{prior}}\left (T\right )$$ , which is consistent with the rate coefficient-temperature function.

     
In accordance with the three tasks above, in Sects. 3.2, 3.3, and 3.4 all background information and definitions needed to solve the problems are provided. A possible solution is presented using the example of elementary reaction OH+H2=H2O+H. This chapter is based on the publications of Nagy et al. [1–3]. All program codes used at the preparation of the examples are freely available from Web site [2, 4].

3.2 Calculation of the Uncertainty Parameters f original in Selected Temperature Points

The temperature dependence of rate coefficient k of a chemical reaction is usually described by the modified Arrhenius equation
$$\displaystyle \begin{aligned} k(T)=A \{T\}^n \exp\left( - \frac{E}{RT}\right), \end{aligned} $$
(3.1)
where R is the universal gas constant, T is temperature in Kelvin and parameters A, n, E are the so-called Arrhenius parameters. In accordance with the recommendations [5], curly brackets are used to denote the numerical value of the enclosed physical quantity at the predefined units, which are in this case cm, K, s, mol. Introducing transformed parameters $$\kappa (T):= \ln \lbrace k(T)\rbrace ,\; \alpha := \ln \lbrace A\rbrace $$ and ε := E∕R, the linearized form of the modified Arrhenius equation is
$$\displaystyle \begin{aligned} \kappa(T)=\alpha +n\ln \{T\}-\varepsilon T^{-1}. {} \end{aligned} $$
(3.2)
To show the temperature dependence of κ, the logarithm of the rate coefficient is plotted as a function of the reciprocal of temperature. Such a κ(T −1) curve is called an Arrhenius curve.
Figure 3.1 shows the Arrhenius curves of measured and theoretically determined rate coefficients for reaction OH+H2=H2O+H. The corresponding Arrhenius parameters and their temperature intervals are available in the Supplementary Material of article [2].
../images/475466_1_En_3_Chapter/475466_1_En_3_Fig1_HTML.png
Fig. 3.1

The logarithm of measured and theoretically determined rate coefficients as a function of the reciprocal of temperature belonging to reaction OH+H2=H2O+H

The rate coefficient of an elementary reaction can be determined by various experimental methods. If several measurements were carried out in different laboratories (maybe using different methods) at similar temperatures, then the uncertainty of the rate coefficient can be well assessed at a given temperature or in a narrow temperature interval. If the uncertainty of a rate coefficient is determined from literature data at different temperatures, then these uncertainties can be very different from each other even at nearby temperatures. However, since the measured rate coefficients are interrelated by a common Arrhenius expression, therefore the uncertainties determined at different temperatures are also related. Taking into account the temperature dependence of the rate coefficient, the uncertainty at a given temperature cannot be high if it is low at nearby temperatures.

The next pages discuss the determination of an Arrhenius-equation-consistent uncertainty function from the uncertainties of a rate coefficient valid at given temperatures (or in given temperature intervals) and the features of the corresponding uncertainty domain of the Arrhenius parameters.

Definition 3.1
The temperature dependent uncertainty function $$f\left (T\right )$$ is defined by Baulch et al. [6]:
$$\displaystyle \begin{aligned} f=\log_{10}\left(k^0(T)/k_{\text{min}}(T)\right)=\log_{10} \left( k_{\text{max}}(T)/k^0(T) \right), \end{aligned} $$
(3.3)
where k 0(T) is the recommended value of the rate coefficient at temperature T, k min(T) and k max(T) are the extreme, but still not excludable, physically realistic values at given temperature.

This definition of uncertainty is related to the limits and does not necessarily have a probabilistic inference. According to Eq. (3.3), the upper and lower extreme values differ from the recommended value by a multiplication factor, which means that, on a logarithmic scale, the extreme values are located symmetrically around the recommended value. The uncertainty bounds for rate coefficient k can be converted to the uncertainty bounds of $$\kappa (T):= \ln \{k(T)\}$$ .

Definition 3.2

The uncertainty bound of the logarithm of the rate coefficient, κ(T), is $$\left [\kappa (T)-\frac {f(T)}{\ln 10} ,\kappa (T)+\frac {f(T)}{\ln 10} \right ]$$ .

Our aim is to define a temperature dependent uncertainty function f(T), and thus defining a symmetrical bound for the recommended rate coefficient κ 0(T). This bound will indicate the extreme, but physically realistic rate coefficient values at every temperature.

According to Eq. (3.3), the uncertainty of the investigated rate coefficient at a given temperature point, T i, can be defined by
$$\displaystyle \begin{aligned} \begin{array}{r l} f_{\text{original}}\left(T_i\right)&=\max\left( |\log_{10}k^0\left(T_i\right)-\log_{10}k_j\left(T_i\right)| \right) \\ j&=1,2,\ldots,m_i;\; i=1,2,\ldots, n_T \end{array} \end{aligned} $$
(3.4)
where k j is the j-th experimentally or theoretically determined value of the rate coefficient at temperature T i, m i is the number of all available values of the rate coefficient at temperature T i and n T is the number of temperature points.
Figure 3.2 shows the calculation of the $$f_{\text{original}}\left (T_i\right )$$ values for reaction OH+H2=H2O+H at T =  1000, 1100, 1200 and 1300 K. The investigated measurements and theoretical calculations for the values of the Arrhenius parameters and the corresponding temperature intervals are available in the Supplementary Material of [2].
../images/475466_1_En_3_Chapter/475466_1_En_3_Fig2_HTML.png
Fig. 3.2

The figure shows the Arrhenius curves of the rate coefficients that had been presented in Fig. 3.1, now in a narrower temperature interval (dashed lines). The mean value of the rate coefficient of reaction OH+H2=H2O+H was added (solid red line in the middle). The black circles and the sections with arrowheads (↔) show the extent of the $$f_{\text{original}} \left (T_i \right )$$ values at T =  1000, 1100, 1200 and 1300 K

Chemical kinetic databases contain suggested values of Arrhenius parameters A, n, E (or transformed Arrhenius parameters α, n, ε) and the temperature range where these values are relevant. Such information is available from the NIST Chemical Kinetics Database [7], the evaluations of Warnatz [8], Tsang et al. (see e.g. [9–11]), Baulch et al. (see e.g.[6, 12, 13]) and the review of Konnov [14].

Task 1

Collect all suggested values for the rate coefficient of reaction H+O2=O+OH and create a user-friendly program code to visualize them in an Arrhenius plot. Using Arrhenius parameters A = 2.07E+14 cm3mol−1s−1, n = −0.097, $$\frac {E}{R}= 7560$$ K (the mean rate expression recommended by Baulch et al. [6]), calculate uncertainty parameter $$f_{\text{original}}\left (T_i\right )$$ at temperature points T = 800, 900, …, 2700 K based on Eq. (3.4).

3.3 Calculation of the Uncertainty Values f extreme at Selected Temperature Points

After the calculation of the $$(T_i, f_{\text{original}}\left (T_i\right ))$$ point pairs, it is a natural idea fitting a polynomial to these points and using the fitted curve as a temperature dependent uncertainty function. However, the uncertainty bound defined by this function may contain physically not meaningful values of the rate coefficient. For this reason, a refined approach for the determination of the temperature-uncertainty point pairs is needed to obtain a physically relevant bound.

The procedure described here determines the uncertainty domain of Arrhenius parameters (p = (α, n, ε)T) from the uncertainty information for the rate coefficients. In several cases the temperature dependence of the rate coefficient can be described by two Arrhenius parameters (α, ε) or (α, n). In this case the third Arrhenius parameter is set to zero.

Assume that a central set of Arrhenius parameters p 0 is available and the symmetric uncertainty of the rate coefficient is estimated at several temperatures by uncertainty parameters $$f_{\text{original}}\left (T_i\right ), i = 1, \dots , n_T$$ . It is possible to generate all Arrhenius curves κ(T, p) that lie between the uncertainty limits, fulfilling the following 2n T inequalities:
$$\displaystyle \begin{aligned} -f_{\text{original}}(T_i) \leq \frac{\kappa(T_i;\mathbf{p})-\kappa(T_i;{\mathbf{p}}^0)}{\ln 10}\leq +f_{\text{original}}(T_i), \quad i=1,2,\dots,n_T. \end{aligned} $$
(3.5)
These curves are located symmetrically around the mean rate coefficient curve $$\kappa \left ( T;{\mathbf {p}}^0 \right )$$ , since Arrhenius equation (3.2) is a linear function of parameters α, n, ε, and equation (3.3) defines symmetric linear constraints. A systematic procedure is proposed here for determining the extreme Arrhenius curves, which touch either the lower or the upper uncertainty limit at least at 2 or 3 temperatures for the 2- and the 3-parameter cases, respectively, and also go within the upper and lower uncertainty limits at all other temperatures. Formally, these criteria correspond to Arrhenius functions that fulfil at least 2 or 3 equality relations in Eqs. (3.5) and for the remaining 2n T − 2 or 2n T − 3 cases, respectively, either the equality or the inequality is fulfilled. The minimum and maximum values of these curves at a given temperature define the edges of the band of all possible Arrhenius curves.
In the case of the 3-parameter Arrhenius expression, term $$n\ln {T}$$ usually has a smaller contribution to the temperature dependence of the rate coefficient than − ε∕T, since $$\ln \lbrace T \rbrace $$ changes more slowly than 1∕T at combustion temperatures (800 K–2700 K). The effect of a change in the temperature exponent n on the rate coefficient at high temperatures can be well compensated by adjusting the pre-exponential factor α, leading to a very strong anti-correlation between α and n in most determinations. This implies that values of n, which significantly deviate (i.e. by ± 10) from the central n 0, can also fulfil all the inequality requirements in Eq. (3.5) if the initial uncertainty limits are not too tight. Both theoretical considerations [15] and the typical range of values of n in kinetic databases [7] show that the temperature exponent n of elementary chemical reactions should take values of small negative or positive numbers. Therefore, we recommend confining the range of n values to a narrow (i.e. Δn = 2) symmetric interval around the central value n 0 when the band of possible Arrhenius curves is determined through finding extreme Arrhenius curves.
$$\displaystyle \begin{aligned} -\varDelta n \leq n-n^0 \leq +\varDelta n \end{aligned} $$
(3.6)
The extreme Arrhenius curves are those which fulfil at least 2 or 3 equality relations in Eqs. (3.5) and (3.6) for the two-parameter and the three-parameter cases, respectively. To determine the extreme Arrhenius curves, uncertainty values need to be known at least at 2 temperatures, since in the three-parameter case a constraint is given for parameter n.
Definition 3.3
The minimum and maximum values of the extreme Arrhenius curves (κ min(T)) and (κ max(T)) define new uncertainty limits at any temperature, which are symmetrically located around the mean $$\kappa \left (T;{\mathbf {p}}^0\right )$$ curve. These new limits, obtained from a set of uncertainty values f and a user-defined Δn, uniquely define a new, continuous uncertainty function $$f_{\text{extreme}}\left (T\right )$$ :
$$\displaystyle \begin{aligned} f_{\text{extreme}}\left(T\right)=\frac{\kappa\left(T;{\mathbf{p}}^0\right)-\kappa_{\text{min}}\left(T\right)}{\ln 10}\equiv\frac{\kappa_{\text{max}}\left(T\right)-\kappa\left(T;{\mathbf{p}}^0\right)}{\ln 10}. \end{aligned} $$
(3.7)
By definition, this Arrhenius-equation-consistent uncertainty $$f_{\text{extreme}}\left (T_i\right )$$ is always less than or equal to the original uncertainty $$ f_{\text{original}} \left (T_i\right )$$ , at every temperature T i (i = 1, .., n T).
Task 2

Create a program code to calculate $$f_{\text{extreme}}\left (T_i\right )$$ uncertainty values at temperature points T = 800, 900, …, 2700 K based on the calculated $$\left (T_i, f_{\text{original}}\left (T_i\right )\right )$$ points for elementary reaction H+O2=O+OH.

3.4 Calculation of Uncertainty Function $$f_{\text{prior}}\left (T\right )$$

The temperature dependent rate coefficient k(T) (and its natural logarithm κ(T)) can be considered as a random variable deduced from measurements and calculations. Therefore, a probabilistic meaning may be attributed to f extreme. According to this interpretation, if f extreme corresponds to 3 standard deviations (3σ) [16–21] or 2 standard deviations (2σ) [22–24] of the rate coefficient on a decimal logarithmic scale, the uncertainty parameter f can be converted [17] to the standard deviation of the natural logarithm of the rate coefficient (σ κ)) at a given temperature T:
$$\displaystyle \begin{aligned} \sigma_{\kappa}\left(T\right)=\frac{\ln 10}{\mu} f\left(T\right), \end{aligned} $$
(3.8)
where parameter μ is 3 or 2, respectively.

If the rate coefficients are considered as random variables, then Arrhenius parameters α, n and ε are also random variables, since these can be calculated from the random values of κ(T) at three given temperatures using the linearized Arrhenius equation (Eq. (3.2)). The joint probability density function of the Arrhenius parameters is independent of temperature. This means that all central moments are also independent of temperature, including their expected values ( $$\bar {\alpha },\; \bar {n},\; \bar {\varepsilon }$$ ), variances $$\left ( \sigma _\alpha ^2, \; \sigma _n^2,\; \sigma _\varepsilon ^2\right )$$ and correlations $$\left ( r_{\alpha n},\; r_{\alpha \varepsilon },\; r_{\varepsilon n}\right )$$ .

The following relation was deduced by Nagy [1] between the variance of κ(T) and the elements of the covariance matrix of the Arrhenius parameters ( Σ p):
$$\displaystyle \begin{aligned} \sigma^2_{\kappa}(T)&={\boldsymbol\Theta}^T (\boldsymbol\Sigma_{\mathbf{p}}) {\boldsymbol\Theta}=\\ &= \sigma^2_\alpha+\sigma_n^2 \ln ^2 T+\sigma_\varepsilon^2 T^{-2}+ 2r_{\alpha n}\sigma_{\alpha}\sigma_n\ln T+ \\ & -2r_{\alpha \varepsilon}\sigma_{\alpha}\sigma_{\varepsilon}T^{-1}-2r_{ n \varepsilon}\sigma_{n}\sigma_{\varepsilon} T^{-1}\ln T \end{aligned} $$
(3.9)

A method was proposed [1] for the determination of the covariance matrix of the Arrhenius parameters using Eqs. (3.8) and (3.9) from uncertainty parameter f of the rate coefficient at various temperatures. To determine the elements of the covariance matrix for the three-parameter Arrhenius expression, the uncertainty of the rate coefficient ( $$f_{\text{extreme}}\left (T_i\right )$$ ) has to be known at least at six different temperatures. In the (α, ε) and (α, n) two-parameter cases, the uncertainty of the corresponding Arrhenius parameters can be handled in a similar way and the uncertainty of the rate coefficient has to be known at least at three temperatures [1]. Having calculated the $$f_{\text{extreme}}\left (T_i\right )$$ values at temperature points T i, the $$\sigma _{\kappa }^2 \left (T_i\right )$$ values are also available. The four or six elements of the covariance matrix have to be determined by minimising the difference between the values $$\sigma _{\kappa }^2 \left (T_i\right )$$ and $$ \sigma ^2_\alpha +\sigma _n^2 \ln ^2 T_i+\sigma _\varepsilon ^2 T_i^{-2}+ 2r_{\alpha n}\sigma _{\alpha }\sigma _n\ln T_i-2r_{\alpha \varepsilon }\sigma _{\alpha }\sigma _{\varepsilon }T_i^{-1}-2r_{n \varepsilon }\sigma _{n}\sigma _{\varepsilon }T_i^{-1} \ln T_i$$ for all T i, (i = 1, …, n T).

Equations (3.8) and (3.9) provide a means for storing the $$f_{\text{extreme}}\left (T\right )$$ function in the form of the covariance matrix of Arrhenius parameters. The uncertainty parameter temperature function calculated from the covariance matrix can be used as a first guess (prior) and conservative estimation of the uncertainty of a rate coefficient and it can be used as a starting point of a more refined estimation [25] of the temperature dependent uncertainty of the rate coefficient.

Definition 3.4

The uncertainty function reconstructed from the covariance matrix is called prior uncertainty and denoted as $$f_{\text{prior}}\left (T\right )$$ .

In Eq. (3.8), the parameter μ defines the proportionality between the uncertainty parameter f and the standard deviation σ κ. When the uncertainty f prior is calculated via σ κ from the covariance matrix Σ p, the same parameter μ has to be used. This means that the value of μ is arbitrary in the storage of the f values in the covariance matrix, and the only important assumption here is that the uncertainty parameter f is proportional to the standard deviation of κ.

Figure 3.3 shows the calculated $$(T_i, f_{\text{original}}\left (T_i\right )$$ ) values and the calculated $$f_{\text{extreme}}\left (T\right )$$ function for reaction OH+H2=H2O+H. The red solid line in Fig. 3.3 represents the $$f_{\text{prior}}\left (T\right )$$ function, which was fitted to points $$(T_i,f_{\text{extreme}}\left (T_i\right )$$ ), where T i =  300 K, 400 K,…, 2500 K.
../images/475466_1_En_3_Chapter/475466_1_En_3_Fig3_HTML.png
Fig. 3.3

The calculated f original(T i) values (black circles), the extreme uncertainty function, $$f_{\text{extreme}}\left (T\right )$$ (blue dashed line) and the fitted f prior(T) function (red solid line) for reaction OH+H2=H2O+H at temperature points T = 300 K, 400 K, …, 2500 K

Task 3

Using the results of Task 2, based on Eqs. (3.8) and (3.9), determine the values of the covariance matrix and define the $$f_{\text{prior}}\left (T\right )$$ function for reaction H+O2=O+OH.

Acknowledgements

Project no. ED_18-1-2019-0030 (Application-specific highly reliable IT solutions) has been implemented with the support provided from the National Research, Development and Innovation Fund of Hungary, financed under the Thematic Excellence Programme funding scheme.