Direct Dynamics Simulation of the Thermal - ACS Publications

Article pubs.acs.org/JPCA

Cite This: J. Phys. Chem. A 2018, 122, 4808−4818

Direct Dynamics Simulation of the Thermal 3CH2 + 3O2 Reaction. Rate Constant and Product Branching Ratios Sandhiya Lakshmanan,† Subha Pratihar,† Francisco B. C. Machado,‡ and William L. Hase*,† †

Department of Chemistry and Biochemistry, Texas Tech University, Lubbock, Texas 79409, United States Departamento de Química, Instituto Tecnológico de Aeronáutica, São José dos Campos, São Paulo, Brazil

‡

S Supporting Information *

ABSTRACT: The reaction of 3CH2 with 3O2 is of fundamental importance in combustion, and the reaction is complex as a result of multiple extremely exothermic product channels. In the present study, direct dynamics simulations were performed to study the reaction on both the singlet and triplet potential energy surfaces (PESs). The simulations were performed at the UM06/6-311+ +G(d,p) level of theory. Trajectories were calculated at a temperature of 300 K, and all reactive trajectories proceeded through the carbonyl oxide Criegee intermediate, CH2OO, on both the singlet and triplet PESs. The triplet surface leads to only one product channel, H2CO + O(3P), while the singlet surface leads to eight product channels with their relative importance as CO + H2O > CO + OH + H ∼ H2CO + O(1D) > HCO + OH ∼ CO2 + H2 ∼ CO + H2 + O(1D) > CO2 + H + H > HCO + O(1D) + H. The reaction on the singlet PES is barrierless, consistent with experiment, and the total rate constant on the singlet surface is (0.93 ± 0.22) × 10−12 cm3 molecule−1 s−1 in comparison to the recommended experimental rate constant of 3.3 × 10−12 cm3 molecule−1 s−1. The simulation product yields for the singlet PES are compared with experiment, and the most significant differences are for H, CO2, and H2O. The reaction on the triplet surface is also barrierless, inconsistent with experiment. A discussion is given of the need for future calculations to address (1) the barrier on the triplet PES for 3CH2 + 3O2 → 3CH2OO, (2) the temperature dependence of the 3CH2 + 3O2 reaction rate constant and product branching ratios, and (3) the possible non-RRKM dynamics of the 1CH2OO Criegee intermediate.

I. INTRODUCTION Methylene (CH2) is an intermediate formed in the combustion of acetylene, ethylene, and methane, and its chemistry is vital for the growth of higher unsaturated hydrocarbons.1 In lean flames, the fate of triplet methylene (3CH2) is determined predominantly by the attack of molecular oxygen, and in rich flames, singlet methylene (1CH2) chemistry dominates. Since 3CH2 is about 9 kcal/mol lower in energy than 1CH2,2 the reaction between 3CH2 and 3O2 is of fundamental importance in flame chemistry and has been a topic of interest for several experimental3−14 and very few theoretical studies.15,16 The challenging issue for this reaction is that there are at least eight extremely exothermic product channels, as shown below:

The reaction enthalpies given are calculated from the enthalpies of formation of the reactive species available in Active Thermochemical Tables (ATcT).17 The 3CH2 and 3O2 reactants interact on quintet, triplet, and singlet PESs. Due to the repulsive nature of the quintet surface, this surface is less relevant for combustion kinetics. On the singlet and triplet PESs, the 3CH2 + 3O2 reaction proceeds through a Criegee intermediate (CH2OO), which is an important intermediate in the ozonolysis of alkenes.18 Earlier ab initio calculations showed that the carbonyl oxide rearranges into dioxirane and methylenebis(oxy) isomers before dissociation into stable products.15,16,19−21 A CASSCF/CASPT2 investigation of the 3CH2 + 3 O2 reaction obtained an exothermicity of 52 and 32 kcal/mol for the formation of singlet CH2OO and triplet CH2OO in the singlet and triplet PESs, respectively.15 The study also determined the unimolecular dissociation pathways for CH2OO. CH2OO isomerizes to HCOOH which further dissociates into HCO + OH or isomerizes to dioxirane, which again isomerizes to methylene(bis)oxy and further decomposes to CO2 + H2, HCOOH, and HCO + OH. Reaction may occur on both the singlet and triplet PESs. Using transition state theory (TST) and CASSCF/CASPT2 Received: January 29, 2018 Revised: April 23, 2018 Published: April 26, 2018

© 2018 American Chemical Society

4808

DOI: 10.1021/acs.jpca.8b01002 J. Phys. Chem. A 2018, 122, 4808−4818

Article

The Journal of Physical Chemistry A calculations,15 reaction on the triplet PES was estimated to be negligible at 250 K, but at 1800 K, 26% of the reaction was on the triplet surface. Another CASSCF study of the triplet PES and TST calculation of the 3CH2 + 3O2 rate constant16 found that at high temperatures of 1000−1800 K the rate constant has a negative temperature dependence, contrary to the result at low temperatures. Quantum chemical calculations, combined with a kinetic analysis, have been used to model dissociation of the Criegee intermediate.21 Direct chemical dynamics simulations have been used to study unimolecular dynamics associated with the Criegee intermediate on the singlet PES.22,23 In a simulation of O3 + propene ozonolysis, formation of and energy partitioning to the Criegee intermediate were studied.22 Isomerization and decomposition of the Criegee intermediate were studied by initiating trajectories at the CH2OO → CH2O2 transition state (TS).23 Many experimental studies have focused on determining the rate constant and product yields of the 3CH2 + 3O2 reaction. The measured rate constant is (3.2 ± 0.3) × 10−12 cm3 molecule−1 s−1.13,14 Molecular beam experiments of the 3CH2 + 3 O2 reaction, in the temperature range from 295 to 600 K, indicated an activation energy of 1.5 ± 0.3 kcal/mol for association of 3CH2 and 3O2.5 Baulch et al. evaluated kinetic data for the 3CH2 + 3O2 reaction and recommended a rate constant of 3.3 × 10−12 cm3 molecule−1 s−1.24 Further evaluation of the experimental data6 suggested that the overall rate constant is independent of temperature, suggesting that formation of the CH2OO intermediate takes place without an energy barrier. Shock wave experiments at high temperatures, of the 3CH2 + 3 O2 reaction by Dombrowsky and Wagner, found that a quarter of the reaction proceeds by direct formation of radicals and the remaining leads to stable products.9 The major product was CO. Using flash kinetic spectrometry, Alvarez and Moore determined the absolute yields of CO, CO2, H2CO, and OH formed in the reaction as 0.34, 0.40, 0.16, and 0.30, respectively.3 Blitz et al. determined the yield of H atoms as 0.8 ± 0.2.11 The yield of the OH radical was examined by Blitz et al. using pulsed laser photolysis experiments.6 At room temperature, the rate constants for forming OH in υ = 0 and 1 are 2.39 × 10−12 and 2.33 × 10−12 cm3 molecule−1 s−1, respectively. Using LMR and ESR techniques, Bley et al. obtained upper limits for the yields of O, H, and OH as 0.1, 0.2, and 0.3, respectively, for the 3CH2 + 3O2 reaction.8 The rate constant ratio for reactions R1 and R3 was determined as 0.9 by means of time-resolved IR chemiluminescence.10 At the high temperatures of 1850−2050 K, the sum of the H and O atom yields was determined to contribute 80−90% of the reaction.7 Very scattered data indicates a large uncertainty for the 3CH2 + 3 O2 reaction rate constant at high temperatures.24 Despite many kinetic studies of the 3CH2 + 3O2 reaction, the following interesting and challenging questions remain. The 3 CH2 + 3O2 reaction is highly exothermic, and the initial, primary products may contain sufficient energy to undergo secondary decomposition reactions. These decompositions would compete with collisional stabilization of the energized products. To analyze competition between decomposition and stabilization of these products, it is crucial to know their internal energies. By knowing these energies and branching between the different 3 CH2 + 3O2 reaction pathways, it should then be possible to determine relative product yields. Here, results of direct dynamics classical trajectory simulations are presented for the 3CH2 + 3 O2 reaction on both singlet and triplet PESs. The simulations are used to obtain an atomistic understanding of the 3CH2 + 3O2 reaction dynamics and interpret experimental results.

II. COMPUTATIONAL METHODS The computational methodology used for the present study is summarized below. II.A. Electronic Structure Calculations. Properties of stationary points along the singlet and triplet PESs of the 3CH2 + 3 O2 reaction were calculated with the DFT-M06 functional25 and the 6-311++G(d,p) basis set. Previous work26 has indicated that M06 may provide an appropriate description of the multireference nature of the wave function. Minima on the PES were confirmed with all positive frequencies, and each transition state had only one imaginary frequency, confirming its maxima along the reaction coordinate. The unrestricted (UM06) wave function was used for these calculations, and its stability was checked at each stationary point. In order to verify the reliability of the UM06/6-311++G(d,p) level of theory, results for the stationary points are compared with those for previous higher level theoretical calculations.15,16,19 The NWChem computer program27 was used for the electronic structure calculations. II.B. Direct Dynamics Simulations. Direct dynamics simulations for the 3CH2 + 3O2 reaction were performed using UM06/6-311++G(d,p) theory. The VENUS chemical dynamics computer program28,29 interfaced with the NWChem electronic structure program30 was used for the calculations. The simulations were performed for a temperature of 300 K. The 3CH2 + 3 O2 relative translational energy is 0.9 kcal/mol, equivalent to 3RT/2 at 300 K temperature. The vibrational and rotational energies for 3CH2 and 3O2 were sampled from their 300 K Boltzmann distributions. These energies were transformed into the Cartesian coordinates and momenta required for trajectory calculations using quasiclassical sampling.31 Both 3O2 and 3CH2 were randomly rotated with an initial center-of-mass separation of 8 Å. The velocity-Verlet algorithm32 was used to integrate the trajectories with an integration time step of 2 fs. The total integration time for each trajectory was 2 ps. From the trajectory calculations, various dynamical properties were analyzed. The reactive cross section is calculated from the reaction probability using the expression σrxn =

∫0

bmax

Pr(b)2πb db

(1)

where b is the impact parameter at which the reactions are sampled and bmax is the maximum impact parameter at which reaction was observed. Pr(b) is the reaction probability corresponding to each value of b. A total of 100 trajectories were calculated for each b. Using the reactive cross section and the relative velocity (vrel) of the reactants, the reaction rate constant of the reaction (k) was calculated using the expression

k = σrxnvrel

(2)

which is sufficiently accurate, since the reaction cross section is not strongly energy dependent.33,34 Reactive trajectories were categorized by visualizing and determining their atomic level mechanisms. The rate constant and the product branching ratios were compared with available experimental results. The reported uncertainties in the calculated cross sections and product yields are standard deviations.

III. RESULTS OF ELECTRONIC STRUCTURE CALCULATIONS Reaction between 3CH2 and 3O2 initially forms the CH2OO intermediate on both the singlet and triplet surfaces. Contrary to previous calculations,15,16 CH2OO formation occurs without an 4809


Article

The Journal of Physical Chemistry A

Figure 1. Energy profiles of pathways for the 3CH2 + 3O2 reaction on the singlet and triplet PESs, calculated at the UM06/6-311++G(d,p) level of theory.

The pathway for INT2, involving TS11, was found previously by ́ Martinez-Nú ñez.35,36 This pathway leads to cyclic-CO2 along with H2. INT3 forms formic acid (HCOOH) by a highly exothermic reaction with an enthalpy of −182.5 kcal/mol and a small energy barrier of 1.8 kcal/mol. Formic acid decomposes through two pathways, forming the CO + H2O and CO2 + H2 products. The PES for these two formic acid decomposition pathways was previously studied by Goddard et al.37 and Francisco.38 The CO + H2O products are formed through TS5, with elongations of the C−O and C−H bonds and hydrogen bonding between the cleaved hydroxyl group and migrating H atom from the C−H bond (see Figure 2). The UM06/6-311++G(d,p) structure of this transition state is exactly the same as that calculated by Goddard et al.37 The UM06 energy barrier for this reaction is 61 kcal/mol, and the reaction enthalpy is −174.7 kcal/mol; the latter is in good agreement with the value of 177.8 kcal/mol calculated using the Active Thermochemical Tables (ATcT).17 The CO + H2O products lie below the 1CH2OO intermediate with a lower energy of 109.3 kcal/mol, which is comparable with the corresponding MR-PT2/cc-pVDZ value of 104.37 kcal/mol.22 Formation of the CO2 + H2 products involves a tight fourcentered transition state TS6, for which the two H atoms attached to the C and O atoms of formic acid migrate to form H2.

energy barrier. Formation of CH2OO without a barrier is in accord with recommended kinetic data by Baulch et al.,24 suggesting that the overall rate constant is independent of temperature. In the following, the UM06/6-311++G(d,p) singlet and triplet PESs are discussed, as well as possible intersystem crossing (ISC) between these surfaces. The energetics are those for the UM06/6-311++G(d,p) PESs, unless stated otherwise. III.A. Singlet PES. The stationary point potential energy profile for all of the reaction pathways on the singlet PES is shown in Figure 1. Optimized structures of transition states and intermediates on the singlet surface are shown in Figure 2. 1CH2OO unimolecular dissociation occurs via several reaction pathways, of which its fate is mainly determined by the dioxirane (INT2) and formic acid HCOOH (INT4) intermediates.15,19−21 1CH2OO isomerizes to dioxirane through transition state TS2 with an energy barrier of 21.6 kcal/mol which is in good agreement with CASSCF/6-311+G(3df,2p) 15 and high level CCSD(T)/ (4s3p2d1f/3s2p1d) results, i.e., respective energies of 20.1 and 20.3 kcal/mol.19 The UM06/6-311++G(d,p) energy barrier for TS2 calculated in the present study is also in excellent agreement with the MR-PT2/cc-pVDZ energy barrier of 22.78 kcal/mol for isomerization of 1CH2OO into dioxirane computed by Kalinowski et al.22 INT2 further isomerizes into the methylenebis(oxy) isomer biradical intermediate (INT3) by O−O bond breakage. 4810


Article


Figure 2. Optimized structures of the stationary points along the singlet and triplet PESs of the 3CH2 + 3O2 reaction.

transition state TS9 forming the CO + OH + H products. As shown in Figure 2, in transition state TS9 for this reaction, the O−O bond of INT4 is elongated, which weakens the C−H bond and the OH, CO, and H products are then formed exothermically with an enthalpy of −58.4 kcal/mol. This exothermicity is in excellent agreement with the literature value of −58.9 kcal/mol.17 1 CH2OO may climb up a large energy barrier of 99.6 kcal/mol, via a transition state TS10, to form the CO + H2 + O(1D) products. This reaction is exothermic by −34 kcal/mol. Given the large energy barrier, the formation of CO + H2 + O(1D) is relatively unimportant for the 3CH2 + 3O2 reaction. A complete comparison of the calculated reaction and TS energies reported here, and those determined previously with CASSCF/6-311+G(3df,2p) theory, is given in the Supporting Information. III.B. Triplet PES. 3CH2 and 3O2 react on the triplet PES to initially form the 3CH2OO intermediate (INTT) on both singlet and triplet surfaces. The 3CH2OO intermediate is about 30 kcal/mol less stable than 1CH2OO. 3CH2OO decomposes homolytically via TS1 to H2CO and O(3P) in an exothermic reaction with a reaction enthalpy of −62.6 kcal/mol. The energy barrier for this

This TS is moderately similar to the decarboxylation transition state calculated by Goddard37 and Francisco.38 This reaction has an energy barrier of 67.6 kcal/mol and is extremely exothermic with a ΔH298 value of −193.7 kcal/mol. On comparing the barrier and the reaction enthalpy of the two dissociation pathways for formic acid, it is obvious that both pathways are feasible. The dioxirane intermediate (INT2) decomposes by O−O and C−O bond cleavages simultaneously forming formaldehyde and O(1D) with a reaction enthalpy of −37.1 kcal/mol. The transition state, TS7, through which this reaction occurs, is a late TS in which the C−O bond is already broken and the resulting CH2OO structure is linear. This linear structure decomposes into H2CO and O(1D). Aside from isomerising into dioxirane, 1 CH2OO rearranges into a relatively unstable intermediate HCOOH (INT4) through transition state TS8, with an energy barrier which is around 13 kcal/mol less than that required for isomerizing into dioxirane through TS2. INT4 dissociates by O−O bond breakage and forms CHO and OH in a nearly barrierless reaction, in good agreement with the result of Gutbrod et al.39 Furthermore, INT4 crosses a small energy barrier of 5.4 kcal/mol through 4811


Article

The Journal of Physical Chemistry A reaction is 4.7 kcal/mol, in close agreement with the CASSCF activation barrier of 4.1 kcal/mol.15 Unlike transition state TS7 on the singlet surface, which has a linear C−O−O structure and leads to H2CO and O(1D), for TS1, the C−O−O group is nonlinear. The reaction energy for 3CH2 + 3O2 → H2CO + O(3P) was calculated with different levels of theory and basis sets, and the results are given in the Supporting Information. The UM06 theory used here gives a reaction energy in excellent agreement with that of CCSD(T) and the complete basis set (CBS) limit. III.C. Intersystem Crossing. Though intersystem crossing (ISC) between the triplet and singlet PESs of the 3CH2 + 3O2 reaction is possible, Anglada et al. found through spin−orbit coupling calculations that ISC is unimportant.40 This is consistent with direct reaction on the triplet PES, not allowing long times for coupling of the triplet surface with the singlet surface. In addition, the Criegee intermediate is substantially higher in energy on the triplet than the singlet PES, i.e., 31.7 kcal/mol. As a further test of ISC, energies for all of the INT and TS stationary point structures on the singlet PES were calculated for the triplet PES and found to be substantially higher than for the singlet PES. This finding indicates that the triplet and singlet surfaces do not cross. Crossings would facilitate ISC.41,42

IV. RESULTS OF DIRECT DYNAMICS SIMULATIONS As discussed above, structures of stationary points, barrier heights, and thermochemistry for the 3CH2 + 3O2 reaction, calculated with UM06/6-311++G(d,p), are in close agreement with earlier high level ab initio calculations and experimental results. Given this agreement and that UM06/6-311++G(d,p) is computationally feasible, it was used for the direct dynamics simulations of the 3 CH2 + 3O2 reaction. IV.A. Singlet Surface. 1. Reaction Pathways, Probabilities, and Atomistic Dynamics. Reaction dynamics on the singlet surface 3CH2 + 3O2 reaction were studied versus collision impact parameter b at a temperature of 300 K, for which the collision energy Erel is 0.9 kcal/mol. The simulations were performed for b = 0, 1, 2, 3, 4, 5, and 6 Å. At the largest b of 6 Å, no reactions were observed, identifying bmax. The following pathways were identified in the simulations:

Figure 3. Reaction probability as a function of the impact parameter (b) for the 3CH2 + 3O2 reaction on the singlet surface at 300 K, with 0.9 kcal/mol relative collision energy: (a) total reaction probability; (b) reaction probability of individual product channels. 1

CH2OO complex. In a few cases, 1CH2OO is followed by formation of the dioxirane intermediate before dissociation into final products. Snapshots from reactive trajectories representative of the above eight mechanisms are illustrated in Figure 4. For each trajectory, the 1CH2OO complex further dissociates into the final stable products. Trajectory 1: In formation of the CO + H2O products (R2), dioxirane is formed immediately after 1CH2OO and there is simultaneous stretching of the C−H as well as the O−O bonds of dioxirane (1.05 ps). The CO + H2O products are then formed with O−O and C−H bond ruptures. Trajectory 2: The CO2 + H2 products (R1) are formed from the dioxirane intermediate (0.79 ps). This intermediate undergoes O−O stretching and subsequent bond breakage, resulting in the methylenebis(oxy) intermediate, followed by C−H bond stretching and rupture, and then formation of the H−H bond and the CO2 and H2 products. Trajectory 3: The CO + OH + H products (R8) are formed from HCOOH. 1CH2OO undergoes a rearrangement where one of the H atoms migrates and binds with the terminal O atom, thereby forming HCOOH (1.12 ps). The unstable HCOOH

The total reaction probability and those for the above pathways, versus b, are shown in Figure 3. The total reaction probability is around 0.68−0.55 for b = 0−3 Å and then rapidly drops with increasing b. The probability for formation of CO2 + H2O (R2) product is high for b = 0−4 Å, and the least probable channel is HCO + O + H (R11). The maximum probability for the extremely exothermic reaction R1 is at b = 0, and the probability decreases substantially for larger b. At b = 5 Å, only reactions R5 and R10 are observed. On examining the trajectories, it was seen that all of the trajectories proceeded through the 4812


Article


Figure 4. Representative snapshots of trajectories for the eight product channels for the 3CH2 + 3O2 reaction on the singlet surface.

(1.09 ps). Trajectory 5: Formation of the H2CO + O(1D) products (R10) involves rearrangement of 1CH2OO to dioxirane, followed by C−O and O−O bond stretching before forming the final products. Trajectory 5 events are short-lived with all reactions occurring within 1 ps. Trajectory 6: The HCO + OH products (R5) are formed by unimolecular rearrangement of

complex then undergoes concomitant O−O and C−H bond stretching and after a few oscillations forms the CO + OH + H products. Trajectory 4: For formation of the CO + H2 + O(1D) products (R9), 1CH2OO immediately undergoes C−H and O−O bond stretching, thereby breaking the O−O bond and forming a H−H bond. The products are formed nearly instantaneously 4813


Article

The Journal of Physical Chemistry A 1

The second most contributing product channel is CO + OH + H with ⟨Pr⟩ = 0.07, closely followed by CO2 + H2 and CH2O + O(1D) with respective ⟨Pr⟩ of 0.05 in both cases. The lowest contributions to reaction are from HCO + OH, CO + H2 + O(1D), CO2 + H + H, and HCO + O(1D) + H with ⟨Pr⟩ of 0.03, 0.02, 0.006, and 0.002, respectively. The HCO + OH contribution is small because the weakly bound HCO primarily dissociates into CO + H. The small contribution of CO + H2 + O(1D) may be attributed to the large energy barrier to form CO + H2 + O(1D) from 1CH2OO, as shown in Figure 1. The CO2 + H + H product channel was observed only at small impact parameters. It only occurs if H2 in the CO2 + H2 product channel dissociates into two H-atoms, which is unlikely. The least important product channel is HCO + O(1D) + H, which results if OH for the HCO + OH channel dissociates to O(1D) + H. 3. Singlet Product Yields. Following the procedure used to analyze experimental product yields for the 3CH2 + 3O2 reaction,3 reactive cross sections for the current simulation may be used to determine product yields for reaction on the singlet PES. These results may be compared with experiment if the singlet PES dominates the 3CH2 + 3O2 reaction and reaction on the triplet PES is negligible. Two different analyses were performed, i.e., one which does not include the nonreactive trajectories which remain trapped on the PES and another which accounts for these trajectories. The former is considered first, and to illustrate the calculation of a singlet PES product yield, the H2CO product yield is the cross section for R10 on the singlet surface divided by the total singlet cross section, i.e., 4.6/25.5 = 0.18. Using this approach, the yields for products O(1D), H, and OH are 0.26 ± 0.07, 0.22 ± 0.06, and 0.30 ± 0.08, respectively. These may be compared to the respective upper limits of 0.07, 0.10, and 0.10 determined by laser magnetic resonance.7 The CO, CO2, and H2CO yields calculated from the present study are 0.60 ± 0.20, 0.10 ± 0.03, and 0.18 ± 0.04, respectively, as compared to those of 0.34 ± 0.06, 0.40 ± 0.08, and 0.16 ± 0.04 obtained by Alvarez and Moore3 using pulsed laser photolysis experiments. Their work provided indirect evidence for OH yield of 0.30 ± 0.08, in good agreement with the simulation results. From this analysis of the simulations, the respective yields for the remaining products H2O, H2, and HCO are 0.32 ± 0.08, 0.17 ± 0.04, and 0.10 ± 0.02. The product yields obtained from the simulations are summarized in Table 2. The discrepancy between the simulations and experiments for product yields may result in part from the trajectories which do not react and remain as either 1CH2OO or dioxirane when the trajectories are terminated. As shown in Table 1, the cross section for these trajectories is 5.4 Å2, for which 92% is for 1CH2OO and 8% is for dioxirane. Though an approximation, if these trajectories are assumed to form products in accord with the product branching ratios found here, the cross sections for forming products increase and are given in Table 1. However, as shown in Table 2, there are not any significant changes in product yields. IV.B. Triplet Surface. The 3CH2 + 3O2 reaction on the triplet surface has only one product channel, H2CO + O (3P). Trajectories on the triplet surface were calculated for b of 0, 1, 2, 3, 3.5, 4, 4.5, and 5 Å. At b = 5 Å, no reactions were observed. Figure 5 shows the reaction probability versus b. The probability is a maximum at b of 0−2 Å and decreases with increasing b. All trajectories passing through the 3CH2OO potential minimum resulted in H2CO + O(3P) formation. Figure 6 shows snapshots of a representative reaction event on the triplet surface. All reactive trajectories proceeded through 3 CH2OO (0.63 ps). This was followed by O−O bond stretching,

CH2OO into HCOOH, followed by O−O bond homolysis into HCO + OH. Trajectory 7: Unlike the CO2 + H2O product channel (R2), formation of the CO2 + H + H products (R4) only involves stretching of the 1CH2OO O−O bond, which oscillates for some time. Then, the terminal O atom binds with the C atom, simultaneously breaking the C−H bonds. The two H atoms attempt to form a bond but eventually separate and form the CO2 + H + H products. Trajectory 8: The CHO + O(1D) + H (R11) products are formed from the dissociation of 1CH2OO, where the stretching of the O−O bond is associated with stretching of one C−H bond, leading to the breaking of O−O and one C−H bond forming CHO + O(1D) + H product channels. From the above descriptions, only the CO2 + H2 (R1), CO + H2O (R2), and H2CO + O(1D) (R10) product channels involve the dioxirane intermediate and all other product channels follow 1 CH2OO dissociation. Also, the CO + H2O products (R2) contribute significantly for reaction on the 3CH2 + 3O2 singlet surface. 2. Reaction Cross Sections and Rate Constant. Using the reaction probabilities in Figure 3 and eq 1, the total cross section and cross sections for individual product channels were calculated and are summarized in Table 1. The CO + H2O product Table 1. Calculated 3CH2 + 3O2 Reactive Cross Sections σr on the Singlet PES product channel

σr (Å2)

CO + H2O (R2) CO2 + H2 (R1) CO + OH + H (R8) CO + H2 + O(1D) (R9) H2CO + O(1D) (R10) HCO + OH (R5) HCO + O(1D) + H (R11) CO2 + H + H (R4) trapped trajs.a total cross section with trapped trajs. total cross section without trapped trajs.

8.2 ± 2.1 (9.9)b 2.4 ± 0.6 (2.9) 5.3 ± 1.4 (6.4) 2.0 ± 0.6 (2.4) 4.6 ± 1.1 (5.6) 2.6 ± 0.6 (3.2) 0.1 ± 0.08 (0.12) 0.3 ± 0.1 (0.4) 5.4 (0.0) 30.9 ± 7.5 25.5 ± 6.2

a

Trajectories that are trapped as CH2OO or dioxirane when the trajectories are halted. bCross sections if the trapped trajectories are assumed to form products in accord with the product branching ratios found here.

channel has the largest cross section. The overall rate constant was calculated at 300 K from eq 2 with a correction for spin degeneracies. The collisions are on the quintet, triplet, and singlet surfaces, and with these dynamics included, 1/9 of the collisions are on the singlet surface and the total cross section should be multiplied by this factor.42,43 The resulting rate constant at 300 K is (0.76 ± 0.18) × 10−12 cm3 molecule−1 s−1 and smaller than the recommended value of 3.3 × 10−12 cm3 molecule−1 s−1. A number of trajectories remained as 1CH2OO or dioxirane in Figure 1, when the trajectories were terminated, and did not undergo further dissociation into stable products. The cross section for these trajectories is given in Table 1, and from the cross sections, these trajectories constituted the fraction [5.4/(5.4 + 25.5)] = 0.17 of the trajectory events on the singlet surface. If the cross section for these trapped trajectories is included, the calculated rate constant increases by 21% and becomes (0.93 ± 0.22) × 10−12 cm3 molecule−1 s−1. Average reaction probabilities, ⟨Pr⟩, for the individual product channels may be determined by dividing their reaction cross sections in Table 1 by πbmax2. bmax = 5 Å was used for this analysis, and these ⟨Pr⟩ do not include the 1/9 spin degeneracy factor. The CO + H2O product channel is most significant with ⟨Pr⟩ = 0.13. 4814


Article


in substantial disagreement with experiment. In the next section, the difficulty in establishing the relative importance of reaction on the singlet and triplet PESs is discussed.

Table 2. Calculated Product Yields yield product

reactive trajectories

reactive + trapped trajectories

1

0.26 ± 0.07 0.22 ± 0.06 0.30 ± 0.08 0.60 ± 0.20 0.10 ± 0.03 0.18 ± 0.04 0.32 ± 0.08 0.17 ± 0.04 0.10 ± 0.02

0.26 0.22 0.31 0.60 0.11 0.18 0.32 0.17 0.11

O( D) H OH CO CO2 H2CO H2O H2 HCO

V. RELATIVE IMPORTANCE OF REACTION ON THE SINGLET AND TRIPLET PESS Reaction on the triplet PES to form 3CH2OO is barrierless as for reaction on the singlet PES to form 1CH2OO. With spin degeneracies included, the direct dynamics rate constant on the singlet PES is (0.93 ± 0.22) × 10−12 and (1.8 ± 0.3) × 10−12 cm3 molecule−1 s−1 on the triplet PES. Their sum is 2.7 × 10−12 cm3 molecule−1 s−1, which is nearly the same as the temperature independent recommended experimental value of 3.3 × 10−12 cm3 molecule−1 s−1.24 With these rate constants, reaction on the triplet surface will dominate, and as discussed below, the H2CO yield becomes 0.71, much different than experiment.3,6 CASSCF and CASPT2 calculations by Chen et al.15 find an enthalpy of activation for reaction on the triplet surface 3.5 kcal/mol higher than that for the singlet surface. If this finding is correct, reaction on the triplet surface will be unimportant at 300 K and the reaction rate constant given by reaction on the singlet surface. As discussed above in section IV.A.3, the simulation product yields for reaction on the singlet PES are in qualitative agreement with experiment. In contrast, if the simulation results for reaction on the triplet PES are included, there are substantial differences between the simulation and experimental product yields, since H2CO becomes the dominant simulation product. With reaction on both PESs included, the product yields given in parentheses are O(1D) (0.09), H (0.08), OH (0.11), CO (0.21), CO2 (0.04), H2CO (0.71), H2O (0.11), H2 (0.06), HCO (0.04), and O(3P) (0.65). These yields are substantially different from the experimental yields and those for the singlet PES simulation, discussed in section IV.A.3. If there was efficient intersystem crossing (ISC) from the triplet to singlet PES, with nearly all of the trajectories on the triplet PES undergoing ISC, product formation may be assumed to occur on the singlet PES. However, the dynamics of these ISC trajectories are expected to be different from those for the 3CH2 + 3 O2 → 1CH2OO collision trajectories, since the ISC and collision trajectories are expected to sample different regions of the singlet PES. Though this is an interesting conjecture, ISC is not expected to be important for the 3CH2 + 3O2 reaction dynamics, since as discussed in section III.C triplet and singlet potential energy curves do not cross. Even if ISC did occur, only a fraction of the trajectories on the triplet PES would follow these dynamics.41,42

Figure 5. Total reaction probability as a function of the impact parameter (b) for the 3CH2 + 3O2 reaction on the triplet surface at 300 K, with 0.9 kcal/mol relative collision energy.

immediately leading to O−O bond breakage and formation of H2CO + O(3P). Using the reaction probabilities versus b, the cross section for this reaction is 19.2 ± 4.6 Å2. With collisions on the quintet, triplet, and singlet surfaces and the 3/9 spin degeneracy factor included,42,43 the overall rate constant for reaction on the triplet surface is (1.8 ± 0.3) × 10−12 cm3 molecule−1 s−1. Using bmax = 4.5 Å, the average reaction probability is then 0.12. Above, product yields were calculated assuming reaction occurred on the singlet PES, with negligible reaction on the triplet PES. However, the calculated cross section for reaction on the triplet PES is similar to that for reaction on the singlet PES, and if spin degeneracy factors were included, reaction on the triplet surface would dominate. In comparison to the singlet surface, 3/4 of the collisions are on the triplet surface.42,43 Such a model would predict high yields for H2CO and O(3P) formation,

VI. CONCLUSIONS AND FUTURE DIRECTIONS The electronic structure calculations and direct dynamics simulations presented here provide detailed atomic-level mechanisms

Figure 6. Representative snapshots of trajectories for the product channel for the 3CH2 + 3O2 reaction on the triplet surface. 4815


Article

The Journal of Physical Chemistry A Table 3. Experimental Measurements of 3CH2 + 3O2 Product Yields reference product

8

1

0.10 0.20 0.30

O( D) H OH CO CO2 H2CO H2O H2 HCO

3

11 0.8 ± 0.2

0.30 ± 0.05 0.34 ± 0.06 0.40 ± 0.08 0.16 ± 0.04

of the 3CH2 + 3O2 reaction on the singlet and triplet PESs. The simulations were performed for a temperature of 300 K using UM06/6-311++G(d,p) theory, which gives reaction energies and barrier heights in overall good agreement with previous higher-level ab initio calculations.15,16,19−21 As shown in Figure 1, reaction products are formed via highly exothermic pathways, with multiple potential energy minima and barriers. The barriers are substantially lower than the 3CH2 + 3O2 reactant energy. The 3 CH2 + 3O2 interaction to form CH2OO is purely attractive and barrierless, seemingly in agreement with the recommendation that the overall reaction rate is independent of temperature.24 However, the transition state to form CH2OO is variational44 and the resulting rate constant may have either a small negative or positive temperature dependence, dependent on properties of the attractive PES.45 The direct dynamics simulations show that all of the trajectories initially proceed through the CH2OO potential energy minimum on both the singlet and triplet PES. The triplet surface leads to only the H2CO + O(3P) products. In contrast, the singlet surface leads to the eight product channels listed in Table 1 and illustrated in Figure 1. The 1CH2OO intermediate may isomerize to the dioxirane intermediate and further to the methylene(bis) oxy isomer, leading to formic acid, which dissociates to CO + H2O and CO2 + H2. The dioxirane intermediate further dissociates into H2CO + O(1D). Only these three product channels involve the dioxirane intermediate, and the remaining five follow directly from 1CH2OO. From the simulations at 300 K, the dominant product channel on the singlet PES is CO + H2O, followed by CO + OH + H and H2CO + O(1D). The two least important product channels are CO2 + H + H and HCO + O(1D) + H, with the latter the lesser of the two. At room temperature, reaction is thought to occur on the singlet PES and the singlet product yields determined from the simulations may be compared with the experimental product yields1,3,6,8,11,46 summarized in Table 3. The two most extensive studies are those by Alvarez and Moore3 and Bley et al.8 In addition, using a composite of experimental results, Blitz et al.6 suggested product yields. The most significant disagreements between the product yields of the experimental studies and the simulation results are for H, CO2, and H2O. From the Blitz et al. study and two other studies, the H atom yield is ∼0.80, while the simulation result is 0.22 ± 0.06. However, the studies by Bley et al.8 and Lindstedt and Skevis1 give a yield close to the simulation result. The CO2 yield from the simulations is much lower than that from the experiments, i.e., 0.10 ± 0.03 versus ∼0.40− 0.45. The H2O yield from the simulations is 0.32 ± 0.08, while it is 0.06 from the experiments of Blitz et al.6 In comparing the simulation and experimental yields, it is noteworthy to point out that Blitz et al. assumed the yields for the CO + H2 + O and

6

46

1

0.18 ± 0.04 0.78 0.43 0.38 0.44 0.18 0.06 0.21 0

0.26 0.78

0.10 0.25

0.37 0.37 0.26 0.11 0.11

0.72 0.18 0.10 0.08

HCO + OH channels are zero, while from the simulation cross sections on the singlet PES the yields for the product channels are 0.08 and 0.10. An interesting and important dynamical question is to what extent dissociation of the CH2OO Criegee intermediate is statistical on the singlet PES. 3CH2 + 3O2 association to form 1CH2OO is expected to deposit energy nonstatistically within the intermediate.23 For statistical unimolecular dissociation, intramolecular vibrational energy redistribution (IVR)47 within 1CH2OO must be faster than its dissociation dynamics. Tests for such dynamics may be probed by future direct dynamics simulations, in which 1 CH2OO is excited randomly with a microcanonical ensemble of states. Branching between the product channels for this statistical excitation may be compared with that for the 3CH2 + 3O2 → 1 CH2OO nonstatistical excitation. In addition to IVR, intramolecular dynamics associated with rotational excitation of 1CH2OO and K-mixing48−50 may also be important. Figure 1 shows that branching between the different product channels on the singlet PES depends on the collision orbital angular momentum, l = μbvrel. The orbital angular momentum is not expected to add statistically to the 1CH2OO rotation axes, with the 1CH2OO total angular momentum J the vector sum J = jCH2 + jO2 + l, where jCH2 and jO2 are the 3CH2 and 3O2 rotation angular momentum vectors. The component of J about the 1CH2OO symmetry axis gives rise to the K rotation quantum number, which is not expected to be statistically populated by 3CH2 + 3O2 → 1CH2OO association. However, at longer times, as a result of vibration− rotation coupling,48,49 a statistical population of the K-levels may occur. In future direct dynamics simulations, it would be of interest to identify the origin of the orbital angular dependence for the product branching ratios and the time-scale for statistical K-mixing.50 High-level ab initio calculations and additional DFT calculations to investigate the 3CH2 + 3O2 → CH2OO potential energy curves on the singlet and triplet PESs would be useful. UM06/6-311++G(d,p) theory used here gives a barrierless 3CH2 + 3 O2 → CH2OO reaction on the singlet PES, consistent with the suggested kinetics,24 but also a barrierless reaction on the triplet surface. In previous CASSCF/CASPT2 calculations,15 activation enthalpies of 1.9 and 5.4 kcal/mol at 298 K were found for 3CH2 + 3 O2 → CH2OO reaction on the singlet and triplet PESs, respectively. An activation barrier for the singlet PES seems incorrect, but an activation barrier for the triplet PES is consistent with experiment. It will be of interest to see if the highest level electronic structure calculations give a barrierless reaction on the singlet surface but a barrier for the triplet surface. Simulations of the 3CH2 + 3O2 reaction versus temperature will be important. For example, the overall rate constant may be 4816


Article


(4) Martinez, R. I.; Huie, R. E.; Herron, J. T. The Role of the Criegee Intermediate in the Matrix Thermoluminescence Study of the CH2 + O2 Reaction. J. Chem. Phys. 1981, 75, 5975−5977. (5) Vinckier, C.; Debruyn, W. Temperature Dependence of the Reactions of Methylene with Oxygen Atoms, Oxygen, and Nitric Oxide. J. Phys. Chem. 1979, 83, 2057−2062. (6) Blitz, M. A.; Kappler, C.; Pilling, M. J.; Seakins, P. W. 3CH2 + O2: Kinetics and Product Channel Branching Ratios. Z. Phys. Chem. 2011, 225, 957−967. (7) Lee, P.-F.; Matsui, H.; Chen, W.-Y.; Wang, N.-S. Production of H and O(3P) Atoms in the Reaction of CH2 with O2. J. Phys. Chem. A 2012, 116, 9245−9254. (8) Bley, U.; Temps, F.; Wagner, H. Gg.; Wolf, M. Investigations of the Reaction Between CH2(X̃ 3B1) and O2 in the Temperature Range 233 K ≤ T ≤ 433 K. Ber. Bunsen-Ges. Phys. Chem. 1992, 96, 1043−1048. (9) Dombrowsky, Ch.; Wagner, H. Gg Investigation of the 3CH2 + O2 Reaction in Shock Waves. Ber. Bunsen-Ges. Phys. Chem. 1992, 96, 1048− 1055. (10) Hancock, G.; Haverd, V. A Time-resolved FTIR Emission Study of the Gas Phase Removal Processes of CH2(X̃ 3B1) and CH2(ã1A1) in Collisions with O2. Chem. Phys. Lett. 2003, 372, 288−294. (11) Blitz, M. A.; McKee, K. W.; Pilling, M. J.; Seakins, P. W. Evidence for the Dominance of Collision-induced Intersystem Crossing in Collisions of 1CH2 with O2 and a Determination of the H Atom Yields from 3CH2+O2, Using Time-resolved Detection of H Formation by vuvLIF. Chem. Phys. Lett. 2003, 372, 295−299. (12) Dombrowsky, Ch.; Hwang, S. M.; Röhrig, M.; Wagner, H. Gg The Formation of O and H Atoms in the Reaction of CH2 with O2 at High Temperatures. Ber. Bunsen-Ges. Phys. Chem. 1992, 96, 194−198. (13) Darwin, D. C.; Young, A. T.; Johnston, H. S.; Moore, C. B. Rate Constants for Triplet Methylene (X̃ 3B1) Removal by Oxygen, Nitric oxide and Acetylene from Infrared Diode Laser Flash Kinetic Spectroscopy. J. Phys. Chem. 1989, 93, 1074−1078. (14) Böhland, T.; Temps, F.; Wagner, H. Gg. Direct Determination of the Rate Constant for the Reaction CH2 + O2 with a LMRSpectrometer. Ber. Bunsen-Ges. Phys. Chem. 1984, 88, 455−458. (15) Chen, B. − Z.; Anglada, J. M.; Huang, M. − B.; Kong, F. The Reaction of CH2 (X3B1) with O2 (X3Σg−): A Theoretical CASSCF/ CASPT2 Investigation. J. Phys. Chem. A 2002, 106, 1877−1884. (16) Fang, D. − C.; Fu, X. − Y. CASSCF and CAS+1 + 2 Studies on the Potential Energy Surface and the Rate Constants for the Reactions between CH2 and O2. J. Phys. Chem. A 2002, 106, 2988−2993. (17) Ruscic, B. Active Thermochemical Tables. Updated Active Thermochemical Tables (ATcT) Values Based on Version 1.122 of the Thermochemical Network, 2013, https://atct.anl.gov/ Thermochemical%20Data/version%201.122/index.php. (18) Vereecken, L.; Glowacki, D. R.; Pilling, M. J. Theoretical Chemical Kinetics in Tropospheric Chemistry: Methodologies and Applications. Chem. Rev. 2015, 115, 4063−4114. (19) Cremer, D.; Kraka, E.; Szalay, P. G. Decomposition Modes of Dioxirane, Methyldioxirane, and Dimethyldioxirane - A CCSD(T), MRAQCC, and DFT Investigation. Chem. Phys. Lett. 1998, 292, 97−109. (20) Li, J.; Guo, H. Full-Dimensional Potential Energy Surface and Rovibrational Levels of Dioxirane. J. Phys. Chem. A 2016, 120, 2991−2998. (21) Nguyen, T. L.; Lee, H.; Matthews, D. A.; McCarthy, M. C.; Stanton, J. F. Stabilization of the Simplest Criegee Intermediate from the Reaction between Ozone and Ethylene: A High-Level Quantum Chemical and Kinetic Analysis of Ozonolysis. J. Phys. Chem. A 2015, 119, 5524−5533. (22) Kalinowski, J.; Räsänen, M.; Heinonen, P.; Kilpeläinen, I.; Gerber, R. B. Isomerization and Decomposition of a Criegee Intermediate in the Ozonolysis of Alkenes: Dynamics Using a Multireference Potential. Angew. Chem., Int. Ed. 2014, 53, 265−268. (23) Vayner, G.; Addepalli, V.; Song, K.; Hase, W. L. Post-transition State Dynamics for Propene Ozonolysis: Intramolecular and Unimolecular Dynamics of Molozonide. J. Chem. Phys. 2006, 125, 014317-1− 014317-16. (24) 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.; Walker, R. W.;

temperature independent, but the product branching may be temperature dependent. Finding information with respect to the latter would be quite helpful in interpreting experiments. Due to the exothermic nature of the product channels, many products may have sufficient energy to decompose further. For instance, the HCO radicals are expected to dissociate at high temperatures, thereby initializing radical chain reactions. Thus, the products formed at high temperatures are uncertain. Even though a large amount of kinetic data at high temperatures exists for the 3CH2 + 3 O2 reaction,24 the data is very scattered, suggesting large uncertainty limits for the rate constant. High temperature simulations are needed. In future work, the product translational, rotational, and vibrational energies for the different reaction channels will be determined. With this information, rate constants for forming OH in the n = 0 and n = 1 vibrational states may be determined and compared with experiment.6 Finally, it will be important to perform the 3CH2 + 3O2 direct dynamics simulations with additional electronic structure theories.

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpca.8b01002. Enthalpies of the studied reaction pathways calculated in the present study and in earlier studies and activation enthalpies (kcal/mol) for 3CH2 + 3O2 reaction pathways calculated in the present and earlier study (PDF)

■

AUTHOR INFORMATION

ORCID

Sandhiya Lakshmanan: 0000-0001-9055-096X Francisco B. C. Machado: 0000-0002-2064-3463 William L. Hase: 0000-0002-0560-5100 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The research reported here is based upon work supported by the Air Force Office of Scientific Research (AFOSR) under Grant No. FA9550-17-1-0119 and the Robert A. Welch Foundation under Grant No. D-0005. F.B.C.M. thanks Fundacão de Amparo à Pesquisa do Estado de São Paulo (FAPESP) under grant 2017/ ́ 07707-3 and Conselho Nacional de Desenvolvimento Cientificio e Technológico (CNPq) under grants 307052/2016-8 and 404337/2016-3. The simulations were performed on the Robinson computer cluster of the Department of Chemistry and Biochemistry at Texas Tech University, whose purchase was funded by the National Science Foundation under the CRIF-MU Grant CHE-0840493. Hai Wang and Greg Smith are thanked for important discussions.

■

REFERENCES

(1) Lindstedt, R. P.; Skevis, G. Chemistry of Acetylene Flames. Combust. Sci. Technol. 1997, 125, 73−137. (2) Hase, W. L.; Phillips, R. J.; Simons, J. W. Vibrational Deactivation of Singlet Methylene. Chem. Phys. Lett. 1971, 12, 161−165. (3) Alvarez, R. A.; Moore, C. B. Absolute Yields of CO, CO2, and H2CO from the Reaction CH2(X 3B1) + O2 by IR Diode Laser Flash Kinetic Spectroscopy. J. Phys. Chem. 1994, 98, 174−183. 4817


Article


Experimental Association Rate Constants. J. Chem. Phys. 1991, 95, 8073−8081. (44) Hase, W. L. Variational Unimolecular Rate Theory. Acc. Chem. Res. 1983, 16, 258−264. (45) Hu, X.; Hase, W. L. Properties of Canonical Variational Transition State Theory for Association Reactions without Potential Energy Barriers. J. Phys. Chem. 1989, 93, 6029−6038. (46) Leung, K. M.; Lindstedt, R. P. Detailed Kinetic Modeling of C1 − C3 Alkane Diffusion Flames. Combust. Flame 1995, 102, 129−160. (47) Hase, W. L. Unimolecular and Intramolecular Dynamics: Relationship to Potential Energy Surface Properties. J. Phys. Chem. 1986, 90, 365−374. (48) Zhu, L.; Hase, W. L. Comparison of Models for Calculating the RRKM Unimolecular Rate Constant k(E,J). Chem. Phys. Lett. 1990, 175, 117−124. (49) Zhu, L.; Chen, W.; Hase, W. L.; Kaiser, E. W. Comparison of Models for Treating Angular Momentum in RRKM Calculations with Vibrator Transition States. Pressure and Temperature Dependence of Cl + C2H2 Association. J. Phys. Chem. 1993, 97, 311−322. (50) Ma, X.; Tan, X.; Hase, W. L. Effects of Vibrational and Rotational Energies on the Lifetime of the Pre-reaction Complex for the F− + CH3I S N 2 Reaction. Int. J. Mass Spectrom. 2017, DOI: 10.1016/ j.ijms.2017.07.011.

Warnatz, J. Evaluated Kinetic Data for Combustion Modeling: Supplement II. J. Phys. Chem. Ref. Data 2005, 34, 757. (25) Zhao, Y.; Truhlar, D. G. The M06 Suite of Density Functionals for Main Group Thermochemistry, Thermochemical Kinetics, Noncovalent Interactions, Excited States, and Transition Elements: Two New Functionals and Systematic Testing of Four M06-class Functionals and 12 other Functionals. Theor. Chem. Acc. 2008, 120, 215−241. (26) Alecu, I. M.; Truhlar, D. G. Computational Study of the Reactions of Methanol with the Hydroperoxyl and Methyl radicals. 1. Accurate Thermochemistry and Barrier Heights. J. Phys. Chem. A 2011, 115, 2811−2829. (27) Valiev, M.; Bylaska, E. J.; Govind, N.; Kowalski, K.; Straatsma, T. P.; van Dam, H. J. J.; Wang, D.; Nieplocha, J.; Apra, E.; Windus, T. L.; de Jong, W. A. NWChem: a Comprehensive and Scalable Open-Source Solution for Large Scale Molecular Simulations. Comput. Phys. Commun. 2010, 181, 1477−1489. (28) Hase, W. L.; Duchovic, R. J.; Hu, X.; Komornicki, A.; Lim, K. F.; Lu, D. H.; Peslherbe, G. H.; Swamy, S. R.; Vande Linde, S. R.; Varandas, A.; Wang, H.; Wolf, R. J. Quantum Chem. Program Exch. (QCPE). Bull. 1996, 16, 43. https://opus.ipfw.edu/chemistry_facpubs/81. (29) Hu, X.; Hase, W. L.; Pirraglia, T. Vectorization of the General Monte Carlo Classical Trajectory Program VENUS. J. Comput. Chem. 1991, 12, 1014−1024. (30) Lourderaj, U.; Sun, R.; Kohale, S. C.; Barnes, G. L.; de Jong, W. A.; Windus, T. L.; Hase, W. L. The VENUS/NWChem Software Package. Tight Coupling between Chemical Dynamics Simulations and Electronic Structure Theory. Comput. Phys. Commun. 2014, 185, 1074−1080. (31) Peslherbe, G. H.; Wang, H.; Hase, W. L. Monte Carlo Sampling for Classical Trajectory Simulations. Adv. Chem. Phys. 2007, 105, 171− 201. (32) Schlick, T. Molecular Modelling and Simulation; Springer: New York, 2000; pp 400−401. (33) Swamy, K. N.; Hase, W. L. Dynamics of Ion−Molecule Recombination. III. Trends in the Recombination Efficiency. J. Am. Chem. Soc. 1984, 106, 4071−4077. (34) Paul, A. K.; Kolakkandy, S.; Hase, W. L. Dynamics of Na+(Benzene) + Benzene Association and Ensuing Na+(Benzene)2*Dissociation. J. Phys. Chem. A 2015, 119, 7894−7904. (35) Martínez-Núñez, E. An Automated Method to Find Transition States Using Chemical Dynamics Simulations. J. Comput. Chem. 2015, 36, 222−234. (36) Martínez-Núñez, E. An Automated Transition State Search Using Classical Trajectories Initialized at Multiple Minima. Phys. Chem. Chem. Phys. 2015, 17, 14912−14921. (37) Goddard, J. D.; Yamaguchi, Y.; Schaefer, H. F. The Decarboxylation and Dehydration Reactions of Monomeric Formic Acid. J. Chem. Phys. 1992, 96, 1158−1166. (38) Francisco, J. S. A. Comprehensive Theoretical Examination of Primary Dissociation Pathways of Formic Acid. J. Chem. Phys. 1992, 96, 1167−1175. (39) Gutbrod, R.; Schindler, R. N.; Kraka, E.; Cremer, D. Formation of OH Radicals in the Gas Phase Ozonolysis of Alkenes: The Unexpected Role of Carbonyl Oxides. Chem. Phys. Lett. 1996, 252, 221−229. (40) Anglada, J. M.; Bofill, J. M.; Olivella, S.; Sole, A. Unimolecular Isomerizations and Oxygen Atom Loss in Formaldehyde and Acetaldehyde Carbonyl Oxides. A Theoretical Investigation. J. Am. Chem. Soc. 1996, 118, 4636−4647. (41) Park, K.; West, A.; Raheja, E.; Sellner, B.; Lischka, H.; Windus, T. L.; Hase, W. L.; S. inglet and Triplet Potential Energy Surfaces for the O2 + C2H4 Reaction. J. Chem. Phys. 2010, 133, 184306. (42) Balucani, N.; Leonori, F.; Casavecchia, P.; Fu, B.; Bowman, J. M. Crossed Molecular Beams and Quasiclassical Trajectory Surface Hopping Studies of the Multichannel Nonadiabatic O(3P) + Ethylene Reaction at High Collision Energy. J. Phys. Chem. A 2015, 119, 12498− 12511. (43) Hu, X.; Hase, W. L. Modification of the Duchovic−Hase− Schlegel Potential Energy Function for H + CH3 ↔ CH4. Comparison of Canonical Variational Transition State Theory, Trajectory, and 4818


Direct Dynamics Simulation of the Thermal - ACS Publications

Recommend Documents