Quantum Chemical Study of the Initial Step of Ozone Addition to the

Article pubs.acs.org/JPCA

Quantum Chemical Study of the Initial Step of Ozone Addition to the Double Bond of Ethylene Oleg B. Gadzhiev,*,†,⊥ Stanislav K. Ignatov,† Boris E. Krisyuk,*,‡,⊥ Alexey V. Maiorov,§ Shruba Gangopadhyay,∥,⊥ and Artem ̈ E. Masunov∥ †

Department of Chemistry, N.I. Lobachevsky State University of Nizhny Novgorod, National Research University, 23 Gagarin Avenue, Nizhny Novgorod 603950, Russia ‡ Institute of Problems of Chemical Physics, Russian Academy of Sciences, 1 Academician Semenov Av., Chernogolovka, Moscow Region 142432, Russia § N. M. Emanuel Institute of Biochemical Physics, Russian Academy of Sciences, 4 Kosygin St., Moscow 119991, Russia ∥ NanoScience Technology Center, Department of Chemistry and Department of Physics, University of Central Florida, 12424 Research Parkway, Ste 400, Orlando, Florida 32826, United States S Supporting Information *

ABSTRACT: The mechanisms of the initial step in chemical reaction between ozone and ethylene were studied by multireference perturbation theory methods (MRMP2, CASPT2, NEVPT2, and CIPT2) and density functional theory (OPW91, OPBE, and OTPSS functionals). Two possible reaction channels were considered: concerted addition through the symmetric transition state (Criegee mechanism) and stepwise addition by the biradical mechanism (DeMore mechanism). Predicted structures of intermediates and transition states, the energies of elementary steps, and activation barriers were reported. For the rate-determining steps of both mechanisms, the full geometry optimization of stationary points was performed at the CASPT2/cc-pVDZ theory level, and the potential energy surface profiles were constructed at the MRMP2/cc-pVTZ, NEVPT2/cc-pVDZ, and CIPT2/cc-pVDZ theory levels. The rate constants and their ratio for reaction channels calculated for both mechanisms demonstrate that the Criegee mechanism is predominant for this reaction. These results are also in agreement with the experimental data and previous computational results. The structure of DeMore prereactive complex is reported here for the first time at the CCSD(T)/cc-pVTZ and CASPT2/ccpVDZ levels. Relative stability of the complexes and activation energies were refined by single-point energy calculations at the CCSD(T)-F12/VTZ-F12 level. The IR shifts of ozone bands due to formation of complexes are presented and discussed.

1. INTRODUCTION

tion by the ecologically safe (green) method of ozonolysis under pressure.5 The mechanism of ozonolysis (C2H4 + O3 → products) was intensively studied by experimental (the last kinetic experiments,6 regularly generalized data sets of IUPAC,7 JPL,8 matrix isolation,9 and microwave spectroscopy10) and theoretical methods, from very early ones11 to the most recent and sophisticated.12 The critical remarks on the values collected in the IUPAC and JPL databases have been considered in ref 13. Progress in the area of theoretical chemistry and computational techniques made it possible to use the calculation methods that allow one to quantitatively estimate the activation energy and investigate reaction dynamics at least for the gasphase conditions14 and therefore extrapolate some conclusions to reaction dynamics at the interface and in droplets.

The interaction of ozone with the double bond of unsaturated compounds is important for atmospheric chemistry and formation of photochemical smog. Both reactants are important components of troposphere: ozone formed photochemically and unsaturated compounds (in particular, alkenes) are plant vegetation products possessing significant biochemical functions in plants. The reaction is not only a sink for O3 and unsaturated compounds but also a source of oxygen-containing species (aldehydes, ketones, acids, and radicals HO and HO2) and can contribute to the formation of secondary aerosols.1 The generation of hot molecules and radicals resulting in luminescence in an ozone/unsaturated hydrocarbon system is an interesting feature of the ozonolysis reaction,2 used for some practical applications, including ozonometry.3 Oxidation O3 + C2H4 → products, under the conditions close to normal in N2/ O2 mixtures, is a model of chemical degradation of rubbers4 and biochemical activity of ozone. It is also important for conservation of water resources and wastewater decontamina© 2012 American Chemical Society

Received: August 4, 2012 Revised: September 17, 2012 Published: September 19, 2012 10420

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

molozonide, some specific features of ozonolysis via the Criegee mechanism26 indicate that the Criegee and DeMore mechanisms can compete efficiently. This competition includes the switching between the mechanisms, affected by the various substituent effects. It is proposed that the Criegee mechanism predominates in the case of the reaction of ozone with ethylene26,27 and other investigated alkenes, whereas the DeMore mechanism is favorable for tetrafluoroethylene25 and hexafluoropropylene.28 However, reliable data on the energy and activation barriers of the individual elementary steps and on the relative contribution of the mechanisms for alkenes and alkynes with different substituents to the net reaction are lacking so far. This is due to experimental difficulties in the studies of fast gas-phase highly exothermic reactions, including the generation of electron-translational−vibrational−rotational hot molecules and radicals. However, this specific feature of the thermal oxidation of ethylene by ozone makes it possible to choose the predominant mechanism among the considered ones only by the activation energy of one elementary reaction, i.e., by activation energy of the first step. In this case, it is necessary to employ highly accurate computational procedures for the ground state PES, which can treat efficiently relaxation of the structures. Quantum chemical calculations by the Hartree−Fock and MP2 methods,11a,b as well as semiempirical and empirical methods, prevailed in early works. The spin-restricted KS-DFT methods resulted in the detection of molecular TS1, whereas the application of spin-unrestricted KS-DFT gives both TS2 and TS1. Both transition states were described at one level by the B3LYP and PBE0 DFT methods, ab initio QCISD and CCSD methods, and by the MRMP2 method with geometry optimization performed with CASSCF and UMP2 methods.27 These theory levels do not account for static and dynamics electron correlation in a balanced way. Therefore, they cannot serve as a reliable tool to study the mechanisms and make it possible to perform only a qualitative estimation of the energy of elementary reactions. A similar situation takes place for the ozonolysis of C2H4 and its substituted derivatives. The CASSCF method was used29 for geometry optimization and the multireference CASPT2 method for energy refinement of the found stationary points. The main drawback of this study29 was substantial inaccuracy in the optimized geometric parameters with the CASSCF method in the model restricted active space, which resulted in insufficient treatment of the dynamic electron correlation. Later, the BHandHLYP method was used with a broken symmetry approach to optimize the geometry, and the CCSD(T)/6-311+(2d,p) method was used to refine the energy.30 The estimated Ea, ΔH‡, ΔG‡, canonical rate constant k (T = 298.15 K), and parameters EA and A in the Arrhenius equation were obtained in the range T = 235−362 K for the Criegee and DeMore (exo- and endo-TSs) mechanisms. The Criegee mechanism was concluded to be the major reaction channel. The QCISD(T)/6-311++G(d,p)//MP2/6-311++G(d,p) composite method was applied30 to study the reactivity of ethylene and its derivatives (CH2CHF and CH2CHCl). This study showed that, for ethylene, the stabilization energy of the prereactive complex and activation energy Ea relatively to the initial (completely separated) reactants are 0.25 and 9.21 kJ/ mol, respectively. The estimated Ea has recently been presented12a for the C2H4 ozonolysis following the Criegee mechanism on the basis

According to ref 15, the addition of ozone to ethylene proceeds via the 1,3-cycloaddition mechanism through the symmetric transition state TS1 to form five-membered cyclic molozonide (primary ozonide, POZ) in the first elementary step of Criegee mechanism (Scheme 1). Scheme 1. First Step of Criegee Mechanisma

a

Reactants, ethylene and ozone; saddle point TS1 of synchronous [2 + 3], cycloaddition; and product, primary ozonide (POZ).

Later, this reaction scheme was confirmed by numerous direct and indirect experimental observations.16 In the late 1980s, primary molozonide was experimentally identified. Its structure and the structure of the reactive complex preceding its formation17 were studied by the microwave spectroscopy,18 combined with GVB-RCI19 and MP4(SDTQ) theoretical methods in 6-31G(d,p), 6-311G(d,p), and 6-31G(2d,p) basis sets. An alternative reaction mechanism, where ozone reacts with the double bond similarly to the peroxyl radical to form the biradical transition state TS2, was proposed by DeMore20 and usually named the stepwise mechanism or DeMore mechanism (Scheme 2). Originally, it was suggested to be favored in the case of the acetylene reaction.20 Scheme 2. First Step of DeMore Mechanisma

a

Reactants, ethylene and ozone; saddle points, exo-TS2 and endoTS2; and products of the two distinctive channels, ethylene epoxide and primary ozonide (POZ).

There are other hypotheses for the reaction mechanism: electron transfer with the formation of an ion−radical complex in the reactions of electron-donor substituted alkenes21 and H atom abstraction from ethylene.22 There are no data confirming the one-electron transfer for the reaction C2H4 + O3 → products in the gas phase.23 Several reaction channels after POZ have been considered for the gas-phase reaction.24 While DeMore mechanism explains25 the formation of such products as epoxides and aldehydes upon the decomposition of 10421

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

of the focal point analysis method31 with extrapolation to the infinite (complete) basis set limit (CBS/focal point) as 3.43 ± 0.20 kcal/mol (ΔrH‡(0 K) = 5.34 ± 0.20 kcal/mol) and the relative energy of the prereactive complex equal to −1.84 ± 0.20 kcal/mol. The applicability and metrological characteristics of various quantum chemical methods and approximations for the calculation of the activation energy were studied in detail in ref 12b. The calibration values for the relative energy of the prereactive complex and activation energy were found to be −1.94 and 3.50 kcal/mol (CBS CCSDT+(2)Q) and −2.03 and 3.18 kcal/mol (CBS CCSD(T)+Q). Average-weighed values of 0.25 and 0.72 kcal/mol, respectively, were proposed12b as the best estimates. The concerted and stepwise ozonolysis mechanisms were studied12c with the modern DFT variants (UBHandHLYP, M06, LC-ωPBE, and other functionals) for geometry optimizations and with CCSD(T), BD(T), CAS(4,4)MRMP2, and the high-correlated multireference MK-CCSD method to refine the energies of the found stationary points. The conclusion12c was that the concerted mechanism is more favorable for C2H4, whereas no unambiguous conclusion about the preference of this or another channel could be made for C2H3CN. Thus, even the use of the modern multireference coupled cluster method does not allow one to reliably determine the energies of the TSs and intermediates, when no geometry optimization is performed. Incomplete account for nondynamic electron correlation by the DFT functionals can even result in qualitatively incorrect topology of the PES. Indeed, search for stationary points using various DFT functionals and verification by the composite CCSD(T)// DFT approximation can be unreliable, as shown in ref 32 (cf., ref 33), where not only quantitative disagreement between DFT and CCSD(T) results was observed but also qualitatively different PES topology (reaction mechanism) was determined. In our opinion, the stationary point search and f ull geometry optimization with the CCSD(T) and CASSCF methods is more reliable verification method. Similar methodological arguments were proposed in ref 34. Therefore, it is necessary to apply computational schemes, which combine a balanced account for dynamic and nondynamic correlation with a stationary point search and construction of the PES near transition states. The objectives of this work include the structures of the reactants, prereactive complexes, and transition states (TSs) of the initial reaction step for the Criegee and DeMore mechanisms. The goal is to estimate the energies of elementary steps and energy barriers by methods of the multireference perturbation theory with full geometry optimization. The complete calculation of the reaction pathways corresponding to the initial stages of both ozonolysis mechanisms is also attempted. We also compare the single- and multireference quantum chemical methods and different calculation procedures of the activation energy for the elementary reactions.

NEVPT240 methods were carried out in the UCF by means of the MOLPRO 2010.1 program package41 using the computational resources provided by DOE NERSC. Further calculations were performed with CASSCF method in combination with the cc-pVDZ, aug-cc-pVDZ, and cc-pVTZ basis sets applying the geDIIS optimization method as coded in Gaussian09 suite of programs;42 DIIS and quadratic steepest descent following (QSD) methods were as coded in MOLPRO. 2.1.2. Dual-Level Approach and the PES Profiles Construction. The particular choice of the modern variants of multireference perturbation theory (NEVPT2 and CIPT2) in combination with the CASPT2 formalism is capable in balanced account of both static and dynamic electron correlation effects. It is worthwhile to briefly mention specific features of these methods. The NEVPT2 method40 is based on the Dyall modified Hamiltonian,43 which is two-particle for active orbitals,40 unlike the other variants of the multireference perturbation theory that are one-particle based. For this reason, it has an advantage over other multireference perturbation theory formalisms of the both types (diagonalization− perturbation). The CIPT2 method39 is a hybrid (composite) approach in which excitations from the active space are variationally treated by the MRCI method, and others are treated by the perturbation theory. This method is stable to the appearance of intruder states but is not strictly size-extensive, which can be compensated by the Davidson correction. In addition, the conventional state-specific MRMP2 method44 implemented in GAMESS-US was employed to verify the computational results obtained by NEVPT2 and CIPT2 methods. For the CASSCF wave function calculations, the active space for ethylene corresponds to two electrons in two orbitals (2,2), i.e., bonding and antibonding π-orbitals. The active space for ozone was chosen to be CAS(8,7) and CAS(12,9), following ref 34. This translates to the CAS(10,9) and CAS(14,11) active spaces for the prereactive complexes and transition states. After the TS geometry was determined at the CAS(10,9)/631+G(d,p) and CAS(14,11)/cc-pVTZ theory levels, the reaction coordinate was scanned toward the reactants until the complete convergence of the intrinsic reaction coordinate (IRC) calculations with energy gradients below 10−4 a.u. (default “OPTTOL” value of Firefly package). We also employ a dual−level approximation for the quantitative estimation of the activation energy (Ea) without full optimization of the geometric parameters in the vicinity of two stationary points. Following ref 45, we use the energy calculation with the higher level correlated method (for example, by the multireference perturbation theory) for each point of the IRC curve, built with the method that accounts electron correlation to a lesser extent (for instance, CASSCF). Taking into account the shift of the TS, it is possible to partially compensate the incomplete relaxation of the geometric parameters of the structure corresponding to this figurative point on the PES (see ref 45 for details). The recent application of the dual-level methodology is the study of vinyl hydroperoxide decomposition by a combination of MR-CISD and CASSCF methods.46 This approximation is similar to the IRCMax method.47 Several dual-level approaches (construction of the IRC curve followed by the energy refinement using another method) was reviewed in ref 48, including interpolation corrections. In the present study, an appearance of a local minimum on the PES profile built by renewal of the IRC curve was taken into account consequently using the IRC

2. COMPUTATIONAL DETAILS 2.1. Quantum Chemical Calculations. 2.1.1. General Remarks. The calculations were performed at the Computational Center of the Institute of Problems of Chemical Physics (Russian Academy of Sciences, Chernogolovka) using the Gaussian03,35 GAMESS-US,36 and Firefly37 programs. The Gaussian0335 and CFOUR38 program packages were also used at the Nizhny Novgorod State University. The search for stationary points at the CASPT2/cc-pVDZ theory level and the energy calculations by the CIPT239 and 10422

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

lower than the energy of the fully optimized local minimum. This is likely to be an artifact of the CASPT2 method, because of the external (appeared within the perturbation theory) intruder state. The application of vertical shifts (0.1 and 0.3 a.u.) did not correct the problem. In this case, this method of construction of the quasi-isolated system is inapplicable. The activation energy relative to the initial reactants (isolated O3 and C2H4) was determined by activation energy Ea of the reaction O3···C2H4 → TS and the stabilization energy of the prereactive complex O3···C2H4. 2.1.4. Energy Refinement and Further Calculations with Coupled Cluster and DFT Methods. In order to verify the results of multireference methods and the determination of the isolated system, we performed two sets of calculations: (1) single-point calculations at CCSD(T)-F12/VTZ-F12 theory level for quasi-isolated O3 and C2H4, TSs of Criegee and DeMore mechanisms, and corresponding O3···C2H4 complexes optimized at the RS2C/cc-pVDZ level; and (2) CCSD(T)F12/VTZ-F12 calculations on O3 and C2H4 monomers, Criegee TS, and the corresponding O3···C2H4 complex. The methods of explicit treatment of electron correlation, including coupled cluster methods, e.g., CCSD(T)-F12 (CCSD(T)-F12a and CCSD(T)-F12b),52 are characterized by faster convergence to basis set limit and stability where influence of nondynamic correlation is important. While reference values were obtained at the most sophisticated theory levels, different methods of KS-DFT (OPBE, OPW91, and OTPSS) were employed to demonstrate ability of the available exchange-correlation functionals to produce qualitatively correct description the molecular systems with significant influence of nondynamic electron correlation and to benchmark a set of computationally affordable methods that can be applied to the systems of experimental interest. The motivation of the specific choice of DFT methods will be discussed in detail. The ChemCraft,53 MOLTRAN,54 and GaussView55 programs were used for the visualization, data processing, and preparation of the materials. 2.2. Calculation of the Rate Constants. Conventional statistical thermodynamics, molecular partition functions based on the harmonic oscillator and rigid rotor approximations, were used to calculate thermodynamic properties (the activation enthalpies and entropies) as a function of temperature. This procedure was performed automatically by means of the MOLTRAN program.54 The kinetic coefficients (canonical rate constants k(T)) were calculated using the conventional transition state theory (TST) in the temperature range from 200 to 400 K for the reaction O3 + C2H4 → TS → products via the Criegee and DeMore mechanisms if both saddle points were localized (the energy and vibrational frequencies were calculated at the CASPT2/cc-pVDZ level) and via the Criegee mechanism only using the data calculated in ref 12a at the CCSD(T)/cc-pVTZ level were used. In ref 12a, no vibrational frequencies were calculated at the CCSD(T)/cc-pVTZ level, but it was shown for the TS [C2H2···O3]‡ that the extension of the basis set from cc-pVDZ to cc-pVTZ did not change the result substantially. The estimated inaccuracy for O3 + C2H4 → TS does not exceed 0.42 kJ/mol (0.1 kcal/mol). We calculated the vibrational frequencies for the fully optimized structures of O3, C2H4, and Criegee TS reported in ref 12a (cf., Supporting Information) at the CCSD(T)/cc-pVTZ level. The difference for ΔZPVE of the elementary process O3 + C2H4 → TS at the CCSD(T,FC)/

curve following method, which will allow one to localize the prereactive complex, estimate its stabilization energy, and refine the activation energy (Ea) of the reaction by the multireference method when the full geometry optimization cannot be performed. The calculations with NEVPT2 and CIPT2 methods were carried out in two different approximations: (1) NEVPT2/ccpVDZ using strong (SC) and partial (PC) contraction formalisms, respectively; (2) CIPT2/cc-pVDZ and CIPT2+Dav/cc-pVDZ (with Davidson’s correction to attain the size consistency). The change of the basis set from 6-31+G(d,p) to cc-pVDZ is due to the necessity to satisfy the electron−electron cusp condition for a more rigorous account for the dynamic correlation energy by the NEVPT2 and CIPT2 methods. Note that it would be more valid to employ the basis set containing polarization functions with the highest angular moments (from f function) for at least O and C atoms for the particular system. However, this would make the CIPT2 and NEVPT2 calculations very time and resource consuming. It was feasible to employ the largest (nearly size consistent) active space CAS(14,11) in conjunction with the very flexible cc-pVTZ basis set for the MRMP2 energy calculations along the IRC curve calculated at the CAS(14,11)/cc-pVTZ level. These calculations were performed by means of the GAMESSUS program package. 2.1.3. PES Topology Determination by Full Geometry Optimization. The geometric parameters of the structures corresponding to minima and saddle points on the PES were fully optimized at the CASPT2/cc-pVDZ level. We chose the CAS(10,9) reference wave function for the CASPT2 method49 (defined as Rayleigh−Schrödinger perturbation theory), the smallest of the active spaces describing both TSs and elementary reaction pathways. This allowed us to reduce computational time, retaining the correct description of the C2H4 + O3 system on the PES fragment corresponding to the reaction path C2H4 + O3 → [C2H4 + O3] → TS. The contracted variant of the multireference perturbation theory CASPT2, i.e., RS2C formalism,49c was used because of the convergence problems in the PT2 procedure upon the application of the noncontracted CASPT2 method, i.e., RS2 formalism.49a,b To reduce the calculation time, the convergence criterion of the stationary point search was specified by the gradient value equal to 10−4 a.u. (0.26 kJ/mol) being an appropriate value for weakly bound molecular complexes. The QSD method50 did not allow us to find any of the TSs (DeMore and Criegee TSs) since oscillation and even divergence were observed. The DIIS method51 made it possible to localize stationary points. These saddle points corresponded to those found by the TS search at the CAS(10,9)/631+G(d,p) and CAS(14,11)/cc-pVDZ levels. The calculation of vibrational frequencies in the harmonic approximation with the numerical calculation of the Hessian at the RS2C/cc-pVDZ level made it possible to characterize the found stationary points as two local minima (prereactive complexes) and two saddle points, identified as Criegee and DeMore TSs by visualization of nuclear motions corresponding to the single negative eigenvalue of the Hessian. The stabilization energy of the prereactive complexes was determined by the constrained optimization (fixed C−OO2 distance) of the quasi-isolated system (molecules O3 and C2H4 at 4.5 Å distance) at the RS2C/cc-pVDZ level. It was found by the relaxed scan along the internuclear C−O(OO′) distance longer than 4.5 Å that the energy of the disconnected system is 10423

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

cc-pVDZ and CCSD(T,FC)/cc-pVTZ calculation levels is 0.72 kJ/mol, which agrees with the estimation given in ref 12a. In the cases with no full geometry optimization, the Gibbs activation energy, ΔG‡, was calculated from the activation energy determined by the IRC curve following method at the NEVPT2/cc-pVDZ and CIPT2/cc-pVDZ levels and thermodynamic corrections obtained at the RS2C/cc-pVDZ level if the CAS(10,9) active space was employed or with activation energy calculated by the dual-level procedure with MRMP2/cc-pVTZ energy calculation and vibrational frequencies determined at the CAS(14,11)/cc-pVTZ theory level. The A and EA parameters in the Arrhenius eq 1

⎛ E ⎞ k(T ) = A exp⎜ − A ⎟ ⎝ RT ⎠

(1) Figure 1. Structural parameters R, α1, α2, and α3 of intermolecular complex O3···C2H4, determined by microwave spectroscopy experiments.18 R, distance between centers of CC and two terminal O atoms; α1, angle between segement R and bisector of angle ∠OOO; α2 and α3, angles between segment R and local C2 axes of HCCH fragments.

were estimated at T = 298 K on the basis of TST by eqs 2 and 3. 0 Δv ‡ ⎛ ‡⎞ kBT ⎛ pi ⎞ ⎜ ⎟ exp⎜ ΔS ⎟ A= ⎜ ⎟ h ⎝ RT ⎠ ⎝ R ⎠

(2)

EA = ΔH ‡ ,0 + 2RT

(3)

Surprisingly, the DeMore exo-TS was localized at the CCSD/ aug-cc-pVDZ and QCISD/aug-cc-pVDZ theory levels.27 We explain this surprising result by the fact that the less sophisticated theory level provides a qualitative agreement with more rigorous multireference methods by the known artifact of perturbative inclusion of triple excitations in the coupled cluster Ansatz, which produces a wrong topology of potential energy curve. The structure of the O 3 ···C 2 H 4 molecular complex determined in ref 18 by microwave spectroscopy is characterized by parameters α1, α2, α3, and R that determine the orientation of the O3 and C2H4 molecules, whose internal coordinates were accepted to be equal to their equilibrium values in the gas phase (Figure 1). For these parameters, the values are given: α1 = 105.8° and 107.6°, α2 = 105.3° and 110.7°, α3 = 74.77° and 69.3°, and R = 3.061 and 3.216 Å as calculated at the CCSD(T)/cc-pVTZ level and determined experimentally, respectively. The prereactive complex is characterized by the point symmetry group Cs and the structure in which the planes of the O3 and C2H4 molecules are not parallel. The calculated rotation constants Ae, Be, and Ce and the centrifugal distortion constants for the ground vibrational state based on the Watson Hamiltonian in the quaternary approximation of centrifugal distortion58 (representation Ir, Areduction) for the intermolecular complex O3···C2H4 and complexes of O3 with ethylene isotopomers (1,1-d2-dideuterioethylene, cis-1,2-d2-dideuterioethylene, and trans-1,2-d2-dideuterioethelene are presented in Tables 1 and 2 along with the experimental data (rotation constants A, B, and C and centrifugal distortion constants). The centrifugal distortion constants for the S-reduced Hamiltonian are given in Supporting Information. The activation energies of elementary reactions (Figures 2 and 3) and geometric parameters (Figures 4 and 5) of stationary points were obtained with the CASSCF method in the active spaces CAS(10,9) and CAS(14,11) in order to take into account the nondynamic electron correlation, which is important for both the Criegee and DeMore TSs search and construction of the MEPs. For the first one, the weight of the first configuration in the expansion of the CI vector is 0.80,

For the density functionals, a similar procedure was carried out with optimizations and frequency calculations at the same theory level.

3. RESULTS AND DISCUSSION 3.1. Quantum Chemical Methods Validation: Calculated vs Experimental Data. The geometric parameters obtained by the full optimization of the reactants, prereactive complexes, and TS for the Criegee mechanism at the RS2C/ccpVDZ level calculated in the presented study and CCSD(T,FC)/cc-pVTZ level used in ref 12a are reported in the Supporting Information and will be discussed in the corresponding sections. These the highest multi- and singlereference theory levels employed thus far for the study of the PES topology. The coupled cluster results are consistent with the experimental data.56 The values are given in the following order: the calculated at the RS2C/cc-pVDZ, CCSD(T,FC)/ccpVTZ levels and experimental data.56 For O3 molecule, the O− O bond lengths are 1.288, 1.272, and 1.276 Å, and the bond angles are 115.9°, 117.0°, and 116.97°; for C2H4 molecule, r(C−H) = 1.094, 1.083, and 1.081 Å, r(CC) = 1.351, 1.337, and 1.334 Å, ∠(HCH) = 117.2, 117.1, and 121.3; for the prereactive complex C2H4···O3 (Figure 1), r(C−C) = 1.350, 1.338, and (1.339) Å, r(C···O) = 3.051, 3.089, and 3.243 Å, r(O−O) = 1.287, 1.275, and (1.276) Å; for Criegee TS, r(C− C) = 1.386, 1.373 Å, r(C···O) = 2.186, 2.187 Å, r(O−O) = 1.314, 1.304 Å. The presented experimental data belong to the equilibrium structure for O3, ethylene molecules, and vibrationally averaged C2H4···O3 complex (Figure 1). For the discussion of inaccuracies at the CCSD(T)/cc-pVTZ level for calculation of C2H4 geometrical parameters and extrapolation of experimental data to determine re, see ref 57. The CCSD(T)/cc-pVTZ level is the minimal singlereference method that provides correct symmetry of vibrational modes and vibrational frequencies in ozone molecule. It was impossible to localize the DeMore TS at the CCSD(T)/ccpVTZ theory level because the effect of nondynamic correlation is beyond the convergence radius of the CCSD(T) method, as indicated by high (0.05) value of the T1 diagnostic. 10424

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

Table 1. Rotation Constants A, B, and C and Centrifugal Distortion Constants for O3···CH2CH2 and O3···trans-CHDCHD Determined Experimentally18 and Calculated at the CCSD(T,FC)/cc-pVTZ Theory Level O3···CH2CH2

O3···trans-CHDCHD

spectral constant (MHz)

A1

A2

theory

A1

A2

theory

A B C ΔJ ΔJK ΔK δJ δK

8246.841(2) 2518.972(4) 2044.248(5) 0.01596(3) 0.2333(2) −0.0524(6) 0.00316(2) 0.168 (2)

8241.897 (4) 2518.941(9) 2044.287(1) 0.01577(6) 0.2350(4) −0.405(1) 0.00328(3) 0.147(5)

8132.236 2786.515 2211.926 0.01310 0.3238 −0.33268 0.00275 0.106

7673.001(51) 2409.216(54) 1958.598(37) 0.0151(12) 0.2116 (20) −0.127(14) 0.002 98(17) 0.152(19)

7671.214(8) 2409.214(1) 1958.581(1) 0.014768(7) 0.21120(6) −0.2813(2) 0.002955(5) 0.1554(6)

7558.486 2651.046 2111.354 0.01187 0.28326 −0.29097 0.00243 0.0872

Table 2. Rotation Constants A, B, and C and Centrifugal Distortion Constants for O3···CD2CH2 and O3···cis-CHDCHD determined experimentally18 and calculated at the CCSD(T,FC)/cc−pVTZ theory level O3···CD2CH2

O3···cis-CHDCHD

spectral constant (MHz)

low- frequency state

high-frequency state

theory

exptl

theory

A B C ΔJ ΔJK ΔK δJ δK

7697.423(2) 2408.323(3) 1961.439(3) 0.015 22(3) 0.2056(1) −0.2893(5) 0.00313(2) 0.152(1)

7699.243(2) 2408.332(3) 1961.451(4) 0.01555(3) 0.2062(2) −0.1343(5) 0.003 26(1) 0.151(2)

7588.1511 2647.0670 2112.9727 0.01196 0.28300 −0.28943 0.00242 0.08424

7689.330 (2) 2433.481 (3) 1974.958 (5) 0.01373(2) 0.2465(2) −0.2380(5) 0.00281(2) 0.154 (2)

7573.9543 2669.9510 2124.0192 0.01058 0.31299 −0.31942 0.00218 0.09248

Figure 2. PES profiles corresponding to first step of Criegee mechanism calculated by IRC method at the CAS(10,9)/6-31+G(d,p) and by IRC following method with energy calculations at the NEVPT2/cc−pVDZ, CIPT2/cc−pVDZ, and CIPT2/cc−pVDZ with Davidson’s correction levels: 1, reagents valley; TS, saddle point; 2, product valley.

of the endo-TS for the DeMore mechanism has been reported earlier.12c However, the stationary points were optimized for various DFT functionals and for BS-UDFT approximation, and the multireference method was applied to validate the values of activation energy. This methodology is not always applicable32,33 so we additionally performed the search of the first TS for the two-step mechanism with analysis of various algorithms for saddle-point geometry optimization and the strategy for the numerical calculation of the Hessian. All the calculations were carried out with CASSCF method in various active spaces. The endo-TS structure denoted as NTS1 in ref 12c was chosen as the initial approximation.

while for the second case, it is 0.63 for the largest applied active spaces, i.e., active space CAS(14,11). Obviously, the quantitative comparison with experimental data requires the application of a method that covers equivalently the dynamic and nondynamic electron correlations. The geometric parameters of the structures obtained by the method of the IRC following with the NEVPT2 and CIPT2 methods in combination with the cc-pVDZ basis set are similar with those optimized at the RS2C/cc-pVDZ level (Figures 4−8). Both transition states and the corresponding prereactive complexes were found at this theory level. 3.2. PES Topology Verification: How Many TSs for the Initial Stage of the DeMore Mechanism? The localization 10425

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

Figure 3. PES profiles corresponding to first step of DeMore mechanism calculated by IRC method at the CAS(10,9)/6-31+G(d,p) and by IRC following method with energy calculations at the NEVPT2/cc-pVDZ, CIPT2/cc-pVDZ, and CIPT2/cc-pVDZ with Davidson’s correction levels: 1, reagents valley; TS, saddle point; 2, product valley.

In this case, the CAS(14,11) active space has the size consistency. We applied several algorithms for stationary point search: QSD and GDIIS in MOLPRO and geDIIS in Gaussian09 with standard (for these program packages) convergence criteria without any geometry constraints. In addition, the convergence criteria standard for the Gaussian09 program package and the Baker convergence criteria59 were used for MOLPRO. For the Gaussian09 program package, the search of the DeMore endo-TS converged to the Criegee TS characterized by the Cs symmetry point group or the slightly distorted Cs (the difference between the C···O internuclear distances does not exceed 0.005 Å). For the MOLPRO program package, the search resulted in a structure close to the Criegee TS, namely, [3 + 2] cycloaddition when the numerical calculation of the Hessian was performed each for five cycles. However, in this case, the optimization converged to the structure with nonequivalent internuclear C···O distances: for the CAS(10,9) active space, the difference does not exceed 0.05 Å for the DIIS method and 0.056 Å for the QSD method; for the CAS(14,11) active space, the difference is 0.098 Å for the aug-cc-pVDZ basis set and the difference does not exceed 0.027 Å for the cc-pVDZ and cc-pVTZ basis sets. The calculation of the Hessian, each for two optimization cycles, gave no substantial change in the geometric parameters in comparison with that reported above. Various convergence criteria did not change significantly the geometric parameters of the structures. The new optimizations from the found stationary points did not change the C1 symmetry point group. In addition, it was shown that the asymmetric structure was neither a conical intersection nor the degenerate point of the singlet and triplet states S0 and T1 (avoided crossing). Thus, the DeMore endo-TS structure located earlier12c was not confirmed by the stationary point search with the CASSCF calculations. 3.3. Activation Energy of the C2H4 + O3 → TS Reaction. 3.3.1. General Remarks. From the theoretical point of view, the IRC curve is monotonic with its final point in an arbitrarily small vicinity of the local minimum. The topology of the IRC−CASSCF curve differs from MEP−NEVPT2 curve for the both mechanisms and from MEP−CIPT2 curve for the Criegee mechanism (Figures 2 and 3). The curves correspond-

Figure 4. PES profile corresponding to first step of Criegee and DeMore mechanisms calculated by IRC following method with energy calculations at the CAS(14,11)-MRMP2/cc-pVTZ: 1 and 1′, reagents valleys; TS1 and TS2, saddle points; 2 and 2′, product valleys.

Figure 5. Criegee TS structures calculated by full geometry optimization (a) with energy calculations at the CAS(10,9)/631+G(d,p) (without brackets), RS2C/cc-pVDZ (in brackets), and CCSD(T)/cc-pVTZ [in square brackets] from ref 12a and (b) by IRC following method with energy calculations NEVPT2/cc-pVDZ, CIPT2/cc-pVDZ, and CIPT2+D/cc-pVDZ. Bond lengths are given in Å.

The CAS(10,9) and CAS(14,11) active spaces in combination with the cc-pVDZ, aug-cc-pVDZ, and cc-pVTZ basis sets were employed for the DeMore endo-TS search. These theory levels were successfully used for localization of the Criegee TS. 10426

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

When the PES profile is constructed by the IRC curve following, one should take into account the displacement of the TS and the appearance of local minima for the Criegee and DeMore mechanisms (see Figures 2 and 3; cf., the IRC− CASSCF and renewed IRC curves). For the Criegee mechanism, the shifted saddle points TS−NEVPT2 and TS− CIPT2 have the coordinate −0.580 amu·Bohr−1/2. For the DeMore mechanism, TS−NEVPT2 and TS−CIPT2 have the coordinates −0.300 amu·Bohr −1/2 and 0 amu·Bohr −1/2 (nonshifted). Hereinafter, all values are presented for the IRC following method defined by a renewal of IRC curves by multireference correlated methods, unless indicated otherwise. The activation energy at the CASSCF level is 90 and 120 kJ/ mol for Criegee and DeMore mechanisms, respectively (Figures 2 and 3), which are overestimated values. The IRC CAS(14,11) curves were calculated with retention of orbitals for each increment of the CASSCF curve by means of the GAMESS-US program, i.e., orbitals were updated for the next point. These CASSCF wave functions were used for the next calculation of the energy at the MRMP2/cc-pVTZ level. The special procedure (standard implemented functionality option of GAMESS-US and FireFly programs) made it possible to avoid intruder states. It was shown by the PC-NEVPT2/cc-pVDZ approximation that Ea (relatively to quasi-isolated C2H4 and O3) for the Criegee mechanism is 22.2 kJ/mol, which is 2.0 kJ/mol higher than the corresponding value 24.2 kJ/mol at the RS2C/ccpVDZ level (Table 3). For the strongly contracted method (SC-NEVPT2), this value is 27.6 kJ/mol. For the CIPT2 method, Ea agrees well with the results of the RS2C approximation, being 22.8 kJ/mol, while the Davidsoncorrected value is 20.4 kJ/mol. At the MRMP2/cc-pVTZ level, the activation energy 26.7 kJ/mol is in very good agreement with the values calculated by NEVPT2 methods, particularly, for the strongly contracted formalism. For the DeMore TS, Ea at the RS2C/cc-pVDZ level is 48.6 kJ/mol, which is considerably lower than a value of 91.4 kJ/mol determined by the energy calculation in the IRC points using the PC-NEVPT2 method (Table 3). In the SC-NEVPT2 approximation, this value is 96.4 kJ/mol. The PC-NEVPT2/ccpVDZ calculation in a larger active space CAS(14,11) gave 89.9 kJ/mol, which is only 0.5 kJ/mol lower than Ea for the calculation in the CAS(10,9) active space. The negligibly small change in Ea with the increase of active space at the same order of the perturbation theory for the RS2C and NEVPT2 methods indicates an insufficient coverage of the dynamic electron correlation by the NEVPT2 method as a possible reason of different results. Another reason can be a stronger dependence of the NEVPT2 results on the basis incompleteness. At the CAS(14,11)-MRMP2/cc-pVTZ level, the activation energy of 61.5 kJ/mol for the DeMore reaction channel is in reasonably good agreement with the corresponding RS2C/cc-pVDZ results (Figure 4 and Table 3). The activation energy Ea calculated with the CIPT2 and CIPT2+Dav methods is 66.1 and 62.5 kJ/mol, respectively (Table 3). These values are in excellent agreement with CAS(14,11)-MRMP2/cc-pVTZ values and are better consistent with the RS2C results than with SC- and PC-NEVPT2, which agrees with the proposed explanation of different results by an insufficient coverage of electron correlation by the NEVPT2 method.

ing to the indicated PES profiles contain a minimum presumably corresponding to the prereactive complexes, i.e., the cyclic π-complex corresponds to the Criegee mechanism and the open complex corresponds to the DeMore mechanism. It is worthwhile to note that the method of IRC curve following (in the approximation introduced in ref 45) can have a disadvantage: the probability for the appearance of intruder states is high when profile of the PES fragment is built. In our opinion, this problem can easily be solved by the fact that the shape of the lower-level IRC curve is known and that the intruders affecting the shape of renewed (higher-level) curve can be identified. This property is an advantage over the scheme of the activation energy calculation by full stationary point optimizations for reactants and saddle point since the intruder states can be readily identified for the low-symmetry systems, for example, the DeMore TS (Figure 3), if IRC following method is applied. Thus, the method of the IRC curve following provides an internal test and, in this sense, is self-consistent. The recent application of the methodology46 demonstrates its superiority over the standard approximations of the MEP construction. 3.4. Calculation of Activation Energy by the Method of IRC Curve Following. The activation parameters of elementary reactions calculated by the stationary point search at the RS2C/cc-pVDZ level, those obtained by the IRC curve following at the SC- and PC-NEVPT2/cc-pVDZ, CIPT2/ccpVDZ, and CIPT2+Dav/cc-pVDZ levels and the data for MRMP2/cc-pVTZ theory level are presented in Table 3. The IRC−CASSCF curves and corresponding PES profiles for the Criegee and DeMore mechanisms are shown in Figures 2, 3, and 4 to illustrate the TS shift and the appearance of local minima when electron correlation is taken into account. The total energies of figurative points are presented in the Supporting Information for all computational methods. Table 3. Activation Energies (Ea, kJ/mol) of Elementary Reactions O3 + C2H4 → TS Calculated for Criegee and DeMore TSs by Different Quantum Chemical Approximations with Full Geometry Optimizations (RS2C, CCSD(T) and All the Density Functionals), the IRC Curve Following Method (NEVPT2 and CIPT2), and Composite Approaches (CCSD(T)-F12); the CCSD(T) Data Are Given from ref 12a method/basis set RS2C/cc-pVDZ PC-NEVPT2/cc-pVDZ SC-NEVPT2/cc-pVDZ CIPT2/cc-pVDZ CIPT2+Dav/cc-pVDZ MRMP2/cc-pVTZ CCSD(T)/cc-pVTZ CCSD(T)-F12/VTZ-F12// CCSD(T)/ cc-pVTZ CCSD(T)-F12/VTZ-F12// RS2C/cc− pVDZ CCSD(T)/cc-pVTZ // RS2C/cc-pVDZ OPW91/aug-cc-pVDZ (OPW91/aug-ccpVTZ) OPBE/aug-cc-pVDZ (OPBE/aug-ccpVTZ) OTPSS/aug-cc-pVDZ (OTPSS/aug-ccpVTZ)

Criegee mechanism

DeMore mechanism

22.21 24.22 27.60 22.80 20.44 26.75 12.97 11.90

48.64 91.44 96.35 66.11 62.46 61.47

10.57

55.58

14.60 17.08 (25.08)

58.59 49.66 (57.05)

17.08 (25.27)

49.22 (56.78)

19.58 (27.70)

52.00 (59.34)

10427

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

re(C···O1) = 1.803, 1.798, and 1.855 Å, re(O1−O2) = 1.394, 1.410, and 1.389 Å, re(O2−O3) = 1.256, 1.318, and 1.316 Å, ∠(O−O−O) = 114.65°, 110.49°, and 111.05° for TS-RS2C, DNEVPT2, and D-CIPT2, respectively. Next we will establish the influence of the IRC curve following approximation and performed shift along the IRC curve on geometrical parameters of the stationary points. The geometry optimized at the RS2C/cc-pVDZ level is used as reference geometry since the TSs were found for both mechanisms. The deviation of the geometric parameters of the TS determined by the method of IRC curve following with the energy calculation by the NEVPT2 and CIPT2 methods and obtained by the full geometry optimization can be considered as approximate measure of inaccuracy of the geometric parameter calculated by the IRC following method. Here, we compare data for the uniformly chosen cc-pVDZ basis set. To compare the structures after the shift along the IRC curve when the PES profiles are constructed by the NEVPT2 and CIPT2 energy calculations and obtained by the standard search with full optimization (saddle points TS-RS2C and TS− CASSCF), we calculated the root-mean-square distance (rmsd) for the corresponding geometric parameters of the structures using the ChemCraft program for the following pairs of points: (1) For the Criegee mechanism, TS-RS2C and shifted TS Cg-PT2 (on the one hand) and TS-RS2C and Cg-TSCASSCF (on the other hand); (2) For the DeMore mechanism, DeMore TS-RS2C and shifted TS D-NEVPT2, for TS-RS2C and nonshifted DCIPT2, which is equivalent to D-TS-CASSCF in this case. The rmsd values characterize the overall quality of the localized structures and deviation from the reference ones. For the Criegee TS, rmsd values are 0.017 and 0.034 Å, respectively. The TS shift improves the agreement between the geometric parameters obtained by the full geometry optimization and the method of IRC curve following and indicates that the effect of the dynamic electron correction on the geometric parameters should be taken into account. For

3.4.1. Geometric Parameters of the Transition States. The optimized stuctures of TSs are shown in Figures 5 and 6. The Cartesian coordinates and vibrational frequencies calculated at the RS2C/cc-pVDZ level are presented in the Supporting Information.

Figure 6. DeMore TS structures calculated by full geometry optimization (a) with energy calculations at the CAS(10,9)/631+G(d,p) or CIPT2/cc-pVDZ, CIPT2+D/cc-pVDZ (the TS is not shifted), and RS2C/cc-pVDZ (in brackets) and (b) by IRC following method with energy calculations NEVPT2/cc-pVDZ. Bond lengths are given in Å.

The transition states for the Criegee mechanism coincide for both the NEVPT2 and CIPT2 methods and are denoted as CgPT2 (Figure 2). For the DeMore mechanism, the TS found for the energy calculation by the CIPT2 method coincides with that observed on the PES CAS(10,9)/6-31+G(d), whereas for NEVPT2, it is shifted (Figure 3). The shifted TS is labeled DNEVPT2, the nonshifted TS on the CIPT2 curve is labeled DCIPT2, and the transition states found by the search for stationary points at the RS2C/cc-pVDZ level is labeled TSRS2C. The internuclear distances for the TS corresponding to the cleaved or formed bonds and the bond angle for O3−O2− O1···C are as follows: (1) for the saddle point of the Criegee mechanism (Figure 5), re(C−C) = 1.386 and 1.371 Å, re(C···O) = 2.186 and 2.180 Å, re(O−O) = 1.314 and 1.322 Å, ∠(O−O−O) = 111.11° and 109.77° for TS-RS2C and CgPT2, respectively; (2) for the saddle point of the DeMore mechanism (Figure 6), re(C−C) = 1.401, 1.421, and 1.401 Å,

Table 4. Arrhenius Parameters, EA (kJ/mol) and Pre-Exponential Factor A (L/(mol·s)), Canonical Rate Constant (k(T), L·(mol·sec)−1) at the Temperature 298 K for Elementary Reaction O3 + C2H4 Proceeds via Criegee TS (k1) and DeMore TS (k2) Calculated on the Basis of Conventional TST and Quantum Chemical Computations (RS2C, NEVPT2, and CIPT2 Methods in Conjunction with cc-pVDZ Basis Set, MRMP2 and CCSD(T) in Conjunction with cc-pVTZ Basis Set, and Composite Approaches with Single-Point Energy Refinement at the CCSD(T)-F12/VTZ-F12 Level) Criegee mechanism

a

DeMore mechanism

method/basis set

A (10−6)

EA

k1

A (10−6)

EA

RS2C PC-NEVPT2 SC-NEVPT2 CIPT2 CIPT2+Dav MRMP2 CCSD(T) CCSD(T)-F12/CCSD(T) CCSD(T)-F12/RS2C OPW91/AVTZ OPBE/AVTZ OTPSS/AVTZ

1.54 1.54 1.54 1.54 1.54 6.49 0.95 0.95 1.54 2.11 2.18 2.01

31.56 32.57 35.95 31.15 28.79 43.58 21.79 20.72 18.92 33.59 33.77 36.22

50.4 22.38 5.72 39.70 102.92 0.179 1.08 × 103 1.66 × 103 5.54 × 103 511.4 19.47 176.5

1.27 1.27 1.27 1.27 1.27 58.90 a a 1.27 2.97 2.96 2.98

58.55 101.35 106.26 76.02 72.37 71.68

5.19 1.63 2.25 4.49 1.96 1.87

65.49 63.09 62.83 65.35

3.14 × 10−5 3.8 × 10−3 2.12 × 10−4 1.5 × 10−3

k2 × × × × × ×

k1/k2 10−4 10−11 10−12 10−7 10−6 10−5

105 1012 1012 108 107 104

108 105 105 105

Stationary point does not exist at this theory level. 10428

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

DeMore TS, rmsd for the first pair of structures is 0.144 Å, and the second rmsd is 0.144 Å, while in the case of the open structure of the transition state, this improvement is insufficient. 3.5. Multireference vs Single-Reference Methods: The Overall Data Comparison. The results of the TS energy calculations by multireference, coupled cluster, and density functional methods and the rate constants for Criegee and DeMore mechanisms are presented in Tables 3 and 4, respectively. The TS1 structure (Scheme 1) was repeatedly analyzed earlier12 (and references therein). It is characterized by equal C···O distances (RCO) in each pair of atoms. The values of RCO bond lengths are 1.97 (MP2), 2.2 (CCSD, QCISD), and 2.3 (B3LYP) Å (the data are given in completeness from ref 27). Only the results for MP2 method deviate significantly from both high-correlated coupled cluster method and multireference perturbation theory. From the general considerations, two configurations are possible for TS2 with different orientations of the ozone moiety relative to the double bond (Scheme 2): trans-configuration of TS2 and cis-configuration of TS2. It should be mentioned that the existence of exo-TS2 is confirmed at all the computational levels, whereas endo-TS2 is not always observed even for different density functionals. Moreover, it was shown by CASSCF calculations that the endo-TS was an artifact of the density functional. If both exo-TS2 and endo-TS2 are observed, they are rather close in energy and all conclusions for them will be identical, but the cis-configuration of TS2 usually (but not always) is more favorable. The activation energy and rate constants calculated in ref 27 using both single-reference methods (DFT and ab initio WFT) with different basis sets (6-31+G(d,p), 6-311+G(d,p), and augcc-pVDZ) are in reasonable good agreement with the best estimates12a for Ea of the Criegee channel particularly for OLYP density functional. The performance of the OLYP method is consistent with the results of higher level quantum chemical methods.12b One can conjecture that, if density functional describes one reaction channel (the first step of Criegee mechanism) satisfactory, this functional can describe the other reaction channel (the first step of DeMore mechanism) as well. Indeed, comparison of the data collected in ref 27 and in Table 4, particularly those based on RS2C/cc-pVDZ optimized geometry, leads us to conclude that OLYP density functional describes both channels with similar quality of high-correlated or multireference approaches. From DFT functional assessment,12b two functionals (HCTH60 and O3LYP61) are characterized by deviation of activation energy from the best estimation, which is smaller than this value for other DFT functionals including OLYP functional. Encouraged by good performance of density functional based on Handy’s OPTX61 modification of Becke’s exchange functional, namely, OLYP and O3LYP, and Handy’s HCTH (called HCTH/407 in Gaussian03 and Gaussian09), we carried out calculations with the OPTX exchange and PW91, PBE, or TPSS correlation functionals, i.e., OPW91, OPBE, and OTPSS density functionals, in conjunction with aug-cc-pVDZ and augcc-pVTZ basis sets (Table 4). All of them provide very close values of activation energies, i.e., about 17 (25) kJ/mol for the Criegee TS and 50 (60) kJ/mol for the DeMore TS; data for the aug-cc-pVTZ basis set are given in parentheses. The geometric parameters of the structures are consistent (cf., Supporting Information). Hence, the nonobvious and nonunique combination of exchange and correlation functionals

provided surprisingly good agreement (near to chemical accuracyof ∼4 kJ/mol) with the best estimate of activation energy value for the Criegee mechanism and calculated by means of high-correlated methods in the presented study activation energy for the DeMore mechanism if the moderate aug-cc-pVDZ basis set was employed. The observed good agreement is based on fortuitous compensation of basis set superposition error (BSSE) and approximate description of the high correlated system. With the more flexible aug-cc-pVTZ basis set, the activation energy of the Criegee reaction channel is overestimated by about 12 kJ/mol; this value is a more common error of density functionals. The rate constant ratios are in the range between 104 and 108, with most values being close to the upper boundary. Let us now consider the data in Table 3 and focus our attention on the results of ab initio methods. First of all, one may note that, for the Criegee mechanism, the calculated activation energies for all multireference methods range from 26 to 38 kJ/mol, i.e., they are fairly close to each other but noticeably exceed the experimental values and the best available estimations by CBS/focal point methodology.12a,b The MRMP2(14,11)/cc-pVTZ data, the results of the RS2C calculations, and the both CIPT2 approximations correspond best of all to the latter. To resolve the observed disagreement between different formalisms of the multireference perturbation theory, which can suffer from nonbalanced description of the electron correlation effects and show high sensitivity to the noncomplete one-electronic basis set, we carried out additional energy calculations at the CCSD(T)-F12/VTZ-F12 theory levels for the structures optimized at the RS2C/cc-pVDZ (both mechanisms) and CCSD(T)/cc-pVTZ (the Criegee mechanism) theory levels. For completeness, the conventional CCSD(T) energy calculations were performed at the RS2C geometries. Activation energies are 10.6 (11.9) kJ/mol for Criegee mechanism and 55.6 kJ/mol for DeMore mechanism, where data in brackets are results of CCSD(T) geometry optimization (Table 3). The agreement between activation energies calculated by RS2C, CIPT2, and MRMP2 (multireference approaches) and CCSD(T)-F12//RSC2 (single-reference approach) validates the single-reference explicitly correlated CCSD(T)-F12 method for the highly correlated system and reveals the energy overestimation at the NEVPT2/cc-pVDZ level. The CCSD(T)-F12/VTZ-F12//RS2C approach combines advantages of explicitly correlated method with very flexible basis set and feasible optimization at the RS2C/ccpVDZ calculations. The CCSD(T)-F12/VTZ-F12//RS2C is the best available approach for the small systems, whereas the results of single-point energy calculation by the CCSD(T)// RS2C method are unreliable due to high influence of nondynamic electron correlation in the system considered. The results of calculations by the employed density functionals are in reasonably good agreement with RS2C, MRMP2, and CIPT2 methods. Slightly worse agreement is observed in comparison between the CCSD(T)-F12//CCSD(T) and CCSD(T)-F12//RSC2 composite methods for the Criegee channel; however, for the particularly difficult for single-reference method DeMore TS this agreement is encouraging. The OPW91, OPBE, and OTPSS density functionals are competitive with the best available estimates (Table 3). The RS2C and CIPT2 (with frequencies calculated at the CAS(10,9)-RS2C/cc-pVDZ level) methods provide rate 10429

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

Table 5. Stabilization Energies Ec (kJ/mol) of the Prereactive Complexes of Criegee and DeMore Mechanisms Calculated by Constraint Optimization at RS2C/cc-pVDZ, estimated by CCSD(T)/aug-cc-pVTZ//CCSD(T)/cc-pVTZ (1) and CCSD(T)F12/VTZ-F12/RS2C//cc-pVDZ (2) Compound Approaches, IRC Following Method with Energy Calculation by PC- and SCNEVPT2, CIPT2, and CIPT2+Dav Methods with cc-pVDZ Basis Set and Active Space CAS(10,9), and at the MRMP2(14,11)/ cc-pVTZ Theory Level stabilization energies (Ec, kJ/mol) of the prereactive Criegee and DeMore complexes calculated by different quantum chemical methods CCSD(T)−F12

a

NEVPT2

complex

RS2C

(1)

(2)

PC

SC

CIPT2

CIPT2+Dav

MRMP2

Criegee mechanism DeMore mechanism

3.19 1.96

8.50 8.03

12.01 6.16

2.83 0.96

2.31 0.73

3.79 2.62

4.01 1.72

1.8 a

Stationary point does not exist at this theory level.

Figure 7. Structures of the prereactive (a) Criegee complex optimized at the RS2C/cc-pVDZ and CCSD(T)/cc-pVTZ [in square brackets] theory levels and (b) the DeMore complex optimized at the RS2C/cc-pVDZ and (c) CCSD(T)/cc-pVTZ theory levels. Bond lengths are given in Å. For the Criegee complex, the CCSD(T) data are from ref 12a.

constant values for the Criegee mechanism (Table 4), which are in reasonable agreement with the reviewed7 data (k(T) = 9.1 × 10−15 [cm3/molecule s] exp(−21.45 [±0.83 kJ/mol]/ RT) = 5.48 × 106 [L·(mol·s)−1] exp(−21.45 [±0.83 kJ/mol]/ RT), and k(298.0) = 1.58 × 10−18 [cm3/molecule s] = 9.51 × 102 L·(mol·s)−1), experimental data62 k(T)= 5.1 × 10−15 [±1.01 × 10−15 cm3/molecule s] exp(−20.37 [±0.83 kJ/ mol]/RT), and the most recent63 k(298 K) = (1.45 ± 0.25) × 10−18 [cm3/molecule s] = (873 ± 145 L·(mol·s)−1) and differ from the last one by 20 times due to the higher activation energy. For the DeMore mechanism, activation energies at the NEVPT2 levels are overestimated, and hence, the ratio of the reaction rate constants for the Criegee and DeMore channels is very high (up to ∼1012), whereas for the RS2C and CIPT2 methods it is not more than ∼108 (Table 4). The MRMP2(14,11)/cc-pVTZ theory level gives reasonable values of Ea (61.5 kJ/mol), and the k1/k2 ratio is about 104. The results of the single-reference approaches show that the coupled cluster method in different variants (CCSD(T) and CCSD(T)-F12) in conjunction with the flexible basis sets provides the results closely corresponding to the experimental data:7,62,63 for the Criegee mechanism, the calculated activation energy is 10−12 kJ/mol and k1 = ∼103 L·(mol·s)−1 (Tables 3 and 4). The best estimated value for the activation energy of the DeMore reaction channel is about 55 kJ/mol (not more than 66.1 kJ/mol determined by CIPT2/cc-pVDZ approach) as obtained by compound approach CCSD(T)-F12/VTZ-F12// RS2C/cc-pVDZ (Table 3), which combines rigorous geometry optimization with RS2C method and state-of-the-art explicitly correlated CCSD(T)-F12 method near the basis set limit for energy estimation. Hence, the CCSD/aug-cc-pVDZ and QCISD/aug-cc-pVDZ methods applied in ref 27 also satisfactory describe both reaction channels but, in our opinion,

underestimate significantly activation energy for the DeMore mechanism. For the reason of high calculated value of the T1diagnosticis (0.05), the agreement between experimental data and those calculated with CCSD and especially with QCISD methods is not reliable. The composite approach CCSD(T)-F12/VTZ-F12//RS2C/ cc-pVDZ provides the best estimate for the activation energy (Ea) on DeMore reaction pathway 55.58 ± 2.40 kJ/mol. The uncertainty (2.40 kJ/mol) was estimated by difference of Ea for the Criegee mechanism at the CCSD(T)/cc-pVTZ and CCSD(T)-F12/VTZ-F12//RS2C/cc-pVDZ and aroused by structures located at the RS2C/cc-pVDZ level nonoptimized at the CCSD(T)-F12/VTZ-F12 level (Table 3). For the Criegee mechanism, the best estimation of the central value of activation energy is 11.90 kJ/mol, i.e., the value obtained by the CCSD(T)-F12/VTZ-F12//CCSD(T)/cc-pVTZ composite approach since the CCSD(T)/cc-pVTZ theory level provides geometry of stationary points with higher precision than the RS2C/cc-pVDZ level. Thus, the upper limit of the uncertainty aroused by the lack of geometry optimization at the CCSD(T)F12/VTZ-F12 theory level can be estimated as the absolute value of the difference (2.37 kJ/mol) of activation energies (12.97 and 14.60 kJ/mol) calculated by CCSD(T) method at the CCSD(T) and RS2C optimized geometries (Table 3). Therefore, the best estimation obtained in the presented study for Ea of the Criegee mechanism is (11.90 ± 2.37) kJ/mol; that is in agreement with the intervals (14.26 ± 0.84) kJ/mol or (14.63 ± 3.00) and (13.30 ± 3.00) kJ/mol reported previously.12a,b Two ways of error estimation employed in the presented study provide very close values. The most affordable DFT can be used for the description of this reaction as the first approximation when special care is taken on the choice of functionals. 10430

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

3.6. Prereactive Complexes. 3.6.1. Geometrical Parameters and Stabilization Energy. Let us consider the appearance of the local minimum on the PES profiles constructed by the method of IRC curve following (Figures 2 and 3) and estimate the stabilization energy of the prereactive complexes for each reaction mechanism (Table 5). For the NEVPT2 and CIPT2 calculations, the IRC following method reveals complexes for both mechanisms; for the MRMP2 calculation, the complex is observed only for Criegee mechanism. We choose as a reference values those geometrical parameters optimized at RS2C/cc-pVDZ (see Figure 7). The structures determined at the NEVPT2 and CIPT2 theory levels are introduced in Figure 8.

stabilization energy does not exceed 5 kJ/mol, and for the second one, it is about 3 kJ/mol. The relative energies of the Criegee intermediate are −12.0 and −8.5 kJ/mol for RS2C and CCSD(T) optimized structures with CCSD(T)-F12/VTZ-F12 single-point energy calculation, respectively. For the prereactive Criegee complex, the relative energy is satisfactorily consistent with a reference value of −8.5 kJ/mol (−2.03 kcal/mol) calculated by the focal point analysis method with extrapolation to the infinite (complete) basis set limit (CBS CCSD(T)+Q).12a It is interesting to note that the extrapolated (focal point/CBS) value is in excellent agreement with the results of single-point calculation. To the best of out knowledge, it is a first comparison of the explicitly correlated coupled cluster method to the focal point analysis with conventional coupled cluster approach. It illustrates the dramatic improvement of the results for the F12 coupled cluster approach recently introduced52 and efficiently implemented in MOLPRO 2010. 3.6.2. Prereactive Complex. We performed full geometry optimization for the DeMore intermolecular complex at the CCSD(T,FC)/cc-pVTZ theory level with 10−6 a.u. convergence criteria of the stationary point search by means of CFOUR program and optimization at the CCSD/cc-pVTZ level performed with the Gaussian09 suite of programs where the “Tight” convergence criteria was used. The converged optimization was accomplished by vibrational frequency calculation. The structure (see Figure 7) was determined as a true local minimum (see Supporting Information for CCSD data and frequencies at the CCSD(T,FC)/cc-pVTZ theory level). The calculated stabilization energies of Criegee and DeMore complexes (Table 5) at the CCSD(T)/cc-pVTZ and CCSD(T)/aug-cc-pVTZ (in parentheses) theory levels are −9.61 (−10.16) and −9.47 (−9.74) kJ/mol, respectively. The relative stability of the complexes is particularly close for the CCSD(T)-F12/CCSD(T) compound approach. The ZPVE corrections (CCSD(T)/cc-pVTZ) estimated in the harmonic oscillator approximation are very close (153.86 vs 153.93 kJ/ mol) for the complexes and have no effect on the relative stability. Calculated shifts for the strongest ν2 ozone band determined by relatively weak O3···C2H4 interaction are −7.3 and 1.0 cm−1 for Criegee complex and −8.2 and −4.3 cm−1 for DeMore complex. Experimentally determined9f shifts are −2.6 and −1.3 cm−1. Calculated shifts for ν1 vibrational mode of O3 are 1.8 and 10.0 cm −1 for Criegee and DeMore complexes, respectively. The reported9f absence of new features of the weak ν1 vibrational frequency is in better agreement with the 1.8 cm−1 shift for Criegee complex. Thus, the experimentally observed IR shifts are in the slightly better agreement with the Criegee complex rather than the complex of the DeMore mechanism. Both complexes exhibit very close shifts for vibrational modes of ethylene. Barrierless formation and very close calculated relative energies (Table 5) can determine a probability to form a mixture. The matrix isolation study,9f frequencies assignment and comparison of calculated and experimental patterns pointed out the existence of the reported here DeMore complex. We also performed the geometry optimization of the located stationary structures at the B3LYP/ cc-pVTZ and B3LYP/aug-cc-pVTZ levels. All the optimizations converged to the Criegee complex. On the basis of the above results, we conclude that, although the values of stabilization energies and IR shifts are comparable with typical uncertainty

Figure 8. Structures of the prereactive Criegee (a) and DeMore (b) complexes optimized by IRC following method at the NEVPT2/ccpVDZ and CIPT2/cc-pVDZ [in square brackets] theory levels. Bond lengths are given in Å.

The shortest internuclear distances re(C···O) are as follows: (1) for the prereactive complex of the Criegee mechanism, re(C···O) = 3.105 and 3.026 Å for Cg-NEVPT2 and Cg-CIPT2, respectively; (2) for the prereactive complex of the DeMore mechanism, re(C···O) = 3.258 and 3.008 Å for D-NEVPT2 and D-CIPT2, respectively. The rmsd values calculated for a comparison of the obtained structures of the Criegee prereactive complex are as follows: 0.106 and 0.111 Å for the RS2C-NEVPT2 and RS2C-CIPT2 pairs, respectively. The data show somewhat better agreement between the geometric parameters obtained by full optimization at the RS2C/cc-pVDZ level and by the method of IRC curve following with the energy calculation by the NEVPT2 method rather than by the CIPT2 method. For the prereactive DeMore complex, rmsd is 0.857 Å in the RS2C-NEVPT2 pair, and in the pair RS2C−CIPT2, it is 0.823 Å. Both are substantially worse than the corresponding values for the prereactive Criegee complex at these theory levels. Thus, the complex of the DeMore mechanism is considerably more loose, and its structural parameters cannot be determined by the IRC curve following method with high precision. However, its existence can be claimed at the different theory levels, and its structure is reported here for the first time. We demonstrated an important advantage of IRC curve following method for pre-equilibrium stage determination. The stabilization energies of the prereactive complexes are consistent with each other (see Table 5). This result can be explained by small gradients on the PES fragments corresponding to intermolecular complexes. The qualitative dependence is as follows: the prereactive complex of the Criegee mechanism is energetically more favorable than the complex for the DeMore mechanism. For the first intermolecular complex, the 10431

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

of the quantum chemical theory and the unambiguous assignment of the experimentally registered structure is still disputable, the new high level calculations indicate the prevalence of the Criegee complex.

A new 1:1 complex of ozone−ethylene was located and characterized by the IR shifts calculated at the CCSD(T)/ccpVTZ theory level. The relative stability of the Criegee and DeMore complexes was estimated by the CCSD(T)/aug-ccpVTZ//CCSD(T)/cc-pVTZ and CCSD(T)-F12/VTZ-F12// CCSD(T)/cc-pVTZ composite approaches. The CCSD(T) calculations show that the experimentally determined IR shifts in the inert matrix isolation study are in a somewhat better agreement with the structure of the Criegee complex. Since the experimental data discussed here were obtained in 80ths, further matrix investigations of complexes formed in the O3/ C2H4 mixture are highly desirable.

4. CONCLUSIONS Multireference calculations of the PES for the C2H4 + O3 system confirm the presence of two competing reaction channels. The active space for this problem should be at least CAS(10,9), while CAS(14,11) is better. All the data obtained here indicate that the concerted [2 + 3]-cycloaddition mechanism prevails for ethylene. The ratio of the reaction rates for two channels in these calculations are in range ∼104− 108 (with preference of the upper bound of the interval) in favor of the Criegee mechanism. The results of calculations at the highest multireference theory level with the full geometry optimization employed here for a first time (CASPT2 method in the RS2C formalism in conjunction with cc-pVDZ basis set) and the results of duallevel approximation denoted as the IRC curve following method with the energy calculation by various methods (MRMP2, NEVPT2, and CIPT2) along IRC curve reveal consistency of two methodologies. The CIPT2 and MPMP2 results are much closer to the RS2C results than the NEVPT2 data for DeMore mechanism. For the Criegee mechanism, the assessed multireference methods provide the results close to the RS2C/cc-pVDZ level. The CIPT2 approach is in somewhat better agreement with the last one. As a whole, the MPMP2 approximation (the least expensive one, even with the flexible cc-pVTZ basis set and the largest active space CAS(14,11)) gives results of a similar quality to the RS2C method applied with moderate cc-pVDZ basis set, i.e., the IRC following method in combination with MRMP2/cc-pVTZ is reliable and a sufficiently accurate approach for quantitative estimations. This conclusion was supported by the results of the CCSD(T)F12//RS2C and CCSD(T)-F12//CCSD(T) composite methods. The IRC curve following method is suitable for weak prereactive complex localization, e.g., for complexes of O3 molecule which is notoriously difficult system to provide adequate results by application of quantum chemical methods. It was evidenced that this dual-level approximation is a rather economical tool to perform activation energy calculations for high-correlated systems with internal test of the procedure and ability to increase significantly the size of systems. The explicitly correlated coupled cluster method (CCSD(T)-F12) provides a new standard of accuracy in quantum chemical calculations. The composite approach CCSD(T)-F12/VTZ-F12//RS2C/ccpVDZ provides the best estimation (55.58 ± 2.40 kJ/mol) of activation energy (Ea) of the DeMore reaction pathway. The estimated Ea for the Criegee mechanism is (11.90 ± 2.37) kJ/ mol; that is in agreement with the results of previous high-level studies. Multireference methods are more rigorous for systems where nondynamic correlation is important; however, we notice that even single-reference methods are quite adequate for this problem. Of the single-reference methods, the coupled cluster approach gives the best results that close to the results of multireference approximation results, whereas the DFT methods (OPW91, OPBE, and OTPSS) can be used for the preliminary estimates or for systems of experimental interest after the presented calibration on the model system.

■

ASSOCIATED CONTENT

S Supporting Information *

Calculated Cartesian coordinates, IRC curves, total energies, vibrational frequencies, and complete citations for references 35, 38, 41, and 42. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected] (O.B.G.); [email protected] (B.E.K). Present Address ⊥

IBM Almaden Research Center, San Jose, California 951206099, USA.

Author Contributions ⊥

These authors contributed equally to this work.

Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work was supported in part by the U.S. National Science Foundation (CHE-0832622). O.B.G., A.E.M., and S.G. are grateful to DOE NERSC, I2lab, and Institute for Simulation and Training (IST) for computer time provided. S.K.I. acknowledges the RFBR project 11-03-00085. O.B.G. is thankful to Dr. Prof. Hans Lischka (Institut für Theoretische Chemie und Strahlenchemie, Universität Wien, A-1090 Vienna, Austria) and Dr. Prof. Thomas Müller (Research Centre Jülich GmbH, D-52425 Jülich, Germany) for helpful discussion of methodological questions.

■

REFERENCES

(1) (a) McGillen, M. R.; Ghalaieny, M.; Percival, C. J. Phys. Chem. Chem. Phys. 2011, 13, 10965−10969. (b) Bateman, A. P.; Nizkorodov, S. A.; Laskin, J.; Laskin, A. Phys. Chem. Chem. Phys. 2009, 11, 7931− 7942. (2) (a) Finlayson, B. J.; Pitts, J. N., Jr.; Akimoto, H. Chem. Phys. Lett. 1972, 12, 495−498. (b) Arora, P. K.; Vohra, K. G. Chem. Phys. Lett. 1983, 96, 161−166. (3) (a) Ohira, S.-I.; Li, J.; Lonneman, W. A.; Dasgupta, P. K.; Toda, K. Anal. Chem. 2007, 79, 2641−2649. (b) Marley, N. A.; Gaffney, J. S. Atmos. Environ. 1998, 32, 1435−1444. (4) George, G.; Goss, B. Elastomers in Oxygen and Ozone, Degradation of. In Encyclopedia of Materials: Science and Technology, 2nd ed.; Buschow, K. H. J., Robert, W. C., Merton, C. F., Bernard, I., Edward, J. K., Subhash, M., Patrick, V., Eds.; Elsevier: Oxford, U.K., 2001; pp 2466−2469. (5) (a) Wang, L. K.; Shamas, N. K. Pressurized Ozonation. In Biosolids Treatment Processes; Wang, L. K., Shammas, N. K., Hung, Y.T., Eds.; Humana Press: New York, 2007; Vol. 6, pp 243−297. (b) von Gunten, U. Water Res. 2003, 37, 1443−1467. (c) von Gunten, U. Water Res. 2003, 37, 1469−1487. 10432

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

(26) Chan, W.-T.; Hamilton, I. P. J. Chem. Phys. 2003, 118, 1688− 1701. (27) Krisyuk, B.; Maiorov, A. Russ. J. Phys. Chem. B 2011, 5, 790− 796. (28) Maiorov, A. V.; Krisyuk, B. E.; Popov, A. A. Russ. J. Phys. Chem. B 2008, 2, 707−710. (29) Ljubić, I.; Sabljić, A. J. Phys. Chem. A 2002, 106, 4745−4757. (30) Li, L.-C.; Deng, P.; Xu, M.-H.; Wong, N.-B. Int. J. Quantum Chem. 2004, 98, 309−316. (31) East, A. L. L.; Allen, W. D. J. Chem. Phys. 1993, 99, 4638−4650. (32) Gadzhiev, O. B.; Ignatov, S. K.; Gangopadhyay, S.; Masunov, A. E.; Petrov, A. I. J. Chem. Theory Comput. 2011, 7, 2021−2024. (33) (a) Gadzhiev, O. B.; Ignatov, S. K.; Razuvaev, A. G.; Masunov, A. E. J. Phys. Chem. A 2009, 113, 9092−9101. (b) Zakharov, I.; Minaev, B. J. Struct. Chem. 2012, 53, 1−11. (34) Müller, T.; Xantheas, S. S.; Dachsel, H.; Harrison, R. J.; Nieplocha, J.; Shepard, R.; Kedziora, G. S.; Lischka, H. Chem. Phys. Lett. 1998, 293, 72−80. (35) Frisch, M. J.; et al. Gaussian 03, version D.02; Gaussian, Inc.: Wallingford, CT, 2007. (36) 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.; Windus, T. L.; Dupuis, M.; Montgomery, J. A. J. Comput. Chem. 1993, 14, 1347−1363. (37) Granovsky, A. A. Firefly, version 7.1. G; see http://classic.chem. msu.su/gran/firefly/index.html (accessed July 2012). (38) Stanton, J. F.; et al. CFOUR. Integral packages MOLECULE (Almlöf, J.; Taylor, P. R.), PROPS (Taylor, P. R.), ABACUS (Helgaker, T.; Jensen, H. J. Aa.; Jørgensen, P.; Olsen, J.), and ECP routines by Mitin, A. V.; van Wüllen, C. For the current version, see http://www.cfour.de (accessed June 2012). (39) Celani, P.; Stoll, H.; Werner, H.-J.; Knowles, P. J. Mol. Phys. 2004, 102, 2369−2379. (40) (a) Angeli, C.; Cimiraglia, R.; Evangelisti, S.; Leininger, T.; Malrieu, J. P. J. Chem. Phys. 2001, 114, 10252−10264. (b) Angeli, C.; Evangelisti, S.; Cimiraglia, R.; Maynau, D. J. Chem. Phys. 2002, 117, 10525−10533. (c) Havenith, R. W. A.; Taylor, P. R.; Angeli, C.; Cimiraglia, R.; Ruud, K. J. Chem. Phys. 2004, 120, 4619−4625. (41) Werner, H.-J.; et al. MOLPRO, version 2010.1, a package of ab initio programs; see http://www.molpro.net (accessed June 2012). (42) Frisch, M. J.; et al. Gaussian 09, version C.0.1; Gaussian, Inc.: Wallingford CT, 2011. (43) Dyall, K. G. J. Chem. Phys. 1995, 102, 4909−4918. (44) (a) Hirao, K. Chem. Phys. Lett. 1992, 190, 374−380. (b) Hirao, K. Chem. Phys. Lett. 1992, 196, 397−403. (45) West, A. C.; Kretchmer, J. S.; Sellner, B.; Park, K.; Hase, W. L.; Lischka, H.; Windus, T. L. J. Phys. Chem. A 2009, 113, 12663−12674. (46) Kurtén, T.; Donahue, N. M. J. Phys. Chem. A 2012, 116, 6823− 6830. (47) Malick, D. K.; Petersson, G. A.; Montgomery, J. J. A. J. Chem. Phys. 1998, 108, 5704−5713. (48) Chuang, Y.-Y.; Corchado, J. C.; Truhlar, D. G. J. Phys. Chem. A 1999, 103, 1140−1149. (49) (a) Celani, P.; Werner, H.-J. J. Chem. Phys. 2003, 119, 5044− 5057. (b) Celani, P.; Werner, H.-J. J. Chem. Phys. 2000, 112, 5546− 5557. (c) Werner, H.-J. Mol. Phys. 1996, 89, 645−661. (50) Sun, J.-Q.; Ruedenberg, K. J. Chem. Phys. 1994, 101, 2157− 2167. (51) Császár, P.; Pulay, P. J. Mol. Struct. 1984, 114, 31−34. (52) (a) Adler, T. B.; Knizia, G.; Werner, H.-J. J. Chem. Phys. 2007, 127, 221106−221104. (b) Knizia, G.; Adler, T. B.; Werner, H.-J. J. Chem. Phys. 2009, 130, 054104−054120. (53) Zhurko, G. ChemCraft, version 1.6 (build 348); see http://www. chemcraftprog.com (accessed March 20, 2012) (54) Ignatov, S. K. MOLTRAN v.2.5, Program for Molecular Visualization and Thermodynamic Calculations; University of Nizhny Novgorod: Nizhny Novgorod, Russia, 2007; see http://ichem.unn.ru/ Moltran (accessed March 20, 2012).

(6) (a) Søndergaard, R.; Nielsen, O. J.; Hurley, M. D.; Wallington, T. J.; Singh, R. Chem. Phys. Lett. 2007, 443, 199−204. (b) Wegener, R.; Brauers, T.; Koppmann, R.; Rodríguez Bares, S.; Rohrer, F.; Tillmann, R.; Wahner, A.; Hansel, A.; Wisthaler, A. J. Geophys. Res. 2007, 112, D13301. (c) Wang, H.; Zhang, Y.; Mu, Y. J. Environ. Sci. 2007, 19, 641−643. (7) Atkinson, R.; Baulch, D. L.; Cox, R. A.; Crowley, J. N.; Hampson, R. F.; Hynes, R. G.; Jenkin, M. E.; Rossi, M. J.; Troe, J. Atmos. Chem. Phys. 2004, 4, 1461−1738. (8) Sander, S. P.; Abbatt, J.; Barker, J. R.; Burkholder, J. B.; Friedl, R. R.; Golden, D. M.; Huie, R. E.; Kolb, C. E.; Kurylo, M. J.; Moortgat, G. K.; Orkin, V. L.; Wine, P. H. Chemical Kinetics and Photochemical Data for Use in Atmospheric Studies, Evaluation No. 17, JPL Publication 10−6,; Jet Propulsion Laboratory: Pasadena, CA, 2011; see http://jpldataeval.jpl.nasa.gov (accessed July 2012). (9) (a) Kuehne, H.; Guenthard, H. H. J. Phys. Chem. 1976, 80, 1238− 1247. (b) Kohlmiller, C. K.; Andrews, L. J. Am. Chem. Soc. 1981, 103, 2578−2583. (c) Andrews, L.; Kohlmiller, C. K. J. Phys. Chem. 1982, 86, 4548−4557. (d) Feltham, E. J.; Almond, M. J.; Marston, G.; Ly, V. P.; Wiltshire, K. S. Spectrochim. Acta, Part A 2000, 56, 2605−2616. (e) Bariseviciute, R.; Ceponkus, J.; Sablinskas, V. Cent. Eur. J. Chem. 2007, 5, 71−86. (f) Nord, L. J. Mol. Struct. 1982, 96, 37−44. (g) Nelander, B.; Nord, L. Tetrahedron Lett. 1977, 18, 2821−2822. (h) Samuni, U.; Haas, Y.; Fajgar, R.; Pola, J. J. Mol. Struct. 1998, 449, 177−201. (10) (a) Kühne, H.; Vaccani, S.; Ha, T. K.; Bauder, A.; Günthard, H. H. Chem. Phys. Lett. 1976, 38, 449−455. (b) Kühne, H.; Vaccani, S.; Bauder, A.; Günthard, H. H. Chem. Phys. 1978, 28, 11−29. (c) Zozom, J.; Gillies, C. W.; Suenram, R. D.; Lovas, F. J. Chem. Phys. Lett. 1987, 140, 64−70. (11) (a) Cremer, D. J. Chem. Phys. 1979, 70, 1898−1910. (b) Cremer, D. J. Chem. Phys. 1979, 70, 1911−1927. (c) McKee, M. L.; Rohlfing, C. M. J. Am. Chem. Soc. 1989, 111, 2497−2500. (d) Hull, L. A.; Hisatsune, I. C.; Heicklen, J. J. Am. Chem. Soc. 1972, 94, 4856−4864. (12) (a) Wheeler, S. E.; Ess, D. H.; Houk, K. N. J. Phys. Chem. A 2008, 112, 1798−1807. (b) Zhao, Y.; Tishchenko, O.; Gour, J. R.; Li, W.; Lutz, J. J.; Piecuch, P.; Truhlar, D. G. J. Phys. Chem. A 2009, 113, 5786−5799. (c) Saito, T.; Nishihara, S.; Kataoka, Y.; Nakanishi, Y.; Kitagawa, Y.; Kawakami, T.; Yamanaka, S.; Okumura, M.; Yamaguchi, K. J. Phys. Chem. A 2010, 114, 12116−12123. (13) Nagy, T.; Turányi, T. Int. J. Chem. Kinet. 2011, 43, 359−378. (14) Copeland, G.; Ghosh, M. V.; Shallcross, D. E.; Percival, C. J.; Dyke, J. M. Phys. Chem. Chem. Phys. 2011, 13, 14839−14847. (15) Criegee, R. Angew. Chem., Int. Ed. 1975, 14, 745−752. (16) (a) Razumovskii, S. D.; Zaikov, G. E. Ozone and Its Reactions with Organic Compounds (Studies in Organic Chemistry); Elsevier Science Ltd: New York, 1984. (b) Razumovskii, S. D.; Gennadii, E. Z. Russ. Chem. Rev. 1980, 49, 1163−1180. (c) Atkinson, R.; Carter, W. P. L. Chem. Rev. 1984, 84, 437−470. (17) Gillies, J. Z.; Gillies, C. W.; Suenram, R. D.; Lovas, F. J.; Stahl, W. J. Am. Chem. Soc. 1989, 111, 3073−3074. (18) Gillies, C. W.; Gillies, J. Z.; Suenram, R. D.; Lovas, F. J.; Kraka, E.; Cremer, D. J. Am. Chem. Soc. 1991, 113, 2412−2421. (19) Wadt, W. R.; Goddard, W. A. J. Am. Chem. Soc. 1975, 97, 3004− 3021. (20) DeMore, W. B. Int. J. Chem. Kinet. 1969, 1, 209−220. (21) Zhang, X.-M.; Zhu, Q. J. Org. Chem. 1997, 62, 5934−5938. (22) Tyupalo, N.; Zakharov, I.; Kolbasina, O.; Zhidomirov, G.; Avdeev, V. J. Struct. Chem. 2000, 41, 194−199. (23) (a) Kuczkowski, R. L. Chem. Soc. Rev. 1992, 21, 79−83. (b) Jonnalagadda, S.; Chan, S.; Garrido, J.; Bond, J.; Singmaster, K. A. J. Am. Chem. Soc. 1995, 117, 562−563. (24) O’Neal, H. E.; Blumstein, C. Int. J. Chem. Kinet. 1973, 5, 397− 413. (25) Krisyuk, B. E.; Popov, A. A. Addition of Ozone to Multiple Bonds. Competition of the Reaction Pathways. In Chemical and Biochemical Physics, Kinetics and Thermodynamics: New Perspectives; Stott, P. E., Zaikov, G. E., Kablov, V. F., Eds.; Nova Science Publisers Inc.: Hauppauge, NY , 2008; pp 256. 10433

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434

The Journal of Physical Chemistry A

Article

(55) Dennington, R.; Keith, T.; Millam, J. GaussView; Semichem, Inc.: Shawnee Mission, KS, 2009. (56) (a) Tanaka, T.; Morino, Y. J. Mol. Spectrosc. 1970, 33, 538−551. (b) Barbe, A.; Secroun, C.; Jouve, P. J. Mol. Spectrosc. 1974, 49, 171− 182. (c) Duncan, J. L.; Wright, I. J.; Van Lerberghe, D. J. Mol. Spectrosc. 1972, 42, 463−477. (57) Martin, J. M. L.; Lee, T. J.; Taylor, P. R.; Francois, J.-P. J. Chem. Phys. 1995, 103, 2589−2602. (58) Gordy, W. C.; Cook, R. L. Microwave Molecular Spectra, Techniques of Chemistry, 3rd ed.; John Wiley and Sons: New York, 1984; Vol. XVIII, p 944. (59) Baker, J. J. Comput. Chem. 1993, 14, 1085−1100. (60) (a) Hamprecht, F. A.; Cohen, A. J.; Tozer, D. J.; Handy, N. C. J. Chem. Phys. 1998, 109, 6264−6271. (b) Boese, A. D.; Doltsinis, N. L.; Handy, N. C.; Sprik, M. J. Chem. Phys. 2000, 112, 1670−1678. (c) Boese, A. D.; Handy, N. C. J. Chem. Phys. 2001, 114, 5497−5503. (61) (a) Cohen, A. J.; Handy, N. C. Mol. Phys. 2001, 99, 607−615. (b) Hoe, W.-M.; Cohen, A. J.; Handy, N. C. Chem. Phys. Lett. 2001, 341, 319−328. (62) Treacy, J.; Hag, M. E.; O’Farrell, D.; Sidebottom, H. Ber. Bunsenges.Phys. Chem. 1992, 96, 422−427. (63) Alam, M. S.; Camredon, M.; Rickard, A. R.; Carr, T.; Wyche, K. P.; Hornsby, K. E.; Monks, P. S.; Bloss, W. J. Phys. Chem. Chem. Phys. 2011, 13, 11002−11015.

10434

dx.doi.org/10.1021/jp307738p | J. Phys. Chem. A 2012, 116, 10420−10434