Reaction of H2 with O2 in Excited Electronic States: Reaction

Nov 27, 2017 - Comprehensive quantum chemical analysis with the use of the multireference state-averaged complete active space self-consistent field ...
2 downloads 4 Views 1MB Size
Article Cite This: J. Phys. Chem. A 2017, 121, 9599−9611

pubs.acs.org/JPCA

Reaction of H2 with O2 in Excited Electronic States: Reaction Pathways and Rate Constants Alexey V. Pelevkin,†,‡ Boris I. Loukhovitski,† and Alexander S. Sharipov*,† †

Central Institute of Aviation Motors, Moscow 111116, Russia Moscow Institute of Physics and Technology, Dolgoprudny, Moscow Region 141700, Russia



S Supporting Information *

ABSTRACT: Comprehensive quantum chemical analysis with the use of the multireference state-averaged complete active space self-consistent field approach was carried out to study the reactions of H2 with O2 in a1Δg, b1Σ+g , c1Σ−u , and A′3Δu electronically excited states. The energetically favorable reaction pathways and possible intersystem crossings have been revealed. The energy barriers were refined employing the extended multiconfiguration quasi-degenerate second-order perturbation theory. It has been shown that the interaction of O2(a1Δg) and O2(A′3Δu) with H2 occurs through the Habstraction process with relatively low activation barriers that resulted in the formation of the HO2 molecule in A″ and A′ electronic states, respectively. Meanwhile, molecular oxygen in singlet sigma states (b1Σ+g and c1Σ−u ) was proved to be nonreactive with respect to the molecular hydrogen. Appropriate rate constants for revealed reaction and quenching channels have been estimated using variational transition-state theory including corrections for the tunneling effect, possible nonadiabatic transitions, and anharmonicity of vibrations for transition states and reactants. It was demonstrated that the calculated reaction rate constant for the H2 + O2(a1Δg) process is in reasonable agreement with known experimental data. The Arrhenius approximations for these processes have been proposed for the temperature range T = 300−3000 K.



only to the reaction kinetics of the H2 + O2(a1Δg) process, including reactive31−37 and quenching38,39 channels. The rate constant for the H2+O2(b1Σ+g ) reaction was estimated in the past, to our best knowledge, based on semiempirical schemes only,36,40 though the collisional quenching of O2(b1Σ+g ) by H2 was studied experimentally elsewhere.38,41−45 In addition, the measured branching ratios for the formation of H2 + {O2(X3Σ−g ), O2(a1Δg), O2(b1Σ+g )} products in the course of the H + HO2 reaction were reported at room temperature.42,46,47 Nevertheless, there is a pressing need to investigate these reaction pathways involving singlet oxygen molecules in a1Δg and b1Σ+g states more thoroughly, especially with an allowance for possible nonadiabatic transitions between different potential energy surfaces (PESs).7,27,48 Meanwhile, the reactions of higher electronic states of O2 with H2 were considered neither theoretically nor experimentally until now, with the exception of some qualitative observations.49 These facts dictate the necessity of extensive quantum chemical analysis of possible kinetic processes in the H2 + O2(a1Δg), H2 + O2(b1Σ+g ), H2 + O2(c1Σ−u ), and H2 + O2(A′3Δu) reacting systems. The present work deals with precisely this subject.

INTRODUCTION In the past, a great deal of effort has been devoted to the study of elementary processes involving electronically excited oxygen molecules owing to their importance in atmospheric,1−4 combustion,5−10 biological,11−13 and electric discharge14−16 chemistry. The fact is that the electronic excitation of oxygen molecules can substantially modify their reactivity.7,12,13,17,18 Excited oxygen in different electronic states also participates in the physicochemical processes behind the strong shock waves in oxygen (air)19−21 and in the active media of lasers.22,23 The reactions of electronically excited oxygen in the a1Δg and 1 + b Σg states (excitation energy values Te are 0.98 and 1.63 eV, respectively24) with molecular hydrogen are assumed to be the critical chain initiation channels under the conditions of laserinduced and plasma-assisted combustion of H2-containing mixtures.6,7,25−28 The reactions with higher (Herzberg) electronic states of O2, such as c1Σ−u (Te = 4.05 eV29) and A′3Δu (Te = 4.20 eV29), can also be potentially important in the latter case.7,9,14 Note that the oxygen molecules in these states have rather high lifetimes for spontaneous radiative transitions: 3900, 12, 6, and 1 s for a1Δg, b1Σ+g , c1Σ−u , and A′3Δu states, respectively.30 This stimulates the research of reaction kinetics in the system comprising electronically excited oxygen and molecular hydrogen. However, the available kinetic data on these processes are rather scarce. In the past, relatively much attention was paid © 2017 American Chemical Society

Received: October 10, 2017 Revised: November 22, 2017 Published: November 27, 2017 9599

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A



COMPUTATIONAL DETAILS Exploration of Potential Energy Surfaces. In the present work, ab initio quantum chemical calculations for reaction of ground state H2 with an O2 molecule, lying in one of the five lowest electronic states (X3Σ−g , a1Δg, b1Σ+g , c1Σ−u , and A′3Δu), were conducted to explore the PESs and to analyze the kinetic processes that take place in this tetra-atomic system. In this connection, the state-averaged complete active space selfconsistent field (SA-CASSCF) method,50 making it possible to explore both the ground and excited PESs within a single frame, was chosen as a basic computational tool. In doing so, stateaveraged orbitals with equal weighting for the states under investigation were specifically employed.51,52 As far as the CASSCF level of theory does not treat the dynamic electron correlation, it has a tendency to overestimate the energy barriers of metathesis reactions.36,53 Therefore, in order to refine the values of activation barriers, we employed pointwise electronic energy calculations with the use of extended multiconfiguration quasi-degenerate second-order perturbation theory (XMCQDPT2).54 Thereby, the CASSCF wave function was taken as the zeroth-order function for calculating energy at the second order of perturbation theory. Following the related studies published in the past,33,37,55,56 the full valence active space, composed of 14 electrons and 10 orbitals [denoted as (14,10)], was applied to the system under study. Dunning’s correlation-consistent basis sets with diffuse functions (aug-cc-pVXZ, X = 2.4)57 were used throughout the work. All ab initio calculations were performed by using the Firefly QC program package,58 which is partially based on the GAMESS(US) source code.59 For all considered PESs, the geometry of the reactants (R), transition states (TS), intermediates (IM), and possible products of the reaction (P) was optimized at the CASSCF(14,10)/aug-cc-pVDZ level of theory. For each critical point, vibrational frequency analysis at the same level of theory was performed in order to check the number of imaginary vibrational frequencies and to assess the zero-point energy (ZPE) contribution to the total energy. It should be emphasized that in the present work the ZPE scaling was not applied to the computed CASSCF(14,10) vibrational frequencies. The special analysis showed (see the Supporting Information) that the incorporation of a single scaling factor for the system under study cannot provide systematic improvement for vibrational frequencies (calculated at the CASSCF(14,10)/aug-cc-pVDZ level of theory) compared to spectroscopic data.24,29,60 For each located TS, the corresponding minimum-energy path (MEP) was followed applying the fourth-order Runge− Kutta method61 (at the CASSCF(14,10)/aug-cc-pVDZ level of theory as well). Besides, to study possible intersystem crossings (ICs), the projections of the obtained MEPs on the adjacent PESs were obtained. In order to define the position of TSs more accurately, we employed the version of the well-known IRCMax approach62 in which one selects the maximum from the energy points in the vicinity of the TS along the low-level calculated MEP (CASSCF(14,10)/aug-cc-pVDZ), recalculated at the high level of theory (XMCQDPT2/aug-cc-pVQZ). The commonly applied notation for such calculations is IRCMax[XMCQDPT2/aug-cc-pVQZ]:[CASSCF(14,10)/aug-ccpVDZ]. For R, IM, and P structures, the refinement of energy values by using single-point XMCQDPT2 calculations was employed also (common notation for this type of calculations is

XMCQDPT2/aug-cc-pVQZ//CASSCF(14,10)/aug-ccpVDZ). Note that a similar computational scheme was successfully applied to the excited PESs of the H + O2 reaction system in the past.63 Calculation of Reaction Rate Constants. It is known that for the elementary chemical reactions with nonzero activation barrier, when a single potential maximum along the MEP on the electronic PES can be associated with the TS structure, separating the reactant and product PES regions, the rate constant can be estimated applying transition-state theory.64 In the present work, both conventional transition-state theory (CTST) and canonical variational theory (CVT) with the use of Eckart-type tunneling correction were applied in the same manner as in the previous works of one of the authors.63,65 Therefore, for estimation of the CTST rate constant, the energy values, calculated at the IRCMax[XMCQDPT2/aug-ccpVQZ]:[CASSCF(14,10)/aug-cc-pVDZ] level of theory, were used (this approximation for rate constants is denoted here and hereinafter as kCTST PT2/QZ). In order to evaluate the influence of variational effects on the reaction rate constant, the projected vibrational frequency technique66 was applied at the CASSCF(14,10)/aug-cc-pVDZ level of theory (just this level of theory was used for the MEP following) along the reaction path. Then, the reaction rate constant, estimated by CVT for the reaction profile determined at the CASSCF(14,10)/aug-cc-pVDZ level of theory (denoted as kCVT CAS/DZ), was compared with the rate constant calculated employing CTST with the same activation barrier (kCTST CAS/DZ). The resulting rate constant is governed, in this case, by the expression k(T ) =

CTST kPT2/QZ

× r (T )

r (T ) =

CVT k CAS/DZ (T ) CTST k CAS/DZ (T )

(1)

where the factor r(T) shows the extent of overprediction in the value of the rate constant caused by the use of nonvariational theory against the variational one. It should be noted that in the CVT calculations, the vibrational partition functions were calculated under the harmonic oscillator (HO) approximation. Such practice is often applied for rate constant evaluation elsewhere (especially if the reactants and TS do not contain the floppy vibrational modes).64,67,68 In this case, the cancellation of anharmonicity corrections between the TS and reactants (products) is supposed implicitly. However, neglect of the anharmonicity of the TS and reactants vibrations can potentially lead to appreciable errors in the predicted partition functions and, consequently, in the estimated rate constants.69,70 In order to allow for the effect of vibrational anharmonicity on the predicted rate coefficients, in the course of kCTST PT2/QZ estimations, the vibrational partition functions were calculated with the use of an anharmonic oscillator (AHO) approximation. In doing so, the whole spectrum of vibrational states was represented by the system of independent one-dimensional Morse oscillators (specifying the individual vibrational modes) with the corresponding depths of the potential well Emax,i. The resulting vibrational partition function was calculated by direct summation of the relative populations of vibrational levels for each mode. The potential well depth for each mode of the structure under study was estimated by means of the simple approach,71,72 allowing one to avoid extensive investigation of the PES for the molecules and TS under study. The idea is based on the seminal suggestion of Zavitsas73 that the bond 9600

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A

treatment of anharmonic effects on the vibrational partition function of small polyatomic molecules. The probabilities of nonadiabatic transitions between different PESs were estimated to a first approximation by the use of the simple and computationally inexpensive Landau− Zener model.78,79 In so doing, the probability of surface hopping Psh from one term (i) to another (j) in the IC point in a pseudo-one-dimensional frame is governed by the equation

dissociation energy of the polyatomic molecule correlates with the corresponding frequency of the vibrational mode. Thereby, the Emax,i values, specifying the energy of the highest vibrational level in the ith mode of the considered structure, can be expressed as Emax, i

a b Eat ⎛ ωei ⎞ ⎛ M i ⎞ ⎜ OH ⎟ ⎜ OH ⎟ = n − 1 ⎝ ωe ⎠ ⎝ M ⎠

(2)

⎧0 Etr < E IC Psh = ⎨ ⎩1 − PLZ Etr > E IC ⎪

where Eat is the atomization energy of the considered n-atomic structure, ωei and ωeOH are the characteristic vibrational frequencies of the ith mode of the structure and OH molecule (the properties of the OH molecule were chosen as the reference for the reactive system under study), respectively, Mi and MOH are their reduced masses, and a and b are the adjustable parameters. Note that the assumed character of Emax,i(n) dependence (∼1/(n − 1)) in refs 71 and 72 provides correct extrapolation to the limit of the dissociation energy for diatomic molecules. Meanwhile, the other possible types of correlation between the bond energy of the polyatomic molecule and corresponding vibrational frequency were discussed earlier;71,73−75 however, just eq 2 can ensure flexible approximation of potential well depths in a broad range of Emax,i values (10−2−102 eV).72 The values of a and b parameters were chosen on the basis of PES analysis along the trajectories of normal vibrations for the following ground-state molecules, relevant for the system under study: H2, O2, OH, H2O, HO2, H2O2, O3. Therefore, for stretching modes, the value of Emax,i was identified as the bond dissociation energy. The estimates of Emax,i were performed either (i) from thermochemical analysis or (ii) from quantum chemical calculations of corresponding TSs (at the UB97-2/ aug-cc-pvTZ level of theory). As a result, the following values of adjustable parameters were obtained: a = 1.05, b = 0.54. Details of these calculations can be found in the Supporting Information. Note that the choice of OH parameters (ωOH = e 5408 K and MOH = 1.07 g/mol) as the reference for anharmonicity of the H2 + O2 system is rather arbitrary, and in the case of a more representative training data set, they could also be used as the adjustable ones. In order to assess the performance of the applied variant of the AHO approximation, we compared its predictions for the vibrational partition function with the results of direct summation of the relative populations of vibrational levels for the ground state HO2 and H2O molecules. Table 1 lists the HO values of the QAHO vib /Qvib ratio for these molecules at different temperatures calculated by means of eq 2 and from the wellestablished data on vibrational energy levels.76,77 As is seen, the applied AHO methodology, despite some essential assumptions made, such as disregard for quartic anharmonicity and mode− mode coupling, demonstrates reasonable accuracy in the



⎛ 2πVij 2 PLZ = exp⎜⎜ − ⎝ ℏ|Fj − Fi|

HO2(A″) Wang et al.76

this work

Polyansky et al.77

300 600 900

1.001 1.013 1.032

1.001 1.012 1.027

1.000 1.004 1.013

1.000 1.002 1.006

(3)



RESULTS AND DISCUSSION Calculations of Electronic Terms. The PES diagram for the five lowest electronic states of the H2 + O2 system, obtained at the XMCQDPT2/aug-cc-pVQZ//CASSCF(14,10)/aug-ccpVDZ level of theory for reactants, products, and intermediates and by using the IRCMax[XMCQDPT2/aug-cc-pVQZ]: [CASSCF(14,10)/aug-cc-pVDZ] approach for TSs, is presented in Figure 1. Structures, frequencies of normal vibrations,

Figure 1. Considered electronic terms and reaction pathways of the H2 + O2 system calculated at the XMCQDPT2/aug-cc-pVQZ // CASSCF(14,10)/aug-cc-pVDZ level of theory with ZPE correction.

and rotational constants for saddle points and minima are given as Supporting Information. Note, the ZPE correction was calculated at the CASSCF(14,10)/aug-cc-pVDZ level of theory. Herein and hereafter, the following designations for R, TS, IM, and P critical points are used: the left superscript denotes the multiplicity, and the right subscript denotes the order number of a given critical point among these of the same multiplicity. It should be noted that when calculating the energies of twoparticle reagents and products the distances between separate particles were fixed at a value of about 4−5 Å. Table 2 lists the Te values for electronic states of O2 and HO2 molecules, obtained utilizing the computational scheme of the present work, in comparison with the reference spectroscopic data.24,29,60 One can conclude that though adequate values of

H2O

this work

⎞ ⎟ ⎟ ⎠

where ℏ is the Planck constant, Fi and Fj are the magnitudes of the gradients of two diabatic PESs at the crossing point, Etr is the kinetic energy of the reactants’ nuclear motion, EIC is the potential energy of the IC point with respect to reagents, μ is the reduced mass of the bond (mode) whose direction coincides with the reaction coordinate qr, and Vij is the spin− orbit coupling (SOC) matrix element between two interacting electronic terms i and j.

HO Table 1. Values of QAHO vib /Qvib for HO2(A″) and H2O Molecules, Estimated Using the Approach of the Present Work and Obtained Directly from the Data on Vibrational Energy Levels76,77

T, K

μ 2(Etr − E IC)

9601

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A

the term of H2 + O2(b1Σ+g ) does not correlate with any stable bimolecular products (similar to atomic hydrogen and HO2 in any electronic state), since the corresponding electronic states of HO2 are dissociable.88,89 The analysis exhibited that the corresponding MEP leads to the dissociation of the O2(b1Σ+g ) molecule [pathway (d)] (with a ∼3.4 eV energy barrier), rather than to the H atom abstraction. This fact is visualized in Figure 2, where the respective PES section based on the relaxed PES

Table 2. Values of Te, Obtained at the CASSCF/aug-ccpvDZ and the XMCQDPT2/aug-cc-pVQZ// CASSCF(14,10)/aug-cc-pVDZ Levels of Theory (Assuming ZPE Correction), As Compared to the Spectroscopic Data24,29,60 (in eV) state

CASSCF

XMCQDPT2

ref

O2(a1Δg) O2(b1Σ+g ) O2(c1Σ−u ) O2(A′3Δu) HO2(A′)

0.934 1.436 3.627 3.846 0.825

0.972 1.637 4.044 4.239 0.853

0.97624 1.62724 4.05029 4.20429 0.84560

Te can be obtained even at the CASSCF/aug-cc-pvDZ level of theory, the use of second-order perturbation theory and a relatively large basis set allows one, in the framework of the (14,10) active space, to improve coincidence with experiment for Te predictions substantially. It is notable that the use of a larger number of virtual orbitals for the O2 subsystem can result in deterioration of Te calculation accuracy for the lowest (a1Δg and b1Σ+g ) electronic states of molecular oxygen.80 Let us remark that the main features of the PES, correlating with the H2 + O2(X3Σ−g ) reacting system, were established using quantum chemical methods at the turn of the century33,81−83 (apart from the earliest semiempirical studies31). Refining of the PES for this reaction has been continued in recent work.35,56,84,85 Therefore, from Figure 1, one can see that the reaction of H2 with O2(X3Σ−g ) leads to the formation of H + HO2(A″) products through the saddle point 3TS0 (Cs symmetry) and local minimum 3IM0 (representing the weakly bound van der Waals H···HO2(A″) complex). The 3IM0 complex lies slightly lower on the lowest triplet PES than the products [let us denote this pathway by the (a) index]. This result qualitatively coincides with the findings of other researchers33,35,56,81,82 and with analysis of experimental data on the reaction rate constant.36,86 As for the H2 + O2(a1Δg) reacting system, in view of the 2fold orbital degeneracy of an oxygen molecule in the singlet delta state,13 the interaction of H2 with O2(a1Δg) can occur along two different PESs, via the transition states 1TS0 of C1 symmetry [pathway (b)] and via the 1TS1 of Cs symmetry [pathway (c)] with different energy, resulting in the formation of a HO2 molecule in different electronic states (A″ and A′) and a H atom. Note that two optical isomeric forms of the 1TS0 structure, separated by a low inversion barrier of ∼0.03 eV height, were detected. Similar to the 3TS0 structure, the 1TS1 one lies slightly lower than the corresponding products. Note also that the local minima 3IM0 and 1IM0 exist on the surfaces of pure electronic energy only. They practically vanish if one takes into account the energy of zero-point vibrations (as shown in Figure 1). This result correlates with the rather small value of the Lennard-Jones potential well depth for the H−HO2 pair εH−HO2 recommended for evaluation of transport properties87 (subject to the kind of combining rule, εH−HO2 varies in the range of 10−45 K). The values of energy barriers for pathways (b) and (c) are 1.78 and 2.32 eV, respectively. The energy barrier for pathway (b) is considerably lower than that for pathway (a) (2.48 eV), which agrees with previous expectations.36,40 It is remarkable that pathway (c), arising from the degeneracy of the H2 + O2(a1Δg) term, was not revealed in the previous relevant studies.33−36 As can be seen from Figure 1, the O2(b1Σ+g ) molecule interacts with H2 through the PES of repulsive nature. In fact,

Figure 2. Internal coordinates used for relaxed PES scan calculations and the resulting contour plot of the singlet PES, correlating with the H2 + O2(b1Σ+g ) term. The segment of the crossing seam between this PES and triplet the H2 + O2(A′3Δu) term is outlined by a curve with open squares.

scan is presented. In these calculations, the ROO and ROH distances (see the scheme in Figure 2) were varied, and for each fixed pair of OO and OH distance values, partial optimization of the rest of the coordinates (HH distance and all angles) was conducted. One can see that the H abstraction, leading to H + OH + O(3P) product formation (as well as O abstraction, leading to H 2O + O( 1D) formation), is energetically unfavorable as compared to O2(b1Σ+g ) dissociation. Note also that the possibility of formation of the OH radical and H2O molecule in the course of the H2 + O2(b1Σ+g ) reaction was found negligible in experiments at least, at atmospheric conditions.44 The analogous situation was observed for interaction of O2(c1Σ−u ) with H2 [pathway (e)]. In fact, as is seen, the molecular oxygen in the singlet sigma states proved to be nonreactive with respect to molecular hydrogen, contrary to some previous estimates.6,25 Meanwhile, the triplet term, corresponding to the highest considered electronic state of oxygen (A′3Δu), is reactive again. Interaction of H2 with O2(A′3Δu) proceeds via the 3TS1 of Cs symmetry, resulting in the formation of H and an excited HO2(A′) molecule [pathway (f)]. The value of the energy barrier for the (f) channel is rather small (∼0.17 eV). The reaction pathways, revealed upon analysis of the triplet and singlet PESs for the H2 + O2 reacting system, are listed in Table 3. The respective values of the activation energy in both Table 3. Revealed Reaction Pathways and Corresponding Values of the Activation Energy in Both Directions (E+a and E−a ) and the Reaction Enthalpy (ΔrH°) (in eV) label (a) (b) (c) (d) (e) (f) 9602

pathway H2 H2 H2 H2 H2 H2

+ + + + + +

O2(X3Σ−g ) → H + HO2(A″) O2(a1Δg) → H + HO2(A″) O2(a1Δg) → H + HO2(A′) O2(b1Σ+g ) → H2 + 2O(3P) O2(c1Σ−u ) → H2 + 2O(3P) O2(A′3Δu) → H + HO2(A′)

E+a

E−a

2.478 1.779 2.316

0.230

0.171

1.035

ΔrH° 2.521 1.549 2.402 3.402 0.995 −0.865

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A directions E+a and E−a and the reaction enthalpy ΔrH° are given there also. The pathways involving the four lowest states of O2 were proved to be highly endothermic. The interaction of H2 with the A′3Δu Herzberg state of oxygen is the only exothermic pathway among the ones considered. It is notable that the obtained ΔrH° value for pathway (a) is close to that calculated from the thermochemical data. Therefore, adopting the best currently available ΔfH298 ° values for the ground state HO2 from the literature (12.3−21.4 kJ/mol90−92), one may derive that the ΔrH°(a) values must lie within the range 2.4−2.5 eV. It would be interesting to compare the value of E+a for pathway (b) with the high-level calculations reported elsewhere. The comparison of results obtained in the course of quantum chemical calculations by different researchers33−35 and the calculations of the present work are presented in Table 4. One can see that all predictions, based on the ab initio

the known parameters of the reactants and products. Note that such approaches are frequently used for express evaluation of energy barriers and rate constants of the reactions with electronically excited molecules in the absence of accurate quantum chemical or experimental data.5,7,21,26,36,96 As seen, the semiempirical methods give a somewhat lower value of activation barrier (by ∼0.3 eV) for pathway (b) than the multireference and composite computational schemes. Note also that the same situation was also observed with respect to the other chain initiation reactions involving O2(a1Δg).36 Investigation of Possible Intersystem Crossings. In order to explore the possible ways of nonadiabatic transitions between different PESs under study, projections of the obtained pathways (a)−(f) on the neighboring PESs were calculated at the CASSCF(14,10)/aug-cc-pVDZ level of theory. Thereby, the mutual arrangement of the PESs along the corresponding MEPs was revealed. The calculations showed that the lowest triplet PES, correlating with the H2 + O2(X3Σ−g ) system, does not cross any one, at least along pathway (a) (see Figure 3a). An analogous situation was also observed for pathways (b) and (c) (see Figure 3b,c). Thus, the possibility of ICs cannot lead to the appearance of additional reaction routes for these terms. As for pathway (d), joining of the singlet PES correlating with the H2 + O2(b1Σ+g ) system with other terms was observed near the O2(b1Σ+g ) dissociation limit only (see Figure 3d). This implies that nonadiabatic transitions apparently are not of chemical significance in this case. A similar situation, as evidenced from Figure 3e, was observed for reaction pathway (e). It is remarkable that in the previous work of one of the authors,36 the nonadiabatic process

Table 4. Quantum Chemical Calculations and Semiempirical Estimates of the Activation Barrier for Pathway (b) method

E+a , eV

ref

MR(S)DCI//CASSCF CASPT2//CASSCF CCSD(full) G3 XMCQDPT2//CASSCF XMCQDPT2 BEBO MMVT

1.817 1.831 1.772 1.415 1.779 1.778 1.518 1.505

Filatov et al.33 Sayos et al.34 Mousavipour and Saheb35 Mousavipour and Saheb35 this worka this workb Mayer and Schieler40 Starik and Sharipov36

a

Basic methodology. bXMCQDPT2/aug-cc-pVQZ //XMCQDPT2/ aug-cc-pVDZ.

1

H 2 + O2 (b Σ+g ) → H + HO2 (A′)

multireference and coupled-cluster methods, are very close to each other [E+a (b) is within the range 1.77−1.83 eV], whereas the G3 composite technique provides an appreciably lower E+a (b) value (1.415 eV35). It appears that the main reason for this discrepancy results from the performance of the level of theory, adopted for primary geometry optimization in the G3 method [viz. RHF/6-31G(d)].93 In fact, the RHF/6-31G(d) level of theory fails to predict the proper geometry (and even the symmetry group) of the 1TS0 structure. Therefore, the 1TS0 structure, optimized at this level of theory, has C2v symmetry instead of C1, obtained at the CASSCF(14,10)/aug-cc-pVDZ level. Meanwhile, some other composite methods that are free from the indicated drawback, such as G494 and CBS-QB395 techniques (the B3LYP functional is used for the primary geometry optimization there), provide more consistent values of E+a (1.740 and 1.738 eV, respectively). In order to validate additionally the basic methodology employed throughout the present work, we tried to estimate the measure of uncertainty resulting from the fact that the highlevel (XMCQDPT2/aug-cc-pVQZ) refining calculations are performed at the geometry obtained at the relatively low (CASSCF/aug-cc-pVDZ) level of theory. In doing so, for pathway (b), the single-point XMCQDPT2/aug-cc-pVQZ calculations were performed at the geometry optimized at the XMCQDPT2/aug-cc-pVDZ level of theory. As is seen from Table 4, such methodology provides practically the same value of E+a (b) as the basic one (based on the IRCMax approach). Also shown in Table 4 are the semiempirical estimates of E+a (b) by means of the “bond energy−bond order” (BEBO)40 method and the modified method of vibronic terms (MMVT),36 which involve the model relationships between

(x)

was suggested. In view of the findings of the present work, this elementary process can be feasible only if there is a nonzero probability of surface hopping from the singlet H2 + O2(b1Σ+g ) term to another one, correlating with the H + HO2(A′) system (see Figure 1). However, this is not the case, at least along the MEP (d). Besides, in the course of the past experimental kinetic studies of the H/O2 system at T = 300 K,46,47 the pathway 1

H + HO2 (A″) → H 2 + O2 (b Σ+g )

(y)

O2(b1Σ+g )

was postulated as mainly the one responsible for the appearance. It was reported46,47 that a ∼10−3 fraction of the molecular oxygen formed in the reaction of H with HO2(A″) is in the b1Σ+g state. Meanwhile, the relative probability of this process with respect to the process backward to (b) was found to be as high as 0.019.47 Assuming that E−a (b) = 0.230 eV (see Table 3), one can derive from the elementary analysis based on the Arrhenius formula that the activation energy for the hypothetical process (y) can be as low as 0.33 eV. However, as evidenced from the Figures 1 and 3d, the term H2 + O2(b1Σ+g ) does not correlate with the H + HO2(A″) products. In principle, one can expect the presence of certain IC far from the explored MEP (d) region, making this reaction channel feasible. In order to check processes x and y, appropriate PES exploration was conducted in the regions outside of the MEP (d). In doing so, PES intersection with the H2 + O2(A′3Δu) term only was detected. The resulting segment of the crossing seam is marked at the contour plot of Figure 2. As is seen, the crossing between the PES, correlating with the H2 + O2(b1Σ+g ) 9603

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A

Figure 3. Electronic energy profiles as a function of the reaction coordinate for pathways (a)−(f) (solid curves) and energy profiles resulting from their projections on the other PESs (dashed curves).

system, and any other one (the H2 + O2(A′3Δu) term in this case) lies higher than the O2(b1Σg+) dissociation limit. Therefore, the possibility of the elementary processes x and y with reasonable values of the activation energy seems unlikely. As a possible reason for this discrepancy with experimental studies,46,47 one can assume that some additional processes, responsible for O2(b1Σg+) formation in the reacting H/O2 system, were not taken into account in the course of observed signal interpretation in these studies. For example, the process OH + HO2(A′) → O2(b1Σ+g ) + H2O, suggested on the basis of quantum chemical calculations recently,36 can be an efficient source of O2(b1Σg+), provided that the OH radicals are abundant in the reacting mixture. In this connection, it is remarkable that in the then46,47 current literature (e.g., see the measurements97) the branching ratio for two OH radicals’ production (H + HO2 = 2OH process) in the overall H + HO2 reaction was seriously underrated with respect to later data (by a factor of ∼3 at T = 300 K).91 Thus, the role of the OH + HO2(A′) → O2(b1Σ+g ) + H2O reaction could be unintentionally understated as the additional source of O2(b1Σ+g ) in the kinetic analysis.46,47 In addition, the process O + HO2(A″) → OH + O2(b1Σ+g ) with a supposedly6 low activation barrier is also possible in the

H/O2 system at experimental conditions.46,47 At last, the O2(b1Σ+g ) molecule can be formed simply by recombination of two O(3P) atoms (note that the oxygen atoms are abundant in the H/O2 system due to the H + HO2 = H2O + O(3P) and H + O2(X3Σ−g ) = OH + O(3P) reactions). As a result, the net rate of O2(b1Σ+g ) production from the mentioned sources could be erroneously ascribed to the single hypothetical process y. As for pathway (f), two ICs for the H2 + O2(A′3Δu) term were detected. Therefore, from Figure 3f, one can observe that this pathway intersects the H2 + O2(c1Σ−u ) (before passing through 1TS0) and H2 + O2(b1Σ+g ) (just after 1TS0) terms. The energies of these IC points are by 0.13 and 0.02 eV lower than that of the 1TS0 structure (at the CASSCF(14,10)/aug-ccpVDZ level of theory). This means that the H2 + O2(A′3Δu) system has a probability to leave twice the path (f), leading to H + HO2(A′) products. This fact should be taken into account when estimating the respective rate constants. Therefore, following the methodology applied earlier,63,79,98,99 one can associate these nonadiabatic transitions with collisional nonradiative quenching of an electronically excited O2(A′3Δu) molecule by H2. In order to obtain the locations of IC points more accurately (and to refine the potential energy of IC points with respect to 9604

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A

Thereby, the relatively large SOC between H2 + O2(A′3Δu) and H2 + O2(c1Σ−u ) terms was observed in the IC1 point, and one can expect intense nonadiabatic transitions between these PESs. At the same time, the very small Vij value in the IC2 point suggests that IC is, apparently, forbidden in this case. It is notable that the magnitudes of computed Vij matrix elements are a little sensitive to the accuracy of the IC points’ location. Therefore, the values of the same matrix elements, calculated for IC1 and IC2 points, obtained in course of CASSCF calculations (see Figure 3f), are equal to 100.3 and 0.1 cm−1, respectively. Further, the P1sh(Etr) and P2sh(Etr) dependences, specifying the influence of nonadiabatic transitions in the IC1 and IC2 points on the reaction kinetics in the system, were calculated with the use of eq 3 in a broad range of reagents’ kinetic energies. In so doing, the magnitudes of the gradients of diabatic PESs at the crossing points were taken from Figure 4. Reaction Kinetics. As it follows from the above analysis, the following reactive (Rn) and quenching (Qn) channels can occur in the reactive system under study

reagents EIC), the projections of pathway (f) on the neighboring PESs were recalculated using the XMCQDPT2 method. The resulting sections are depicted in Figure 4. One

Figure 4. Electronic energy profile as a function of the reaction coordinate for pathway (f) (solid curve) and the energy profiles resulting from their projections on the other PESs (dashed curves), recalculated at the XMCQDPT2/aug-cc-pVDZ level of theory. The IC points of interest are marked by circles.

can see that when applying the second-order perturbation theory the mutual arrangement of the different PESs remained the same as that in the case of CASSCF calculations. The first IC point (IC1), in the former case, is a bit closer to the TS region than that in the case of CASSCF projection (Ea − EIC = 0.02 eV), whereas the second IC point (IC2) is slightly shifted downhill toward the product region (Ea − EIC = 0.11 eV). In order to make it possible to evaluate the surface-hopping probabilities in these points, the SOC matrix elements between the interacting triplet and both singlet PESs were computed with the use of the effective nuclear charge approximation100 implemented in the Firefly program package.58 In doing so, a common set of molecular orbitals obtained with the stateaveraged CASSCF wave function calculations was used. The values of effective nuclear charges were chosen in line with the recommendations of Koseki et al.101 Note that the Vij value depends not only on the spin quantum numbers Si(j) of the interacting PESs (for IC1 and IC2 points we have Si = 1 and Sj = 0) but also on the projections of the spin on the preferred Z axis, MSi(j). Table 5 lists the calculated complex magnitudes of required matrix elements for different Ŝz operator eigenvalues.

IC1

IC2

−1 0 1

11.2 + i69.6 0 11.2 − i69.6

0 i0.8 0

S

Vij =

∑ Mi

Si =−Si

S

∑ Mj

Sj =−Sj

(R2)

H 2 + O2 (a1Δg ) → H + HO2 (A′)

(R3)

H 2 + O2 (A′3Δu) → H + HO2 (A′)

(R4)

H 2 + O2 (A′3Δu) → H 2 + O2 (c1Σ−u )

(Q1) (Q2)

Currently, the R1 channel is known to be a main chain initiation reaction in H2−O2(air) and CO−H2−O2(air) mixtures (let us remark that the role of different initiation reactions in hydrogen-containing mixtures was a subject of controversy for a long time, as discussed elsewhere8,82,103,104). The rate constant of this reaction is quite well defined over a wide temperature range, particularly by the extensive theoretical35,81,82 and precise experimental studies.82,105 At present, there is a reasonable consensus among the best data available in the literature on the R1 reaction rate constant. Therefore, the rate constant approximation proposed by Michael et al.,82 elaborated for the temperature range T = 400−2300 K, is supposed to be the most accurate:8,104,106 k(T) = 7.395 × T2.433 exp(26926/T) cm3 mol−1 s−1 (depicted in Figure 5). Therefore, there is no pressing need to refine the rate constant for reaction R1 in the present work. As for the reaction of H2 with singlet delta oxygen, supposed to be responsible for enhancement of chain processes in H2containing mixtures by means of electric discharge or resonance laser radiation,6,7,26 there are two different channels (R2 and R3) with different activation barriers, producing the HO2 molecule in the ground A″ and excited A′ electronic states, respectively. Evidently, due to the difference in activation energies (see Table 3), kR3(T) ≪ kR2(T) (at not very high temperature), and reaction R3 can be neglected during evaluation of the effective rate constant of the overall reaction H2 + O2(a1Δg) → products. The estimations for the temperature-dependent rate constants are also given in Figure 5. Meanwhile, the existence of two optical isomeric forms for the 1TS0 structure was taken into account in its

|⟨Si , MSi|Ĥ SO|Sj , MSj⟩|2

2 min(Si , Sj) + 1

H 2 + O2 (a1Δg ) → H + HO2 (A″)

1

According to Fedorov and co-workers,100,102 the square of the matrix element Vij was expressed via a bra and ket formalism as 2

(R1)

H 2 + O2 (A′3Δu) → H 2 + O2 (b Σ+g )

Table 5. Magnitude of the SOC Matrix Element in cm−1 between Triplet and Singlet Terms in the IC1 and IC2 Points for Different MSi Values (MSj = 0) MSi

H 2 + O2 (X 3Σ−g ) → H + HO2 (A″)

(4)

where Ĥ SO is the operator of SOC. The resulting Vij values for IC1 and IC2 points in Figure 4 are equal to 99.6 and 0.8 cm−1. 9605

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A

As seen from Table 6, the anharmonicity of the TS and reactants vibrations is also substantial for the reaction under study, especially at high temperature. The account of anharmonic effects leads to a notable increase in the channel R2 reaction rate constant (by a factor of 1.1−2.1), as long as these effects are more substantial for the TS structure than those for the reactants. Note that similar behavior was observed for other reactions under study (see the Supporting Information) and elsewhere.70 Note also that, although the applied variant of the AHO approximation provides reasonable accuracy for small polyatomic molecules (see Table 1), the vibrational partition function for the TS apparently can still be underestimated due to neglect of mode−mode coupling.69,109 Thereby, the reported values of an anharmonic correction factor for the rate constants can be treated as the lower estimate. It would be interesting to compare the kR2(T) dependence, obtained in the present work, with the data available elsewhere. Therefore, as far back as 1989, Basevich and Belyaev derived the rate constant for the channel R2 process in the temperature range of T = 300−800 K from fitting the complex mechanism to describe the measurements of the laminar flame velocity in the hydrogen/oxygen mixture with singlet oxygen addition.32 Their rate constant is also presented in Figure 5. One can observe that our calculated kR2 rate constant is somewhat smaller than that obtained from the kinetic modeling of Basevich and Belyaev32 (by a factor of 2 at T = 800 K). Note also that Basevich and Belyaev conceded that in the course of their evaluation the upper estimate of kR2(T) was provided. In this connection, the coincidence observed in Figure 5 is thought to be reasonable. The estimate of kR2(T), performed based on a semiempirical MMVT scheme,36 is also depicted here. As was mentioned above, the MMVT method underestimates the activation barrier for pathway (b). As a result, the MMVT-based rate constant is considerably larger than both the data extracted from the laminar flame experiments32 and the calculations of the present work. In addition, the branching ratios for the formation of molecular oxygen in different electronic states in the course of the H + HO2(A″) reaction were measured at room temperature by Washida et al.47 Washida et al. found that the k−R2/k−R1 ratio is equal to 0.025 ± 0.005 at T = 300 K. Afterward, this result was confirmed experimentally by Michelangeli et al.42 (k−R2/ k−R1 = 0.028). In order to compare the rate constant for the −R2 channel, calculated in the present work, with these measurements, for k−R1(T = 300 K), we adopted a theoretical estimation for the −R1 channel, based on the CCSD(full) calculated barrier35 (the approximation of Michael et al.82 is invalid for this temperature). In order to estimate the uncertainty range for k−R1 at this temperature, we also applied the other k−R1(T = 300 K) values reported elsewhere.81,82,91 Thus, we assumed that k−R1(T = 300 K) = 4.22 × 1012 ± 8.0 × 1011 cm3 mol−1 s−1. Combining these k−R1 and k−R2/k−R1 values, we can obtain the sought rate constant along with the uncertainty range: k−R2(T = 300 K) = 1.1 × 1011 ± 3 × 1010 cm3 mol−1 s−1. Figure 6 shows the comparison of present calculations for k−R2 with other data.32,35,47 One can see that the calculated k−R1(T) dependence coincides rather well with both the experimental data by Washida et al.47 (especially considering the estimated uncertainty) and the values derived from the laminar flame velocity experiments.32 The rate constant for the −R2 channel, estimated on the basis of

Figure 5. Temperature-dependent rate constants for the reaction channels R2 and R3, estimated in the present work, as well as the data for reactions R1 and R2 reported elsewhere.32,36,82

partition function (which is equivalent to the 2-fold reaction path degeneracy107) when calculating the rate constant for channel R2. One can see that even at T = 3000 K, kR2/kR3 ≈ 10. However, the fact that in the course of the reaction of H2 with singlet delta oxygen the excited HO2(A′) can arise is quite interesting from the spectroscopic diagnostics viewpoint. In addition, the reaction backward to eq R3, (−R3), seems to be important with respect to the kinetics of electronically excited HO2 (let us remark that the excited HO2 is believed to be an important actor in the H2/O2 mixture under nonequilibrium conditions27,48,108). The role of various factors, incorporated in the calculated reaction rate constant of channel R2, such as variational r(T), tunneling Γ(T), and anharmonic corrections, are presented for different temperatures in Table 6. One can see that the omitting 3 −1 −1 Table 6. Rate Coefficient kCTST s ) for the PT2/QZ (in cm mol R2 Reaction Channel Computed with the Use of CTST as Well as Values of r(T), Γ(T), and Anharmonic Correction Factors and the Total Rate Constant ktot(T) at Different Temperatures

T, K 300 500 1000 1500 2000 2500 3000

kCTST PT2/QZ 1.10 9.41 1.39 2.36 1.19 1.42 7.96

× × × × × × ×

10−16 10−05 105 108 1010 1011 1011

r(T)

Γ(T)

AHO corr.

0.41 0.44 0.42 0.39 0.41 0.38 0.37

203.46 11.65 2.19 1.48 1.27 1.17 1.12

1.09 1.25 1.72 2.00 2.08 2.04 1.93

ktot(T) 9.90 6.03 2.20 2.72 1.27 1.30 6.34

× × × × × × ×

10−15 10−4 105 108 1010 1011 1011

variational analysis can lead to the considerable overestimation of the channel R2 reaction rate constant (by a factor of ∼2.5). Meanwhile, neglect of tunneling correction can lead to serious underestimation of the reaction rate constant at low temperature (by 2 orders of magnitude at T = 300 K). Regarding the accuracy of the tunneling factor, the method applied in the present work is based on the assumption that the reaction and tunneling paths coincide. This leads to underestimation of the Γ value as long as the corner cutting effects are ignored. Therefore, we do not recommend using the rate constants obtained in the present work at low temperatures, where Γ is far from unity (T < 300 K). Certainly, more advanced and computationally demanding methods for the tunneling effect, based on the inclusion of divergence between the tunneling and the reaction paths,64,69 are desirable for application at such low temperatures. 9606

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A

presented in Figure 7. One can see that, despite a high probability of nonadiabatic transitions in the IC1 point, the

Figure 6. Temperature-dependent rate constant for the reaction channel −R2: calculations of the present work (solid curve) and those of Mousavipour and Saheb35 (dotted curve), measurements47 (with uncertainty bars), and estimations32 (symbols). The k−R1(T) rate constant reported by Michael et al.82 is also presented for comparison (dashed curve).

Figure 7. Temperature-dependent rate constant for the R4, Q1, and Q2 reaction channels, estimated in the present work.

resulting channel R4 rate constant is still rather high mainly due to a low activation barrier of pathway (f). Let us proceed to the nonadiabatic processes Q1 and Q2 that, as was mentioned above, can be attributed to the collisional quenching of an electronically excited O2(A′3Δu) molecule by H2. It should be noted that extensive analysis based on the surface-hopping trajectory calculations112,113 or accurate quantum dynamics114,115 is a preferable way to provide accurate estimates for the rate constants of electronic energy transfer processes. However, as a first approximation, the channel Q1 and Q2 rate constants can be estimated by using the MEP leading to O2(A′3Δu) quenching by H2, depicted in Figure 4. Therefore, the TS for the channel Q1 process is associated with the region of the IC1 point. However, because the IC1 point is located close to the saddle point 3TS1, one can suppose that the partition function of the system at the IC1 point is roughly equal to the partition function of the TS. Therefore, the temperature-dependent rate constant of the channel Q1 process also can be approximately expressed, following refs 79 and 111, in terms of modification of k0R4(T) dependence by the formulas

quantum chemical calculations by Mousavipour and Saheb,35 is also presented there. Let us proceed to the reactions with higher electronic states of O2. As was noted above, the interaction of H2 with molecular oxygen in the singlet sigma states O2(b1Σ+g ) and O2(c1Σ−u ) is not reactive. The respective MEPs [(d) and (e)], in fact, correspond to simple thermal dissociation of oxygen in b1Σ+g and c1Σ−u states 1

O2 (b Σ+g ) + M → O(3P) + O(3P) + M O2 (c1Σ−u ) + M → O(3P) + O(3P) + M

where M = H2. It is interesting that low reactivity of O2(b1Σ+g ) with respect to different substances was also observed in solution.18 However, the reaction of H2 with molecular oxygen in the highest considered electronic state of oxygen (A′3Δu), due to detected ICs, splits into three channels: reactive (R4) and quenching (Q1, Q2) ones. If one neglects the possibility of nonadiabatic transitions, the rate constant of the reaction R4 process, estimated with the use of the present methodology, can be expressed as follows: k0R4(T) = 1.487 × 107T2.042 exp(−346.2/T) cm3 mol−1 s−1. However, the presence of IC1 and IC2 points results in certain probabilities for the reacting system to leave the reaction path (f) in these points. Therefore, the Maxwell-averaged transition probability in the ICi point pi (i = 1, 2) for the molecules capable of overcoming the activation barrier Ea in reaction path (f), in compliance with refs 110 and 111, can be expressed as ∞

pi (T ) =

∫E Pshi (Etr) exp(−Etr /T ) a



∫E exp(−Etr /T ) a

0 k Q1(T ) ≈ kR4 (T )

1 ∞ dE tr dE / exp(− Etr /T ) Ttr T 0



(7)

q1(T ) =

∫0



1 Psh (Etr) exp(− Etr /T )

dEtr / T

∫0



exp(− Etr /T )

dEtr T (8)

where q1(T) is the Maxwell-averaged probability to overcome the EIC1 barrier and to pass through the IC1 point at the same time and the denominator of the quotient given by eq 7 represents the thermally averaged classical probability to overcome the Ea barrier. The rate constant of channel Q2, in turn, can be easily expressed as

(5)

Notice that the denominator of the quotient given by eq 5 equals exp(−Ea/T). Thus, the rate constant of reactive channel R4 can be estimated as 0 kR4(T ) ≈ kR4 (T )(1 − p1 (T ))(1 − p2 (T ))

∫E exp(−Etr /T ) a

dE tr T

dE tr T

q (T )



0 k Q2(T ) ≈ kR4 (T )(1 − p1 (T ))p2 (T )

(9)

The obtained rate constant temperature dependences are also plotted in Figure 7. As seen, the rate of O2(A′3Δu) quenching by H2 via channel Q1 prevails over the rate of the reactive channel R4 in the whole temperature range (by 1−2 orders of magnitude). This is caused mainly by the large SOC between corresponding PESs in the IC1 point. Meanwhile, O2(A′3Δu)

(6)

The calculations revealed that at T = 300−3000 K the p1 probability lies in the range of 0.98−0.83, whereas p2 in only about ∼3 × 10−5. The resulting rate constant dependence is 9607

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A

Table 7. Coefficients of Arrhenius Approximation Obtained for the Rate Constants of Reactions under Study (k in cm3 mol−1 s−1) k−

k+ elementary process R2 R3 R4 Q1 Q2

A 4.929 1.600 4.005 1.630 8.381

× × × ×

103 106 104 108

n

Ea, K

3.084 2.402 2.566 1.723 2.301

17513.0 25957.3 709.2 390.0 779.8

× × × ×

104 107 105 108

n

Ea, K

2.488 1.607 2.313 1.707 2.545

−5.1 −464.9 11199.4 2623.7 31024.6

precise trajectory calculations should be done to estimate them more accurately.

quenching via the IC2 point, located after the TS, is practically negligible. As the final result, all obtained temperature-dependent rate constants for channels R2−R4, Q1, and Q2 were approximated over the temperature range T = 300−3000 K by Arrhenius formula k(T) = ATn exp(−Ea/T) (cm3 mol−1 s−1). The corresponding coefficients are listed in Table 7.



A 5.294 6.705 2.764 1.630 2.589



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpca.7b09964. Tables comprising the geometrical structure and properties of computed critical points, their Cartesian coordinates, comparison of calculated vibrational frequencies with the spectroscopic data, data used to obtain a and b parameters in eq 2, and values of tunneling and anharmonic correction factors for the reactions considered (PDF)

CONCLUSIONS

Ab initio multireference calculations within the framework of the state-averaged complete active space self-consistent field approach were conducted to examine the PESs and to study the reaction kinetics in the H2 + {O2(a1Δg), O2(b1Σ+g ), O2(c1Σ−u ), O2(A′3Δu)} systems. The MEPs were found, and the energy values of critical points were refined employing second-order perturbation theory. Special attention was paid to the search for possible ICs between the considered PESs. It was revealed that the interaction of O2(a1Δg) and O2(A′3Δu) with H2 occurs through H-abstraction with relatively low activation barriers that results in the formation of the HO2 molecule in one of the electronic states (A″ and A′, respectively). Hence, the excitation of molecular oxygen to these states by means of electric discharge or resonance laser radiation allows one, in principle, to accelerate the chain initiation processes in H2-containing combustible mixtures, though the efficiency of O2(A′3Δu) as a promoter of chain initiation is reduced due to the existence of a nonadiabatic channel leading to its fast nonradiative quenching. At the same time, the molecular oxygen in singlet sigma states (b1Σ+g and c1Σ−u ) was proved to be nonreactive with respect to the molecular hydrogen, i.e., despite expectations, the excited O2(b1Σ+g ) and O2(c1Σ−u ) molecules cannot enhance the chain initiation processes in H2/O2 mixture under the conditions of laser-induced and plasma-assisted combustion. The appropriate temperature-dependent reaction and quenching rate constants were determined based on variational transition-state theory, taking into account the tunneling effect, anharmonicity of transitions states and reactants’ vibrations, and existing ICs. It was demonstrated that the rate constant for the reactive channel H2 + O2(a1Δg), obtained theoretically in the present work, is in reasonable agreement with known experimental data. The three-parameter Arrhenius approximations for both reactive and quenching channels were elaborated for the temperature range T = 300−3000 K. Let us remark that the methodology applied in the present work to treat nonadiabatic transitions allows one to make only rough estimates of the rate constant of quenching processes, and the temperature dependences of these rate constants are open for revision in future studies. Specifically, more accurate determination of the PES’s topology, application of more sophisticated models of the surface-hopping probability, and



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Boris I. Loukhovitski: 0000-0003-4092-8914 Alexander S. Sharipov: 0000-0002-1538-8672 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was supported by the Russian Science Foundation (Project No. 16-19-00111) and by the Russian Foundation for Basic Research (Grant No. 17-01-00810) for part of the methodology of reaction rate constant evaluation for the reactions with electronically excited species.



REFERENCES

(1) Slanger, T. G.; Copeland, R. A. Energetic Oxygen in The Upper Atmosphere and The Laboratory. Chem. Rev. 2003, 103, 4731−4765. (2) Yankovsky, V. A.; Manuilova, R. O. Model of Daytime Emissions of Electronically-Vibrationally Excited Products of O3 and O2 Photolysis: Application to Ozone Retrieval. Ann. Geophys. 2006, 24, 2823−2839. (3) Krasnopolsky, V. A. Excitation of The Oxygen Nightglow on The Terrestrial Planets. Planet. Space Sci. 2011, 59, 754−766. (4) Kirillov, A. S. Influence of Electronically Excited N2 and O2 on Vibrational Kinetics of These Molecules in The Lower Thermosphere and Mesosphere during Auroral Electron Precipitation. J. Atmos. Sol.Terr. Phys. 2012, 81−82, 9−19. (5) Ombrello, T.; Won, S. H.; Ju, Y.; Williams, S. Flame Propagation Enhancement by Plasma Excitation of Oxygen. Part II: Effects of O2(a1Δg). Combust. Flame 2010, 157, 1916−1928. (6) Starik, A. M.; Bezgin, L. V.; Kopchenov, V. I.; Loukhovitski, B. I.; Sharipov, A. S.; Titova, N. S. Numerical Study of The Enhancement of Combustion Performance in a Scramjet Combustor due to Injection of Electric-Discharge-Activated Oxygen Molecules. Plasma Sources Sci. Technol. 2013, 22, 065007. 9608

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A

(29) Slanger, T. G.; Cosby, P. C. O2 Spectroscopy below 5.1 eV. J. Phys. Chem. 1988, 92, 267−282. (30) Klotz, R.; Peyerimhoff, S. D. Theoretical Study of the Intensity of the Spin- or Dipole Forbidden Transitions between The c1Σ−u , A′3Δu, A3Σ+u and X3Σ−g , a1Δg, b1Σ+g States in O2. Mol. Phys. 1986, 57, 573−594. (31) Minaev, B. F. Quantum-Chemical Investigation of The Mechanism of Chain-Reaction Intiations during Hydrogen Combustion. Khim. Fiz. [Sov. J. Chem. Phys.] (in Russian) 1984, 3, 983−987. (32) Basevich, V. Y.; Belyaev, A. A. Evaluation of Hydrogen-Oxygen Flame Velocity Increase at Singlet Oxygen Addition. Chem. Phys. Rep. 1989, 8, 1124−1127. (33) Filatov, M.; Reckien, W.; Peyerimhoff, S. D.; Shaik, S. What Are The Reasons for The Kinetic Stability of a Mixture of H2 and O2? J. Phys. Chem. A 2000, 104, 12014−12020. (34) Sayos, R.; Oliva, C.; Gonzalez, M. Ab Initio CASPT2//CASSCF Study of The O(1D)+H2O(X1A1) Reaction. J. Chem. Phys. 2001, 115, 8828−8837. (35) Mousavipour, S. H.; Saheb, V. Theoretical Study on The Kinetic and Mechanism of H+HO2 Reaction. Bull. Chem. Soc. Jpn. 2007, 80, 1901−1913. (36) Starik, A.; Sharipov, A. Theoretical Analysis of Reaction Kinetics with Singlet Oxygen Molecules. Phys. Chem. Chem. Phys. 2011, 13, 16424−16436. (37) Coelho, D. V.; Brandao, J. A Full Dimensional Potential for H2O2 (X1A) Covering All Dissociation Channels. Phys. Chem. Chem. Phys. 2017, 19, 1378−1388. (38) Becker, K. H.; Groth, W.; Schurath, U. The Quenching of Metastable O2(1Δg) and O2(1Σ+g ) Molecules. Chem. Phys. Lett. 1971, 8, 259−262. (39) Du, S.; Leng, J.; Wang, J.; Yang, H.; Sha, G.; Zhang, C. Study on Collisional Deactivation of O2(1Δg) by H2 and D2. Chem. Phys. 2011, 383, 83−85. (40) Mayer, S. W.; Schieler, L. Activation Energies and Rate Constants Computed for Reactions of Oxygen with Hydrocarbons. J. Phys. Chem. 1968, 72, 2628−2631. (41) Stuhl, F.; Welge, K. H. Deactivation of O(1S) and O2(b1Σ+g ). Can. J. Chem. 1969, 47, 1870−1877. (42) Michelangeli, D. V.; Choo, K. Y.; Leu, M. T. Yields of O2(b1Σ+g ) and O2(a1Δg) in the H+O2 Reaction System, and The Quenching of O2(b1Σ+g ) by Atomic Hydrogen. Int. J. Chem. Kinet. 1988, 20, 915− 938. (43) Borrell, P.; Richards, D. S. Quenching of Singlet Molecular Oxygen, O2(a1Δg) and O2(b1Σ+g ), by H2, D2, HCl and HBr. J. Chem. Soc., Faraday Trans. 2 1989, 85, 1401−1411. (44) Talukdar, R. K.; Dunlea, E. J.; Brown, S. S.; Daniel, J. S.; Ravishankara, A. R. Kinetics of O2(b1Σ+g ) Reaction with H2 and an Upper Limit for OH Production. J. Phys. Chem. A 2002, 106, 8461− 8470. (45) Kozlov, D. N.; Kobtsev, V. D.; Stel’makh, O. M.; Smirnov, V. V. Study of Collisional Deactivation of O2(b1Σg+) Molecules in a Hydrogen-Oxygen Mixture at High Temperatures Using LaserInduced Gratings. J. Exp. Theor. Phys. 2013, 117, 36−47. (46) Hislop, R. J.; Wayne, R. P. Production of O2(b1Σ+g ) in The H + O2 System. J. Chem. Soc., Faraday Trans. 2 1977, 73, 506−516. (47) Washida, N.; Akimoto, H.; Okuda, M. Formation of Singlet State Molecular Oxygen in The Reaction of H + O2. J. Phys. Chem. 1978, 82, 18−21. (48) Chukalovsky, A. A.; Klopovsky, K. S.; Palov, A. P.; Mankelevich, Y. A.; Rakhimova, T. V. Reaction of Hydrogen Atoms with Singlet Delta Oxygen (O2(a1Δg)). Is Everything Completely Clear? J. Phys. D: Appl. Phys. 2016, 49, 485202. (49) Harteck, P.; Reeves, R. R., Jr. Formation and Reactions of the Excited O2(A3Σ+u ) Molecules. Discuss. Faraday Soc. 1964, 37, 82−86. (50) Roos, B. O. In Advances in Chemical Physics: Ab Initio Methods in Quantum Chemistry Part 2; Lawley, K. P., Ed.; John Wiley & Sons, Inc.: Hoboken, NJ, 1987; Vol. 69.

(7) Starik, A. M.; Loukhovitski, B. I.; Sharipov, A. S.; Titova, N. S. Physics and Chemistry of The Influence of Excited Molecules on Combustion Enhancement. Philos. Trans. R. Soc., A 2015, 373, 20140341. (8) Konnov, A. A. On The Role of Excited Species in Hydrogen Combustion. Combust. Flame 2015, 162, 3755−3772. (9) Castela, M.; Fiorina, B.; Coussement, A.; Gicquel, O.; Darabiha, N.; Laux, C. O. Modelling The Impact of Non-equilibrium Discharges on Reactive Mixtures for Simulations of Plasma-Assisted Ignition in Turbulent Flows. Combust. Flame 2016, 166, 133−147. (10) Larsson, A.; Zettervall, N.; Hurtig, T.; Nilsson, E. J. K.; Ehn, A.; Petersson, P.; Alden, M.; Larfeldt, J.; Fureby, C. Skeletal Methane-Air Reaction Mechanism for Large Eddy Simulation of Turbulent Microwave-Assisted Combustion. Energy Fuels 2017, 31, 1904−1926. (11) Minaev, B. F.; Agren, H. Spin-Orbit Coupling Induced Chemical Reactivity and Spin-Catalysis Phenomena. Collect. Czech. Chem. Commun. 1995, 60, 339−371. (12) Paterson, M. J.; Christiansen, O.; Jensen, F.; Ogilby, P. R. Overview of Theoretical and Computational Methods Applied to the Oxygen-Organic Molecule Photosystem. Photochem. Photobiol. 2006, 82, 1136−1160. (13) Minaev, B. F.; Murugan, N. A.; Agren, H. Dioxygen Spectra and Bioactivation. Int. J. Quantum Chem. 2013, 113, 1847−1867. (14) Gudmundsson, J. T. Recombination and Detachment in Oxygen Discharges: The Role of Metastable Oxygen Molecules. J. Phys. D: Appl. Phys. 2004, 37, 2073−2081. (15) Ionin, A. A.; Kochetov, I. V.; Napartovich, A. P.; Yuryshev, N. N. Physics and Engineering of Singlet Delta Oxygen Production in LowTemperature Plasma. J. Phys. D: Appl. Phys. 2007, 40, R25−R61. (16) Azyazov, V. N.; Torbin, A. P.; Pershin, A. A.; Mikheyev, P. A.; Heaven, M. C. Kinetics of Oxygen Species in an Electrically Driven Singlet Oxygen Generator. Chem. Phys. 2015, 463, 65−69. (17) Cvetanovic, R. J. Excited State Chemistry in the Stratosphere. Can. J. Chem. 1974, 52, 1452−1464. (18) Bodesheim, M.; Schmidt, R. Chemical Reactivity of Sigma Singlet Oxygen O2(b1Σ+g ). J. Phys. Chem. A 1997, 101, 5672−5677. (19) Bultel, A.; Cheron, B. G.; Bourdon, A.; Motapon, O.; Schneider, I. F. Collisional-Radiative Model in Air for Earth Re-Entry Problems. Phys. Plasmas 2006, 13, 043502. (20) Loukhovitski, B. I.; Starik, A. M. Modeling of VibrationElectronic-Chemistry Coupling in the Atomic-Molecular Oxygen System. Chem. Phys. 2009, 360, 18−26. (21) Starik, A. M.; Titova, N. S.; Arsentiev, I. V. Comprehensive Analysis of The Effect of Atomic and Molecular Metastable State Excitation on Air Plasma Composition behind Strong Shock Waves. Plasma Sources Sci. Technol. 2010, 19, 015007. (22) Yuryshev, N. N. Chemically Pumped Oxygen-Iodine Laser. Quantum Electron. 1996, 26, 567. (23) Azyazov, V. N. Excited States in The Active Media of OxygenIodine Lasers. Quantum Electron. 2009, 39, 989. (24) Huber, K. P.; Herzberg, G. Molecular Spectra and Molecular Structure. Vol. 4. Constants of Diatomic Molecules; Van Nostrand Reinhold: New York, 1979. (25) Starik, A. M.; Titova, N. S. Kinetics of Detonation Initiation in The Supersonic Flow of The H2+O2(air) Mixture in O2 Molecule Excitation by Resonance Laser Radiation. Kinet. Catal. 2003, 44, 28− 39. (26) Sharipov, A. S.; Starik, A. M. Kinetic Mechanism of CO−H2 System Oxidation Promoted by Excited Singlet Oxygen Molecules. Combust. Flame 2012, 159, 16−29. (27) Chukalovsky, A.; Rakhimova, T.; Klopovsky, K.; Popov, N.; Mankelevich, Y.; Proshina, O. Specific Features of The Kinetics of H2O2-O2(a1Δg) Mixtures: I. Formation and Quenching of Electronically and Vibrationally Excited HO2(A′) Molecules in H2-O2-O2(a1Δg) Mixtures at a Temperature of 300 K. Plasma Phys. Rep. 2014, 40, 34− 51. (28) Gerasimov, G. Y.; Shatalov, O. P. Numerical Modeling of the Ignition of Hydrogen-Oxygen Mixtures Under Nonequilibrium Conditions. J. Eng. Phys. Thermophys. 2014, 87, 1063−1070. 9609

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A

Constants in Polyatomic Molecules. J. Phys. Chem. 1987, 91, 5573− 5577. (74) Cremer, D.; Wu, A.; Larsson, A.; Kraka, E. Some Thoughts about Bond Energies, Bond Lengths, and Force Constants. J. Mol. Model. 2000, 6, 396−412. (75) Sharipov, A. S.; Loukhovitski, B. I.; Starik, A. M. Quantum Chemical Study of Small BnCm Cluster Structures and Their Physical Properties. Eur. Phys. J. D 2015, 69, 211. (76) Wang, X.; Ding, S.; Wang, P.; Xie, J.; Zhong, W. Algebraic Approach to the Potential Energy Surface and Vibration Energy of the Transition Molecule−HO2. Chin. J. Phys. 2005, 43, 1051−1057. (77) Polyansky, O. L.; Ovsyannikov, R. I.; Kyuberis, A. A.; Lodi, L.; Tennyson, J.; Zobov, N. F. Calculation of Rotation-Vibration Energy Levels of the Water Molecule with Near-Experimental Accuracy Based on an ab Initio Potential Energy Surface. J. Phys. Chem. A 2013, 117, 9633−9643. (78) Landau, L. D.; Lifshitz, E. M. Quantum Mechanics: NonRelativistic Theory; Pergamon Press: New York, 1977; Vol. 3. (79) Harvey, J. N. Understanding The Kinetics of Spin-Forbidden Chemical Reactions. Phys. Chem. Chem. Phys. 2007, 9, 331−343. (80) Minaev, B. F.; Minaeva, V. A. MCSCF Response Calculations of The Excited States Properties of The O2 Molecule and a Part of Its Spectrum. Phys. Chem. Chem. Phys. 2001, 3, 720−729. (81) Karkach, S. P.; Osherov, V. I. Ab initio Analysis of The Transition States on The Lowest Triplet H2O2 Potential Surface. J. Chem. Phys. 1999, 110, 11918−11928. (82) Michael, J. V.; Sutherland, J. W.; Harding, L. B.; Wagner, A. F. Initiation in H2/O2: Rate Constants for H2+O2 → H+HO2 at High Temperature. Proc. Combust. Inst. 2000, 28, 1471−1478. (83) Thurman, G. M.; Steckler, R. Ab Initio Direct Dynamics Study of OOH + H → H2 + 3O2. J. Phys. Chem. A 2001, 105, 2192−2197. (84) Kalugina, Y.; Alpizar, O. D.; Stoecklin, T.; Lique, F. A New Ab Initio Potential Energy Surface for The Collisional Excitation of O2 by H2. Phys. Chem. Chem. Phys. 2012, 14, 16458−16466. (85) Chefdeville, S.; Kalugina, Y.; van de Meerakker, S. Y. T.; Naulin, C.; Lique, F.; Costes, M. Observation of Partial Wave Resonances in Low-Energy O2-H2 Inelastic Collisions. Science 2013, 341, 1094−1096. (86) Tsang, W.; Hampson, R. F. Chemical Kinetic Data Base for Combustion Chemistry. J. Phys. Chem. Ref. Data 1986, 15, 1087− 1280. (87) Paul, P.; Warnatz, J. A Re-Evaluation of The Means Used to Calculate Transport Properties of Reacting Flows. Symp. Combust., [Proc.] 1998, 27, 495−504. (88) Langhoff, S. R.; Jaffe, R. L. Theoretical Study of The Four Lowest Doublet Electronic States of The Hydroperoxyl Radical: Application to Photodissociation. J. Chem. Phys. 1979, 71, 1475−1486. (89) Kendrick, B.; Pack, R. T. Potential Energy Surfaces for The Low-lying 2A” and 2A’ States of HO2: Use of The Diatomics in Molecules Model to Fit Ab Initio Data. J. Chem. Phys. 1995, 102, 1994−2012. (90) Francisco, J. S.; Zhao, Y. A Re-Examination of The Heat of Formation of HO2 Using Ab Initio Molecular Orbital Theory. Mol. Phys. 1991, 72, 1207−1208. (91) Baulch, D. L.; Bowman, C. T.; Cobos, C. J.; Cox, R. A.; Just, T.; Kerr, J. A.; Pilling, M. J.; Stocker, D.; Troe, J.; Tsang, W.; et al. Evaluated Kinetic Data for Combustion Modeling. J. Phys. Chem. Ref. Data 2005, 34, 757−1397. (92) Ruscic, B.; Pinzon, R. E.; Morton, M. L.; Srinivasan, N. K.; Su, M. C.; Sutherland, J. W.; Michael, J. V. Active Thermochemical Tables: Accurate Enthalpy of Formation of Hydroperoxyl Radical HO2. J. Phys. Chem. A 2006, 110, 6592−6601. (93) Curtiss, L. A.; Raghavachari, K.; Redfern, P. C.; Rassolov, V.; Pople, J. A. Gaussian-3 (G3) Theory for Molecules Containing First and Second-Row Atoms. J. Chem. Phys. 1998, 109, 7764−7777. (94) Curtiss, L. A.; Redfern, P. C.; Raghavachari, K. Gaussian-4 Theory. J. Chem. Phys. 2007, 126, 084108. (95) Montgomery, J. A., Jr; Frisch, M. J.; Ochterski, J. W.; Petersson, G. A. A Complete Basis Set Model Chemistry. VI. Use of Density

(51) Schmidt, M. W.; Gordon, M. S. The Construction and Interpretation of MCSCF Wavefunctions. Annu. Rev. Phys. Chem. 1998, 49, 233−266. (52) Freidzon, A.; Tsybizova, A. CASSCF and Firefly: A Tutorial. 2017; http://photon190573.livejournal.com/ (May 2017). (53) Bertran, J. Some Fundamental Questions in Chemical Reactivity. Theor. Chem. Acc. 1998, 99, 143−150. (54) Granovsky, A. A. Extended Multi-Configuration QuasiDegenerate Perturbation Theory: The New Approach to Multi-state Multi-Reference Perturbation Theory. J. Chem. Phys. 2011, 134, 214113. (55) Bytautas, L.; Ruedenberg, K. Accurate Ab Initio Potential Energy Curve of O2. I. Nonrelativistic Full Configuration Interaction Valence Correlation by The Correlation Energy Extrapolation by Intrinsic Scaling Method. J. Chem. Phys. 2010, 132, 074109. (56) Li, J.; Dawes, R.; Guo, H. An Accurate Multi-Channel MultiReference Full-Dimensional Global Potential Energy Surface for The Lowest Triplet State of H2O2. Phys. Chem. Chem. Phys. 2016, 18, 29825−29835. (57) Kendall, R. A.; Dunning, T. H., Jr.; Harrison, R. J. Electron Affinities of the First-Row Atoms Revisited. Systematic Basis Sets and Wave Functions. J. Chem. Phys. 1992, 96, 6796−6806. (58) Granovsky, A. A. Firefly V. 8.0. http://classic.chem.msu.su/ gran/firefly/index.html (2017). (59) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S. T.; Gordon, M. S.; Jensen, J. H.; Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S.; et al. General Atomic and Molecular Electronic Structure System. J. Comput. Chem. 1993, 14, 1347−1363. (60) Jacox, M. E. Vibrational and Electronic Energy Levels of Polyatomic Transient Molecules. Supplement B. J. Phys. Chem. Ref. Data 2003, 32, 1−437. (61) Baldridge, K. K.; Gordon, M. S.; Steckler, R.; Truhlar, D. G. Ab Initio Reaction Paths and Direct Dynamics Calculations. J. Phys. Chem. 1989, 93, 5107−5119. (62) Petersson, G. A. In Quantum-Mechanical Prediction of Thermochemical Data; Cioslowski, J., Ed.; Kluwer Academic Publishers: New York, 2002; Chapter Complete Basis Set Models for Chemical Reactivity: from the Helium Atom to Enzyme Kinetics, pp 99−130. (63) Sharipov, A. S.; Starik, A. M. Analysis of The Reaction and Quenching Channels in a H + O2(a1Δg) System. Phys. Scr. 2013, 88, 058305. (64) Sumathi, R.; Green, W. H., Jr A Priori Rate Constants for Kinetic Modeling. Theor. Chem. Acc. 2002, 108, 187−213. (65) Sharipov, A. S.; Starik, A. M. Theoretical Study of the Reactions of Ethanol with Aluminum and Aluminum Oxide. J. Phys. Chem. A 2015, 119, 3897−3904. (66) Miller, W. H.; Handy, N. C.; Adams, J. E. Reaction Path Hamiltonian for Polyatomic Molecules. J. Chem. Phys. 1980, 72, 99− 112. (67) Canneaux, S.; Bohr, F.; Henon, E. KiSThelP: A Program to Predict Thermodynamic Properties and Rate Constants from Quantum Chemistry Results. J. Comput. Chem. 2014, 35, 82−93. (68) Euclides, H. O.; Barreto, P. R. P. APUAMA: a Software Tool for Reaction Rate Calculations. J. Mol. Model. 2017, 23, 176. (69) Fernandez-Ramos, A.; Miller, J. A.; Klippenstein, S. J.; Truhlar, D. G. Modeling the Kinetics of Bimolecular Reactions. Chem. Rev. 2006, 106, 4518−4584. (70) Jasper, A. W.; Miller, J. A. Theoretical Unimolecular Kinetics for CH4+M = CH3+H + M in Eight Baths, M = He, Ne, Ar, Kr, H2, N2, CO, and CH4. J. Phys. Chem. A 2011, 115, 6438−6455. (71) Loukhovitski, B. I.; Sharipov, A. S.; Starik, A. M. Theoretical Study of Physical and Thermodynamic Properties of AlnNm Clusters. Eur. Phys. J. D 2016, 70, 250. (72) Loukhovitski, B. I.; Sharipov, A. S.; Starik, A. M. Quantum Chemical Study of Small AlnBm Clusters: Structure and Physical Properties. Chem. Phys. 2017, 493, 61−76. (73) Zavitsas, A. A. Quantitative Relationship between Bond Dissociation Energles, Infrared Stretching Frequencies, and Force 9610

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611

Article

The Journal of Physical Chemistry A Functional Geometries and Frequencies. J. Chem. Phys. 1999, 110, 2822−2828. (96) Carr, R. W., Jr Predictions of The Rates of Hydrogen Abstraction by CH2(3B1) by The Bond-Energy Bond-Order Method. J. Phys. Chem. 1972, 76, 1581−1586. (97) Westenberg, A. A.; De Haas, N. Steady-State Intermediate Concentrations and Rate Constants. HO2 Results. J. Phys. Chem. 1972, 76, 1586−1593. (98) Mebel, A. M.; Hayashi, M.; Kislov, V. V.; Lin, S. H. Theoretical Study of Oxygen Isotope Exchange and Quenching in the O(1D) + CO2 Reaction. J. Phys. Chem. A 2004, 108, 7983−7994. (99) Lehman, J. H.; Lester, M. I.; Klos, J.; Alexander, M. H.; Dagdigian, P. J.; Herráez-Aguilar, D.; Aoiz, F. J.; Brouard, M.; Chadwick, H.; Perkins, T.; et al. Electronic Quenching of OH A2Σ+ Induced by Collisions with Kr Atoms. J. Phys. Chem. A 2013, 117, 13481−13490. (100) Fedorov, D. G.; Koseki, S.; Schmidt, M. W.; Gordon, M. S. Spin-Orbit Coupling in Molecules: Chemistry beyond The Adiabatic Approximation. Int. Rev. Phys. Chem. 2003, 22 (3), 551−592. (101) Koseki, S.; Schmidt, M. W.; Gordon, M. S. Effective Nuclear Charges for The First- through Third-Row Transition Metal Elements in Spin-Orbit Calculations. J. Phys. Chem. A 1998, 102, 10430−10435. (102) Fedorov, D. G. Theoretical Study of Spin-Orbit Coupling in Molecules. Retrospective Theses and Dissertations. Paper 12662, Iowa State University, Ames, IA, 1999. (103) Konnov, A. A. Remaining Uncertainties in The Kinetic Mechanism of Hydrogen Combustion. Combust. Flame 2008, 152, 507−528. (104) Bezgin, L. V.; Kopchenov, V. I.; Sharipov, A. S.; Titova, N. S.; Starik, A. M. Evaluation of Prediction Ability of Detailed Reaction Mechanisms in the Combustion Performance in Hydrogen/Air Supersonic Flows. Combust. Sci. Technol. 2013, 185, 62−94. (105) Koike, T. Shock Tube Studies of The H2-O2 Reaction by Atomic Resonance Absorption Spectroscopy. Bull. Chem. Soc. Jpn. 1989, 62, 2480−2483. (106) Burke, M. P.; Chaos, M.; Ju, Y.; Dryer, F. L.; Klippenstein, S. J. Comprehensive H2/O2 Kinetic Model for High-Pressure Combustion. Int. J. Chem. Kinet. 2012, 44, 444−474. (107) Fernandez-Ramos, A.; Ellingson, B. A.; Meana-Paneda, R.; Marques, J. M. C.; Truhlar, D. G. Symmetry Numbers and Chemical Reaction Rates. Theor. Chem. Acc. 2007, 118, 813−826. (108) Hack, W.; Kurzke, H. The Production of H2(S) Atoms on The Energy-Transfer Reaction of O2(1Δg) with HO2(A″). Chem. Phys. Lett. 1984, 104, 93−96. (109) Christoffel, K. M.; Bowman, J. M. RRKM Theory beyond the Separable Harmonic Approximation: The HCO2 -> H + CO2 Unimolecular Decomposition. J. Phys. Chem. A 1999, 103, 3020−3030. (110) Carr, R. W. In Comprehensive Chemical Kinetics; Green, N. J. B., Ed.; Elsevier, 2007; Vol. 42; Chapter Elements of Chemical Kinetics, pp 43−99. (111) Sharipov, A. S.; Starik, A. M. Theoretical Study of the Reaction of Ethane with Oxygen Molecules in the Ground Triplet and Singlet Delta States. J. Phys. Chem. A 2012, 116, 8444−8454. (112) Kaledin, A. L.; Heaven, M. C.; Morokuma, K. Theoretical Prediction of The Rate Constant for I+O2(a1Δg) Electronic Energy Transfer: A Surface-Hopping Trajectory Study. J. Chem. Phys. 2001, 114, 215−224. (113) Takayanagi, T. Quantum Scattering Calculations of the O(1D) + N2(X1Σg+) -> O(3P) + N2(X1Σg+) Spin-Forbidden Electronic Quenching Collision. J. Phys. Chem. A 2002, 106, 4914−4921. (114) Topaler, M. S.; Allison, T. C.; Schwenke, D. W.; Truhlar, D. G. Test of Trajectory Surface Hopping Against Accurate Quantum Dynamics for an Electronically Nonadiabatic Chemical Reaction. J. Phys. Chem. A 1998, 102, 1666−1673. (115) Tully, J. C. Perspective: Nonadiabatic Dynamics Theory. J. Chem. Phys. 2012, 137, 22A301.

9611

DOI: 10.1021/acs.jpca.7b09964 J. Phys. Chem. A 2017, 121, 9599−9611