Effects of Roaming Trajectories on the Transition ... - ACS Publications

Article pubs.acs.org/JPCA

Effects of Roaming Trajectories on the Transition State Theory Rates of a Reduced-Dimensional Model of Ketene Isomerization Inga S. Ulusoy,† John F. Stanton,‡ and Rigoberto Hernandez*,† †

Center for Computational and Molecular Science and Technology, School of Chemistry and Biochemistry, Georgia Institute of Technology, Atlanta, Georgia 30332-0400, United States ‡ Department of Chemistry and Biochemistry, The University of Texas at Austin, Austin, Texas 78712-0165, United States ABSTRACT: The rates of chemical reactions (or any activated process) are by definition determined by the flux of reactants (or initial states) that end up as products (or final states). The forward flux through any surface that divides reactants from products is a constant as long as only those trajectories that were reactants in the infinite past and products in the infinite future are included in the flux once and only once. Transition state theory (TST) ignores this last clause, thereby overestimating the rate if any of the trajectories recross the dividing surface. However, its advantage is that it replaces a dynamical calculation with a statistical integral over the TST geometry. The recent identification of roaming trajectoriesthose that persist for a long time as neither reactant nor product without ever visiting near the col on the energy landscapeapparently challenges the dogma that TST’s only error lies in the omission of recrossing trajectories. This question is investigated using the isomerization reaction of ketene in which the experimental values are verified to be in reasonable agreement with both the exact and TST values. We have found two trajectories for the ketene isomerization that carry the signature of roaming, but their effect on the calculation of the reaction rate constant using classical transition state theory is small. Indeed, the existence of roaming trajectories is seen to impose a limitation on which dividing surfaces are appropriate for the calculation of either exact or approximate TST rates, but in this case, they do not unseat the existence of dividing surfaces that can be used safely to calculate TST rates.

■

INTRODUCTION The theoretical description of chemical reactions is tightly connected to the concept of reaction coordinates and transition states. In naive transition-state theory (TST), the reaction rate constant is determined by the potential energy landscape in the vicinity of the transition state of a reaction, as well as that describing the reactants. The transition state is a first-order saddle point along the reaction coordinate, and it gives rise to a hyperdimensional dividing surface (DS) separating the reactant and product regions. Within the TST approximation the transition state DS is assumed to be “recrossing-free,” giving rise to a geometric calculation of the rate.1−4 With the recent and important discovery of roaming pathways,5 a challenge has been posed to TST. Roaming can be described as the meandering motion of one atom or a group of atoms around the molecule, exploring large regions of the PES, followed by an “intramolecular” abstraction of an atom or group of atoms, as in the formation of H2 and CO from H roaming about the HCO fragment in formaldehyde.6−12 The roaming fragment thus bypasses the traditional reaction path, and the energy profile of a roaming reaction may not reveal any relevant dynamics in the saddle point region. Roaming therefore presents an alternative to standard reaction mechanisms in which reactive pathways can stray so far from the internal reaction coordinate (IRC),13 that they not only show non-IRC behavior after crossing a © 2013 American Chemical Society

barrier but also circumvent the barrier present in the IRC altogether.14,15 The challenges of roaming to TST are mainly 2fold: First, is the concept of a dividing recrossing-free surface applicable to a very flat potential? Second, if there exists not only a roaming but also a feasible traditional reaction pathway, how do these compete? And how does the competition change with energy? These questions have been addressed in a series of papers, in which relative rates and branching ratios between roaming and traditional pathways have also been computed.12,16−20 Furthermore, attempts have been made to finding a common DS in configuration space for both traditional and roaming pathways.8,21 In the present paper, we explore the first of these questions. We quantify roaming for a very simple model system and investigate the degree to which TST does in fact provide a useful framework for investigating it. To be specific, a reduced-dimensional model of the ketene isomerization serves as model reaction: H 2CC′O ⇌ OCC′H 2

(1)

Special Issue: Joel M. Bowman Festschrift Received: March 7, 2013 Revised: June 12, 2013 Published: June 17, 2013 7553

dx.doi.org/10.1021/jp402322h | J. Phys. Chem. A 2013, 117, 7553−7560

The Journal of Physical Chemistry A

Article 2 ⎛ d1qF 4 ⎞ ⎜ ⎟ Vcoup(qF ,q1) = k1⎜q1 + k1 ⎟⎠ ⎝

The classical reaction rate constants for the isomerization are calculated, using both the exact expression22,23 for the reaction rateas derived from the flux through the DS between reactant and productand the TST expression. The calculated rate constants are compared to experiment.24,25 By visual inspection of the trajectories leading from the reactant to product, two roaming pathways have been identified for the ketene isomerization. We discuss properties of representative trajectories in the context of roaming and investigate the effect these trajectories have on the reaction rate. Reaction rates for different positions of the DS are calculated to illustrate the importance of the correct choice of the DS.22,23,26 Although we recognize that the present reaction does not involve an intermolecular abstraction and results in closed-shell dissociation products as in “typical” roaming reactions, we use the term more broadly here to indicate that most of the relevant dynamics occur very far from the standard transition state region and on parts of the potential energy surface that are quite flat and relatively featureless. The theoretical framework of this article is discussed in the following section. We then present the results for the rates and for the two roaming trajectories that we have found for the ketene isomerization reaction. Perhaps surprisingly, in this simple system, we find that the TST rates can be obtained accurately, if the DS is chosen judiciously, despite the presence of the roaming trajectories.

The potential contains a 1-dimensional term V1D(qF) along the reaction coordinate qF and a 2-dimensional coupling term Vcoup(qF,q1). The parameters are given in Table 1. We use three Table 1. Parameters Contained in the Potential Function,28 with Three Different Sets of Parameters for the Coupling Potential Vcoup(qF,q1), Giving Rise to Three Different Potentials i, ii, and iii parameter a2 a4 a6 c d k1 k1 d1 k1 d1

THEORY Model System. We construct a model system which is reminiscent of roaming. A common feature in the roaming reactions investigated thus far has been the motion of one atom or a group of atoms around the molecule, followed by intramolecular abstraction.5,7,12,27 An example of the former type of motion is the isomerization of ketene (Figure 1), which proceeds via an oxirene intermediate.

−2.3597 × 10−3 (Eh a0−2) 1.0408 × 10−3 (Eh a0−4) −7.5496 × 10−5 (Eh a0−6) 7.7569 × 10−3 (Eh a0−2) 1.9769 (a0−2) 0.0 (Eh a0−2) 1.0074 × 10−2 (Eh a0−2) 0.0 (Eh a0−5) 1.0074 × 10−2 (Eh a0−2) −2.45182 × 10−4 (Eh a0−5)

Figure 2. Model potentials used in the calculations. Top: potential ii contains a harmonic term in q1. Bottom: potential iii contains a coupling term between the reaction coordinate qF and q1.

Figure 1. Isomerization of ketene via an oxirene intermediate. The potential along the reaction coordinate qF is shallow and exhibits minor features.

The coordinates are referenced to the oxirene geometry, with the reaction coordinate qF corresponding to a coupled in-plane motion of the oxygen and hydrogen atoms. The second mode q1, which is explicitly included in our model, corresponds to a local mode of oxirene:28 It is a linear combination of the a2 CH asymmetric out-of-plane mode and the b1 CH symmetric outof-plane mode, resulting in an out-of-plane motion of one of the hydrogen atoms with a frequency of ω1 = 514 cm−1. Reaction Rates. The reaction rate constant for a classical particle moving from the reactant to product regions on these three potentials is calculated. The exact microcanonical rate

A reduced-dimensional potential for this reaction has already been constructed by Miller et al.28 The potential energy is given by

V1D(qF) = a 2qF 2 + a4qF 4 + a6qF6 + cqF 2 exp−dqF

(i) (ii) (ii) (iii) (iii)

value

different sets of parameters for the coupling term, giving rise to distinct potentials: potential i contains no coupling term (k1 = 0, Vcoup(qF,q1) = 0; Figure 1), potential ii contains a harmonic term in q1 with no coupling to qF (d1 = 0, Vcoup(qF,q1) = Vcoup(q1) = k1q12), and potential iii contains a coupling term between qF and q1. Potentials ii and iii are shown in Figure 2.

■

V (qF ,q1) = V1D(qF) + Vcoup(qF , q1)

(4)

(2) 2

(3) 7554



Article

constant is obtained using the flux through the DS (whose points are denoted through the (Ndof − 1)-dimensional vector qs) between the product and reactant regions on the potential:22 1 kFlux(E) = ρR (E)hNdof

a = 1 − βw(Ep) β=

∫ dp dq s(p·ns) χr (p,q s)·

δ(E−H(p,q s))

w=

G (E ) =

∫ dE ∫ dE δ(E−(E + E ))· ∫ dpim dq ims δ(Eim−Him(pim,q sim))· ∫ dpex dq exs (pex ·ns) χr (pex,q exs )· im

im

im

≈ G (E ) F(E−E )

FFlux(E) =

(13)

N im N ex E+ E Ndof Ndof

(14)

∫ dpex dq exs (pex ·ns) χr (pex,q exs )· (15)

where χr(p has to be replaced by χ+(p for a TST calculation, and the unimolecular rate constant becomes ex

(7a)

ex

,qex s )

kFlux(E) =

,qex s ,ns)

G (E ) ρR (E)hN

ex

−1

(16)

which carries the unit of inverse time. The flux integral in TST approximation for one explicit degree of freedom (model i) can be solved directly1 FTST,i(E) =

(7b)

2(E − V (qs))

(17)

and in two dimensions, the momentum integral can be solved and one obtains FTST,ii/iii(E) = 2

∫ dqs

2(E − V (qs))

(18)

For models ii and iii, the integration of the reactive flux is carried out in local polar coordinates, with the origin of the coordinate system at the point of minimal potential energy along the DS; the flux integral becomes FFlux,ii/iii(E) =

im

G‡(Eim) =

(Ep ≥ 1.0)

)) δ(E−H ex(pex ,q ex s

The first term in product 7b Gim‡(E) is the number of states of the activated complex perpendicular to the explicit degrees of freedom.31 The second term F(E) in product 7b is the flux integral (that is, the numerator in eq 5) restricted to the explicit degrees of freedom. The Nim vibrational modes included implicitly in G‡(Eim) are treated as harmonic oscillator modes, with vibrational frequencies taken from ref 32; the number of states is calculated according to the Whitten−Rabinovitch expression,33−35 im (E + aEzim)N N im N im! ∏i = 1 hνi

(0.1 < Ep < 1.0)

The energies are given relative to an effective zero of energy E0, which we found to lie approximately at 3600 cm−1 below the barrier. At the barrier Eim is available as kinetic energy for the implicit modes. The energy in the explicit modes E ex corresponds to that required to activate the reactant above the threshold to reaction. The difference between these energies, Eex − E‡, is the excess energy. The flux integral of eq 5 is now carried out over the Nex − 1 explicit degrees of freedom

)) δ(E ex −H ex(pex , q ex s ‡

exp( −2.4191Ep1/4)

E = Eim + E ex =

(6)

ex

(11)

Within the RRKM approximation, the total energy E is distributed equally between the implicit and explicit degrees of freedom into the correspondingly designated implicit and explicit energies Eim and Eex

∫ dp dq s(p·ns) χ+ (p,q s,ns)·

im

(Ep ≤ 0.1)

(12)

Here, the characteristic function χr(p,qs) is approximated as χ+(p,qs,ns) = 1 for momenta p·ns > 0. Hence kTST(E) should yield an upper bound to the exact classical rate, as is wellknown.31 The integration in eq 5 is in principle carried out over all Ndof degrees of freedom of the molecule. To reduce the computational effort of the calculation, only one (model i) or two (models ii and iii) degrees of freedom are treated explicitly (Nex i = 1 and Nex ii/iii = 2). Assuming approximate separability between the implicit and explicit variables, the reactive phase space can be divided as ex

(10)

w = (5.0Ep + 2.73 Ep + 3.51)−1

(5)

δ(E−H(p,q s))

Ndof − 1 ⟨ω 2⟩ Ndof ⟨ω⟩2

w = −0.477Ep + 0.25222

The function χr(p,qs) is the characteristic function of reactive phase points, χr(p,qs) = 1 if (p, qs) lies on a reactive trajectory, and χr(p,qs) = 0 otherwise. ns is a unit vector perpendicular to the DS, pointing in the direction of product formation. ρR(E) is the reactant density of states, for which we use the experimentally measured value29,30 4.5 × 104 states/cm−1 for a ketene molecule with 28 000 cm−1 internal energy. The expression in the numerator can be interpreted as the cumulative reaction probability, N(E). In the transition state approximation, the rate simplifies to a form that requires no dynamical trajectories, namely 1 k TST(E) = ρR (E)hNdof

(9)

f (ϕ0) =

(8)

∫ f (ϕ0) dqs ∫0

2π

dϕ(p·n s) χ (p,qs)

mF cos2(ϕ0) + m1 sin 2(ϕ0)

(19) (20)

where ϕ0 is the angle of the DS with the Cartesian x-coordinate and corresponds to the angle between the local polar and Cartesian coordinate systems.

This expression contains a weighted zero-point energy N correction, Ez = ∑i=1dofℏωi/2 and Ep = E/Ez. The parameters a and b and the function w(Ep) are given by33−35 7555



■

Article

RESULTS Microcanonical Rate Constants for the Ketene Isomerization. We calculate the microcanonical reaction rates for the ketene isomerization using the methods described in the previous sections. Experimental rates are taken from the literature.24,25 In these experiments, the isomerization rates were measured indirectly by the exchange yield of isotopically labeled carbon atoms after dissociation into CH2 and CO. Thus, by this means, the isomerization rate can only be measured above the dissociation threshold (about 28 000 cm−1), which lies slightly above the barrier for isomerization. The barrier height for isomerization as estimated from quantum mechanical calculations32 is Ebarr = 348.1 kJ mol−1. We adjust the zero of energy to 3600 cm−1 below the barrier, to account for the equipartitioning of energies, where the energies denoted in the plot are the excess energies defined above. In Figure 3, the calculated classical and the experimental24,25 microcanonical reaction rates are shown. We plot rate versus

TST rate constants kTST,iii(E) are higher than the exact rates kFlux,iii(E). Clearly, for this potential, a fair amount of barrier recrossing occurs, which increases with increasing energy, indicative of the effect that the first saddle point of potential iii becomes less controlling of the reaction. To summarize, we obtain reasonable agreement with the experimental values. This agreement is indicative of the physical relevance of the model; the coupling constants between the explicit modes were obtained by fitting to structures and energies from high-level ab initio calculations.28,32 Furthermore, the calculated rates show the expected trends; specifically, kTST(E) = kFlux(E) for potentials i and ii, and kTST,iii(E) > kFlux,iii(E). It is notable that the experimentally measured rates exhibit additional features that are not present in the calculated rates due to the limitations of the simple one- and two-dimensional models. One such feature is the stepwise increase in the microcanonical rate constant, attributable to the opening of additional channels in which uncoupled vibrational modes are excited. These stepwise features cannot be seen as such in the calculated results due to the approximation of the uncoupled modes as classical harmonic oscillators. A second feature is the presence of prominent resonancessuch as tunneling and Feshbach resonances28in the experimentally observed rates for the ketene isomerization that are not captured in the classical calculation. Ketene Isomerization and Roaming Pathways. Roaming is generally identified by inspecting trajectories. Such roaming trajectories cross through regions of coupled internal motions of the molecule, with pathways of high potential and low kinetic energy that result in internally (vibrationally) excited products.5,7,11,12,27 In typical roaming trajectories, one atom or a group of atoms usually moves far out, leading to a partial fragmentation of the molecule with substantial radical character10,20 at typical separation lengths of ≈6 a0 (3 Å) on a relatively flat region of the potential energy surface.18,20 Inspecting the trajectories obtained with the ketene model potential iii (that for which TST overestimates the exact classical rate), we find trajectories matching these criteria (Figure 4). Trajectory 1 passes from reactant to product region A via regions E of high potential energy, completely avoiding the oxirene intermediate D. Trajectory 2 takes a full turn on the plateau region of the PES before it crosses the second barrier to the product site. These two trajectories are thus exemplary of a

Figure 3. Microcanonical reaction rates for different model potentials: i 1-dimensional, ii 2-dimensional, and iii 2-dimensional with coupling, in comparison with experimental rates.24,25

the energy above the barrier, where the experimental onset of isomerization is at Eexc = 0.0 cm−1. The calculated rate constants are obtained for a DS positioned at the first maximum of the potential (at qF = −2.8 a0; Figure 1), where the trajectories are initiated at the DS. In the 1-d model case i, the obtained reaction rates are in the range of the experimental ones, but the slope of the calculated curve is smaller than the slope of the experimental values. In this model, every trajectory that has an energy larger than the potential energy at the barrier passes directly to product, and at each energy, there is only one configuration that is counted explicitly in the flux integral. Furthermore, the reaction rate in the TST approximation is equal to the exact rate, because every trajectory that passes the barrier is reactive in the classical picture. In the 2-d model ii with no explicit coupling between qF and q1, the calculated rates are close to the experimentally measured ones. Again, every trajectory that passes the barrier crosses to the product side with little deviation from the minimum-energy path, so that the reaction occurs quickly and every trajectory passing the DS is reactive. The curvature of the calculated rates is in good agreement with the experimental one; comparing this result with model i, it seems that at least one additional degree of freedom in the flux integral is necessary to reproduce the reactive phase space. In model iii with coupling between qF and q1, the calculated rates are also in satisfactory agreement with experiment, and the

Figure 4. Trajectories that match the criteria for roaming: trajectory 1 (solid line) crosses the barriers B and enters region E while moving from reactant to product A′, whereas trajectory 2 (dotted line) takes a full turn in the plateau region of the PES via the strongly coupled regions of modes qF and q1 (E and F). Both trajectories do not enter region D of the planar oxirene intermediate. 7556



Article

reaction pathway from reactant to product, which largely skirts the planar oxirene intermediate that appears in the conventional reaction path, and are here categorized as roaming pathways to product. The energies of trajectories 1 and 2 are very high, at E1 = 985 kJ mol−1 (0.375 Eh) and E2 = 1218 kJ mol−1 (0.464 Eh) above the barrier. At such high energies, region B/B′ does not represent as tight a bottleneck as at the low energies in Figure 3; the modulations in the potential as experienced by the trajectories passing through region B/B′ are small in comparison with the total energies; the barrier height at B/B′ relative to the plateau region of the potential around the ′ = 24 kJ mol−1 (0.009 Eh)). oxirene intermediate is Ebarr Trajectory 2 in particular crosses the barrier region B/B′ far from the saddle point, and for both trajectories, the potential energy profile along the pathway does not exhibit a maximum at B/B′. The potential energy profile for trajectory 2 is shown in Figure 5, with the energy profile of the minimum energy path (MEP) as dashed line. The regions and corresponding structures are also shown on Figure 5.

Figure 6. Important structures in the ketene isomerization A to A′. The roaming pathways occur via structures E and F instead of the oxirene intermediate D occurring in the traditional reaction path. In both structures B/B′ (the salient transition state) and C/C′ (formylmethylene) the Hα is pointing out of the plane. In structure E, one hydrogen points down while in F one hydrogen points up.

Figure 5. Potential energy profile along the roaming pathway 2 shown as a dashed line in Figure 4 (solid line) and the minimum energy path (dashed line). The corresponding regions of the potential along the minimum energy path are indicated on the plot, with labels for the region C/C′ partially omitted for lack of space. The roaming trajectory does not exhibit a potential energy maximum at the position of the barrier of the conventional reaction pathway. The energy at the roaming trajectory is very high compared to the subtle features of the potential in the plateau region, which lie below 2000 cm−1 (0.009 Eh).

The molecular structures are shown in Figure 6. The structure of the transition state B is nonplanar, with one of the hydrogen atoms pointing out of the plane, and an angle of 125° between the oxygen−carbon and the carbon−carbon bond.32 In the region of structures near D (oxirene), the molecule is planar, with the oxygen connected to both carbon atoms. Region E and region F correspond to structures in which the oxygen and carbon atoms form an equilateral triangle as in oxirene, with one of the hydrogen atoms pointing out of the plane, according to the definition of mode q1 (see above). The kinetic energy projected parallel and perpendicular to the reaction path is shown in Figure 7 for both trajectories. The letters for the corresponding structures are also shown in the figure. For trajectory 1 (top), the kinetic energy parallel to the reaction path (Tpar) decreases when entering the barrier region at t = 1.0 ps. The trajectory passes through the plateau region, with a constant low value of E = 0.075 Eh. After the passage of the second barrier at t = 3.0 ps, Tpar increases again. The kinetic energy perpendicular to the reaction path, Tperp, oscillates with a frequency of ω ≈ 39 cm−1.

Figure 7. Kinetic energy of the roaming trajectories projected along and perpendicular to the reaction path (Tpar and Tperp). Top: decomposed kinetic energy for trajectory 1 with alphabetic representation of the structures of reactant/product and transition state indicated on the plot. Bottom: decomposed kinetic energy for trajectory 2. The inset shows the decomposed kinetic energy after reaching the product state, where the parallel contribution of the kinetic energy roughly corresponds to the translational energy, and the perpendicular contribution, to the internal kinetic (vibrational) energy.

The projected kinetic energies for roaming trajectory 2 (bottom) exhibit similar features, except that the coupling between the degrees of freedom is stronger. The inset shows the kinetic energy contributions after reaching the product state; here, Tpar corresponds to the translational energy in the first approximation, and Tperp to the vibrational energy of the particle, which oscillates strongly. After following the roaming pathway, the product is thus in a vibrationally highly excited state, with the sense that, in roaming, the product fragments 7557



Article

contain high internal energy. We observe the same long-time behavior for roaming trajectory 1. Such behavior of the kinetic energy is usually taken as a signature of a roaming. Note that the time scale of the roaming trajectories found here is in picoseconds. It does not correspond to the faster time scales of the dominant reaction path trajectories on the twodimensional surface because the former trajectories are so rare in the present case that they contribute little to the overall rate. On the other hand, the time scales of the roaming trajectory are faster than the effective reaction times that would be computed from the inverse over the total rateon the order of microsecondsbecause the latter include the time spent in the perpendicular modes. The model employed here is very simple, representing a molecule with nine degrees of freedom by only two explicit degrees of freedom. The putative roaming trajectories found in this model involve atom−fragment separation distances of almost 6 a0 (3 Å), which is essentially consistent with the general definition of roaming. Although we do not calculate the radical character of the resulting structure for these large separation lengths, it is implicitly included in the model via the potential energy surface that was used for the derivation of the analytic potential. However, not present in the 2-d model is an exit channel for radical dissociation, besides the two existing exit channels. In a higher-dimensionality model, such an exit channel would exist, and there would be a competition between the different exit channels. The roaming trajectories shown above are relatively stable and representative of an ensemble of trajectories. In the presence of a third exit channel, presumably some of the putative roaming trajectories would leave through that exit channel, especially at the high energies encountered here. However, inclusion of this channel would also have an influence on the trajectories passing from the traditional pathway. The rate of roaming versus traditional pathway therefore may even increase. Reaction Rate Constants and Roaming Pathways. We calculate the reaction rates using different sets of DSs to illustrate the effect of the choice of the DS. We concentrate on potential iii at the energies of the roaming pathways, to determine the effect of roaming on the rates. For the rate calculations, the actual energy at which the rates are obtained is higher than the energy of the roaming pathway, due to the equipartitioning of the energy between the explicit and implicit degrees of freedom. The deviation between TST and exact result is a direct measure of the quality of the former calculations. The different sets of DSs are shown in Figure 8. The two roaming trajectories are exemplary for an ensemble of trajectories with very similar initial conditions but still correspond to rare events compared to the total number of reactive trajectories at the given energy in this system. All roaming trajectories pass through the DSs qs1 to qs4, but only DS qs1 is crossed only once by both roaming trajectories. In Figure 8, the potential energy contours corresponding to the total energy of the roaming trajectories are highlighted. The roaming trajectories occur at high energies above the barrier (E1 = 0.375 Eh and E2 = 0.464 Eh), where the potential at the first barrier (at structure B in Figure 4) (which represents the bottleneck at low energies), has widened significantly. The rates for the four DSs at the two different roaming tot energies Etot 1 and E2 are shown in Table 2. The rates at these high energies are several orders of magnitude larger than those rates in Figure 3. This is not surprising because the higher

Figure 8. Different positions of DSs qs on potential iii. The contours of potential energy corresponding to the total energy of the roaming trajectories are shown as dashed lines.

Table 2. Reaction Rates Calculated Using Different Positions of the DS qs for the Reaction Rate Calculationa qs1 qs2 qs3 qs4

−1 kFlux(Etot 1 )/s

−1 kTST(Etot 1 )/s

−1 kFlux(Etot 2 )/s

× × × ×

× × × ×

× × × ×

0.94 0.98 1.01 1.01

20

10 1020 1020 1020

1.50 2.12 2.13 2.14

20

10 1020 1020 1020

4.86 4.97 5.15 5.14

20

10 1020 1020 1020

−1 kTST(Etot 2 )/s

7.99 11.16 11.25 11.26

× × × ×

1020 1020 1020 1020

a

Total energies with equipartitioning between the modes: Etot 1 = 1.68311 Eh, Etot 2 = 2.08446 Eh above the barrier.

energies are 700 times larger than the threshold energies displayed in Figure 3. For each energy, the exact microcanonical rate constant is given by kFlux(E), where its value should in principle be independent of the position of the DS. We observe small differences between the values of kFlux(E) at each energy, which can be attributed to numerical issues. Several trends are observable in the comparison between the exact rate constant and the rate in TST approximation: At tot energy Etot 1 , the TST rate at DS qs1, kTST(E1 ), is a factor of 1.6 higher than the exact rate kFlux(Etot ). Positioning the DS further 1 out, the TST rate increases due to increased trajectory recrossings. For energy E2 at qs1, the difference between TST tot and the exact result is slightly higher than at Etot 1 , and kTST(E2 ) is a factor of 1.65 higher than kFlux(Etot ). For the other DSs, the 2 TST rate is again higher than at qs1. The rate constant in TST approximation for DS qs2 − qs4 is more than twice as high as the exact rate. Still, the reaction rates calculated in the TST approximation at qs1 are in good agreement with the exact calculation of the reactive flux. Changing the geometric position of the DS, which more generally is characterized by its position in phase space, increases the obtained reaction rate constants, meaning that they are less representative of the true rate. Generally, an adequate choice of the DS is crucial for a meaningful application of TST. As long as the roaming trajectory crosses a “bottleneck” that is not recrossed during roaming, then its existence will be accurately accounted for within the TST assumption associated with a DS on that bottleneck. This appears to be apperoximately the case here, and perhaps also more generally. Furthermore, the roaming trajectories represent rare events in the total number of possible reactive trajectories at a given energy: In Figure 9, a histogram over the observed values of q1 is displayed, marking at which q1 the trajectories cross the central region of the potential (qF = 0.0 a0). The roaming trajectories cross this region with a high value, far from the oxirene intermediate for which q1 = 0.0 a0. For energy E1, 92% 7558



Article

exhibit much more complicated dynamics that may not admit a global DS.

■

CONCLUSIONS We demonstrate that microcanonical classical isomerization rates for ketene can be calculated using a two-dimensional model with reasonable accuracy in comparison to experiment. For our model potential, we also observe some trajectories for the ketene isomerization that are reminiscent of what is now known as roaming. These trajectories avoid the planar oxirene intermediate and cross from reactant to product through strongly coupled regions of the potential energy landscape via nonplanar structures. These roaming trajectories meet the criteria of low kinetic versus high potential energy, strong coupling between the degrees of freedom, and a high internal energy of the resulting product molecule. We further demonstrate the importance of the choice of the dividing surface on the quality and applicability of TST. We show that the roaming trajectories are rare events and sample regions of the phase space that are not visited by any other trajectories.

Figure 9. Position in q1 at which all trajectories cross qF = 0.0 a0, for all crossings of only the reactive trajectories (black) and of all trajectories (blue) at the two roaming energies E1 and E2, in percent of the total number of crossings. The corresponding values for the roaming trajectories are indicated as dotted lines. Roaming trajectory 1 crosses qF = 0.0 a0 only once at q1 = −5.6 a0. Roaming trajectory 2 crosses qF = 0.0 a0 three times, twice at q1 = −1.9 a0 and once at q1 = 5.9 a0.

■

of all crossings lie in the range of −5.0 a0 ≤ q1 ≤ 4.0 a0; only 5% of all crossings lie at q1 < −5.0, with roaming trajectory 1 crossing at q1 = −5.6 a0. We observe the same features for energy E2: 93% of all crossings occur between −5.0 a0 ≤ q1 ≤ 5.0 a0, where roaming trajectory 2 crosses at q1 = −1.9 a0 and q1 = 5.9 a0. At q1 > 5.0 a0, only 4% of all crossings take place. This illustrates that the roaming trajectories sample regions in phase space that are far from the region sampled by the majority of the reactive events. The occurrence of roaming trajectories is generally energy dependent. Comparisons between the branching ratios of traditional versus roaming pathways have shown a significant contribution by the latter: as much as 18% (formaldehyde decomposition)16,17 and even up to 84 ± 10% (18 ± 10%) (photodissociation of acetaldehyde).12,18,20 This is not the case in the extremely simple model system studied here, where roaming is relatively rare. To our knowledge, the absolute change in reaction rate due to the opening of a roaming pathway has not been calculated so far, though there have been attempts to construct a combined DS which contains both roaming and traditional pathways.8,21 These results underscore the importance of the choice of the DS.22,23,26 With a poor choice of DS, the estimate obtained by TST will be far from the exact rate constant. To the extent that the DS is well-defined in the sense of globally separating reactants from products, roaming pathways will certainly cross the DS. In general, the challenge posed by roaming pathways is that they will cross the DS multiple times far from the naive saddle region of the reaction. In the present case, the family of roaming trajectoriesof which the ones shown in the figures are representativecan be captured by a single global DS through which they essentially do not recross. Moreover, they appear to be relatively infrequent in the simple models explored here and have little effect on the rate in the case of an optimally chosen DS. Thus the quality of TST is exacerbated in reactions that exhibit roaming because the optimization of the DS requires not only local variations in the vicinity of the saddle but also global variations far from it. This is a fairly straightforward thing to do in the present case, due to its low dimensionality; in larger molecules, it is likely to be much more difficult. For example, high-dimensionality systems such as the fragmentation of acetaldehyde to CH4 and CO11,12 clearly

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work has been partially supported by the Air Force Office of Scientific Research through Grant No. FA9550-12-1-0483 (R.H. and J.F.S). I.U. acknowledges the Alexander von Humboldt Foundation, Germany, for support through a Feodor Lynen Fellowship.

■

REFERENCES

(1) Pechukas, P. Transition State Theory. Annu. Rev. Phys. Chem. 1981, 32, 159−177. (2) Truhlar, D. G.; Garrett, B. C.; Klippenstein, S. J. Current Status of Transition-State Theory. J. Phys. Chem. 1996, 100, 12771. (3) Bartsch, T.; Uzer, T.; Hernandez, R. Stochastic Transition States: Reaction Geometry Amidst Noise. J. Chem. Phys. 2005, 123, 204102. (4) Hernandez, R.; Bartsch, T.; Uzer, T. Transition State Theory in Liquids beyond Planar Dividing Surfaces. Chem. Phys. 2010, 370, 270− 276. (5) Townsend, D.; Lahankar, S. A.; Lee, S. K.; Chambreau, S. D.; Suits, A. G.; Zhang, X.; Rheinecker, J.; Harding, L. B.; Bowman, J. M. The Roaming Atom: Straying from the Reaction Path in Formaldehyde Decomposition. Science 2004, 306, 1158. (6) Zhang, X.; Rheinecker, J. L.; Bowman, J. M. Quasiclassical Trajectory Study of Formaldehyde Unimolecular Dissociation: H[sub 2]CO −> H[sub 2] + CO, H + HCO. J. Chem. Phys. 2005, 122, 114313. (7) Pomerantz, A. E.; Camden, J. P.; Chiou, A. S.; Ausfelder, F.; Chawla, N.; Hase, W. L.; Zare, R. N. Reaction Products with Internal Energy beyond the Kinematic Limit Result from Trajectories Far from the Minimum Energy Path: An Example from H + HBr, H2 + Br. J. Am. Chem. Soc. 2005, 127, 16368−16369. (8) Shepler, B. C.; Han, Y.; Bowman, J. M. Are Roaming and Conventional Saddle Points for H2CO and CH3CHO Dissociation to Molecular Products Isolated from Each Other? J. Phys. Chem. Lett. 2011, 2, 834−838. (9) Bowman, J. M.; Shepler, B. C. Roaming Radicals. Annu. Rev. Phys. Chem. 2011, 62, 531−553. (10) Bowman, J. M.; Suits, A. G. Roaming Reactions: The Third Way. Phys. Today 2011, 64, 33−37. 7559



Article

(11) Houston, P. L.; Kable, S. H. Photodissociation of Acetaldehyde as a Second Example of the Roaming Mechanism. Proc. Natl. Acad. Sci. U. S. A. 2006, 103, 16079−16082. (12) Heazlewood, B. R.; Jordan, M. J. T.; Kable, S. H.; Selby, T. M.; Osborn, D. L.; Shepler, B. C.; Braams, B. J.; Bowman, J. M. Roaming is the Dominant Mechanism for Molecular Products in Acetaldehyde Photodissociation. Proc. Natl. Acad. Sci. U. S. A. 2008, 105, 12719− 12724. (13) Fukui, K. Formulation of the Reaction Coordinate. J. Phys. Chem. 1970, 74, 4161−4163. (14) Sun, L.; Song, K.; Hase, W. L. A SN2 Reaction That Avoids Its Deep Potential Energy Minimum. Science 2002, 296, 875−878. (15) Lourderaj, U.; Hase, W. L. Theoretical and Computational Studies of Non-RRKM Unimolecular Dynamics. J. Phys. Chem. A 2009, 113, 2236−2253. (16) Lahankar, S. A.; Chambreau, S. D.; Zhang, X.; Bowman, J. M.; Suits, A. G. Energy Dependence of the Roaming Atom Pathway in Formaldehyde Decomposition. J. Chem. Phys. 2007, 126, 044314(1)− 044314(8). (17) Lahankar, S. A.; Chambreau, S. D.; Townsend, D.; Suits, F.; Farnum, J.; Zhang, X.; Bowman, J. M.; Suits, A. G. The Roaming Atom Pathway in Formaldehyde Decomposition. J. Chem. Phys. 2006, 125, 044303(01)−044303(10). (18) Harding, L. B.; Georgievskii, Y.; Klippenstein, S. J. Roaming Radical Kinetics in the Decomposition of Acetaldehyde. J. Phys. Chem. A 2010, 114, 765−777. (19) Harding, L. B.; Klippenstein, S. J. Roaming Radical Pathways for the Decomposition of Alkanes. J. Phys. Chem. Lett. 2010, 1, 3016− 3020. (20) Klippenstein, S. J.; Georgievskii, Y.; Harding, L. B. Statistical Theory for the Kinetics and Dynamics of Roaming Reactions. J. Phys. Chem. A 2011, 115, 14370−14381. (21) Harding, L. B.; Klippenstein, S. J.; Jasper, A. W. Separability of Tight and Roaming Pathways to Molecular Decomposition. J. Phys. Chem. A 2012, 116, 6967−6982. (22) Miller, W. H. Direct and correct calculation of canonical and microcanonical rate constants for chemical reactions. J. Phys. Chem. A 1998, 102, 793−806. (23) Pechukas, P. In Modern Theoretical Chemistry; Miller, W. H., Ed.; Plenum: New York, 1976; Vol. 2, pp 269−322. (24) Lovejoy, E. R.; Kim, S. K.; Alvarez, R. A.; Moore, C. B. Kinetics of Intramolecular Carbon Atom Exchange in Ketene. J. Chem. Phys. 1991, 95, 4081−4093. (25) Lovejoy, E. R.; Moore, C. B. Structures in the Energy Dependence of the Rate Constant for Ketene Isomerization. J. Chem. Phys. 1993, 98, 7846−7854. (26) Bartsch, T.; Uzer, T.; Moix, J. M.; Hernandez, R. Identifying Reactive Trajectories Using a Moving Transition State. J. Chem. Phys. 2006, 124, 244310-01−13. (27) Shepler, B. C.; Braams, B. J.; Bowman, J. M. Roaming Dynamics in CH3CHO Photodissociation Revealed on a Global Potential Energy Surface. J. Phys. Chem. A 2008, 112, 9344−9351. (28) Gezelter, J. D.; Miller, W. H. Resonant Features in the Energy Dependence of the Rate of Ketene Isomerization. J. Chem. Phys. 1995, 103, 7868−7876. (29) Kim, S. K.; Lovejoy, E. R.; Moore, C. B. Transition State Vibrational Level Thresholds for the Dissociation of Triplet Ketene. J. Chem. Phys. 1995, 102, 3202−3219. (30) Lovejoy, E. R.; Kim, S. J.; Moore, C. B. Observation of Transition-State Vibrational Thresholds in the Rate of Dissociation of Ketene. Science 1992, 256, 1541. (31) Hernandez, R. A Combined Use of Perturbation Theory and Diagonalization: Application to Bound Energy Levels and Semiclassical Rate Theory. J. Chem. Phys. 1994, 101, 9534−9547. (32) Scott, A. P.; Nobes, R. H.; Schaefer, H. F., III; Radom, L. The Wolff Rearrangement: The Relevant Portion of the Oxirene-Ketene Potential Energy Hypersurface. J. Am. Chem. Soc. 1994, 116, 10159− 10164.

(33) Whitten, G. Z.; Rabinovitch, B. S. Accurate and Facile Approximation for Vibrational Energy-Level Sums. J. Chem. Phys. 1963, 38, 2466−2473. (34) Whitten, G. Z.; Rabinovitch, B. S. Approximation for RotationVibration Energy Level Sums. J. Chem. Phys. 1964, 41, 1883−1883. (35) Tardy, D. C.; Rabinovitch, B. S.; Whitten, G. Z. VibrationRotation Energy-Level Density Calculations. J. Chem. Phys. 1968, 48, 1427−1429.

7560


Effects of Roaming Trajectories on the Transition ... - ACS Publications

Recommend Documents