Nonstatistical effects in the unimolecular dissociation ...

Nonstatistical effects in the unimolecular dissociation of the acetyl radical Angeles Peña-Gallego, Emilio Martı́nez-Núñez, and Saulo A. Vázquez Citation: The Journal of Chemical Physics 110, 11323 (1999); doi: 10.1063/1.479073 View online: http://dx.doi.org/10.1063/1.479073 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/110/23?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Quasiclassical trajectory study of formaldehyde unimolecular dissociation: H 2 C O → H 2 + C O , H + H C O J. Chem. Phys. 122, 114313 (2005); 10.1063/1.1872838 Dissociation dynamics of the Ã 2 A ″ state of vinyl radical J. Chem. Phys. 118, 4452 (2003); 10.1063/1.1542878 Competing isomeric product channels in the 193 nm photodissociation of 2-chloropropene and in the unimolecular dissociation of the 2-propenyl radical J. Chem. Phys. 114, 4505 (2001); 10.1063/1.1345877 Theoretical study of an isotope effect on rate constants for the CH 3 + H 2 → CH 4 + H and CD 3 + H 2 → CD 3 H+H reactions using variational transition state theory and the multidimensional semiclassical tunneling method J. Chem. Phys. 110, 10830 (1999); 10.1063/1.479025 Unimolecular dissociation dynamics of highly vibrationally excited DCO (X̃ 2 A) . II. Calculation of resonance energies and widths and comparison with high-resolution spectroscopic data J. Chem. Phys. 106, 5359 (1997); 10.1063/1.473599

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 193.144.81.13 On: Fri, 21 Nov 2014 07:50:10

JOURNAL OF CHEMICAL PHYSICS

VOLUME 110, NUMBER 23

15 JUNE 1999

Nonstatistical effects in the unimolecular dissociation of the acetyl radical Angeles Penã-Gallego, Emilio Martıńez-Nuń˜ez, and Saulo A. Va´zquez Departamento de Quı´mica Fı´sica, Universidad de Santiago de Compostela, Santiago de Compostela E-15706, Spain

~Received 9 February 1999; accepted 24 March 1999! Classical trajectory and statistical variational efficient microcanonical sampling transition state theory calculations were carried out to investigate the dissociation dynamics of the acetyl radical. For this purpose, an analytical potential function was developed based on ab initio and experimental data reported in the literature. This potential function reproduces reasonably well the geometries, frequencies, and energies of the stationary points of the ground state potential energy surface. The dynamics of the reaction was shown to be intrinsically non-Rice–Ramsperger–Kassel–Marcus ~RRKM! at high energies and particularly at 65.9 kcal/mol, at which experimental work showed evidence for nonstatistical behavior. On the other hand, initial excitations of normal modes 507 ~CCO bend!, 1079 ~CC stretch!, 1504 (CH3 umbrella vibration!, and 1939 ~CO stretch! enhance significantly the rate of reaction; specifically, excitation of the CO stretch gives a rate coefficient an order of magnitude higher than the rate obtained under random initial conditions. These mode specific effects are explained in terms of a restricted intramolecular vibrational redistribution ~IVR!. Under statistical initial conditions, the classical trajectory calculations showed a normal isotope effect at the two lowest energies studied, and a slight inverse isotope effect at 65.9 kcal/mol, a result that can be explained with the presence of a methyl free-rotor at the transition state. In contrast, upon initial excitation of the CC and CO stretches and CCO bending at 65.9 kcal/mol, the calculations predicted a normal isotope effect, which agrees with the experimental findings. © 1999 American Institute of Physics. @S0021-9606~99!01223-4#

I. INTRODUCTION

The Rice–Ramsperger–Kassel–Marcus ~RRKM! theory, one of the major accomplishments of chemical kinetics, assumes that all energetically accessible states of the reactant are equally probable ~random lifetime assumption!, that is, the vibrational energy is rapidly randomized on the time scale of reaction throughout all vibrational degrees of freedom so that molecules with a particular energy have a lifetime distribution given by P(t)5k RRKM exp(2kRRKMt), where k RRKM is the RRKM rate constant.1 However, at high energies, both the intramolecular energy redistribution ~IVR! and the reaction take place competitively so that the applicability of the RRKM theory is questionable.2 Moreover, if the barrier to reaction is low in comparison with the initial excitation energy, those molecules prepared in states that are close to the transition state dividing surface can dissociate before the energy randomization occurs, which implies a violation of the RRKM postulates. The acetyl radical (CH3CO) dissociation process is a good candidate to present this type of nonstatistical effects, as the bond scission to form CH3 and CO is a reaction with a dissociation barrier as low as 1761 kcal/mol. 3,4 Several investigations examined the fragmentation of the acetyl radical as an intermediate in the photodissociation of simple ketones.3–14 Recently, Owrutsky and Baronavski12 in their ultrafast deep UV mass-resolved photoionization spectroscopy study of acetone concluded that the second step of the overall decomposition, that is, the acetyl dissociation, may be described by statistical models. Based on RRKM calcula0021-9606/99/110(23)/11323/12/$15.00

tions, they estimated that, under their experimental conditions, the average acetyl internal energy was about 25 kcal/ mol. In contrast, at higher internal energies ~concretely, about 65610 kcal/mol), Kim et al.13 had previously found evidence for nonstatistical dynamics. By femtosecond mass spectrometry, they obtained a time constant of 500 fs for the acetyl radical, which is about ten times higher than that computed by RRKM calculations13 using ab initio vibrational frequencies. They concluded that ‘‘this dramatic difference indicates that the vibrational energy of the acetyl is prevented from reaching the C–C stretch reaction coordinate by a slow or restricted IVR.’’ Furthermore, the isotope effects observed in their study cannot be explained by RRKM calculations. It is interesting to note here that in their experiment the preparation of the vibrationally excited, groundstate acetyl radical involved preferentially the motion of the central carbon atom in the CCO moiety. Therefore, under their experimental conditions, the energy distribution in the nascent acetyl radical is nonrandom and so its dissociation dynamics appears to exhibit, in terms defined by Bunker and Hase,2 apparent non-RRKM behavior. Shibata et al.14 have selected the acetyl radical to present a novel experimental technique to measure the energy-dependent unimolecular dissociation rate k(E) of radical species. They found that, at energies near the threshold, the observed dissociation rates were an order of magnitude smaller than those calculated by RRKM theory, and concluded that IVR and/or vibration– rotation energy exchange prior to C–C bond rupture might be restricted. The above experimental studies motivated us to investi-

11323

© 1999 American Institute of Physics


11324

Penã-Gallego, Martıńez-Nun˜ez, and Va´zquez

J. Chem. Phys., Vol. 110, No. 23, 15 June 1999

gate the dissociation dynamics of the acetyl radical by the classical trajectory15 and the efficient microcanonical sampling–transition state theory ~EMS–TST!16 methods, and the results are presented here. To perform the calculations, we first constructed an analytical potential energy surface ~PES! that reproduces adequately the available experimental3,4,17,18 and ab initio19,20 data. We then calculated rate coefficients k(E) for several ensembles over a wide energy range. Since both methods employed here are classical and both use the same PES, a comparison of rate coefficients evaluated by these two methods provides a straightforward way to explore nonstatistical effects in the dynamics of the acetyl dissociation. In addition, we have analyzed the mechanisms of IVR for ensembles in which a normal mode is initially excited, by incorporating the normal and local mode approximations within the framework of classical trajectories. Finally, we have calculated the rate coefficients for the trideuterated molecule at several energies and for different ensembles. A comparison between the trajectory, RRK, RRKM, and experimentally determined values of k H /k D is also discussed. II. POTENTIAL ENERGY SURFACE

We have developed an analytical PES based on experimental and theoretical data available in the literature.3,4,17–20 The general expression of our PES reads

3

1

( i51

(

i51 3

V ~ u CCH~ i ! ! 1

V ~ R CH~ i ! ! 1V ~ u CCO! 3

( ( i51 j5i11

( V ~ t OCCH~ i ! !

i51

1V ~ R CC , u CCO! ,

~1!

where the subscripts i and j refer to the hydrogen atoms ~1, 2, or 3! attached to the carbon atom. In order to obtain a correct representation of the system at all regions of the PES, many parameters for these functions have to vary significantly as the reactant evolves to products. For this reason, we introduced switching functions with the general expression: SW~ i ! 5exp$ 2C i @ R CC2R CC,eq# 2 % ,

~2!

where R CC is the C–C distance, R CC,eq is the CC equilibrium bond length, and C i is a parameter whose value depends on the particular interaction term. The individual terms of our potential are detailed as follows. The carbon–oxygen stretching interaction V(R CO) is represented by a Morse function: V ~ R CO! 5D CO$ 12exp@ 2 a CO~ R CO2R CO,eq!# % 2D CO , 2

~3!

where D CO is the potential well depth, and a CO the curvature parameter of the CO bond. The CO bond length parameter was varied by SW~1! to alter the CO equilibrium bond length as the acetyl radical dissociates, p p r R CO,eq5R CO 2 ~ R CO 2R CO ! SW~ 1 ! .

The carbon–hydrogen interaction terms are V ~ R CH~ i ! ! 5D CH~ i ! $ 12exp@ 2 a CH~ i ! ~ R CH~ i ! 2R CH,eq!# % 2 2D CH~ i ! .

~6!

For the sake of simplicity, we employed the same parameter value R CH,eq for all the CH bonds in the reactant and the product CH3. Bending interactions are described by harmonic functions: V ~ u ! 50.5K u ~ u 2 u eq! 2 .

~7!

The force constant for the CCO bend was modeled as follows: r SW~ 2 ! . K CCO5K CCO

~8!

Thus, the CCO bending term vanishes as the CC bond dissociates. The CCH(i) bends were also attenuated with SW~2!, r SW~ 2 ! . K CCH~ i ! 5K CCH

p p r u HCH,eq5 u HCH 2 ~ u HCH 2 u HCH ! SW~ 3 ! .

V ~ u H~ i ! CH~ j ! !

3

1V ~ v H~1!CH~2!H~3!! 1

V ~ R CC! 5D CC$ 12exp@ 2 a CC~ R CC2R CC,eq!# % 2 2D CC . ~5!

~9!

For the three H(i)CH( j) bending interactions we employed the switching function SW~3!, which varies the HCH equilibrium bond angle,

3

V5V ~ R CO! 1V ~ R CC! 1

The superscripts p and r denote equilibrium values for product and reactant, respectively. V(R CC) is the carbon–carbon stretching term taken as

~4!

~10!

This switching function changes the equilibrium bond angle u HCH,eq from 109.5° ~the value in the acetyl radical! to 120° ~the value in the product CH3). For simplicity, we have taken the same force constant K HCH for the reactant CH3CO and the product CH3. To model the wagging interaction in the product CH3, we included the following term in the potential function V ~ v H~1!CH~2!H~3!! 50.5K v ~ v ! 2 ,

~11!

where K v 5 $ K vp 1exp@ 2C v ~ R CC2R CC,eq! 2 # % 3 $ 12exp@ 2C v ~ R CC2R CC,eq! 2 # % ,

~12!

K vp

is the wagging force constant of the product CH3. where This expression of K v ensures that the wagging interaction is not present before the CC bond cleavage occurs. In our PES, the torsion is described by means of a cosine series for each OCCH( j) arrangement ~three in all!: 5

V~ t !5

( a i cos~ i t ! ,

i50

~13!

where the torsional terms were attenuated with SW~4!, a i 5a ri SW~ 4 ! .

~14!

The switching function here guaranties that the torsional interaction vanishes as the CC bond dissociates.




11325

TABLE I. Equilibrium geometries,a energetics, and frequencies for the chemical species involved in this study. Internal coordinate CO CC CHa CHb CCO Ha CC Hb CC Ha CHb Hb CHb Ha CCO Hb CCO

CH3CO PESb

ab initioc

PESb

1.187 1.508 1.088 1.088 128.1 109.5 109.5 109.5 109.5 0.0 240.0

1.187 1.508 1.089 1.087 128.1 111.1 108.4

1.147 2.085 1.088 1.088 115.3 94.6 94.6 119.5 119.5 0.0 240.0

0.0 238.2

ab initioc

PESb

1.149 2.125

1.139

ab initioc

exptd

1.139

1.128

1.088 1.088

1.076 1.076 1.076

115.9

117.4 117.4

120.0 120.0

120.0 120.0

Harmonic frequencies/cm21

Assignment (a 9 )torsion (a 8 )CCO bend (a 8 )CH3 sym. rock (a 9 )CH3 asym. rock (a 8 )CC stretche (a 8 )CH3 umbrella vibr.e (a 9 )CH3 asym. scissors (a 8 )CH3 sym. scissors (a 8 )CO stretch (a 8 )CH3 sym. stretch (a 8 )CH3 asym. stretch (a 9 )CH3 asym. stretch

CH31CO

TS

115 507 958 975 1079 1504 1476 1477 1939 2799 2908 2909

101 468 1062 960 884 1361 1477 1478 1928 3188 3193 3193

44 394 815 821 373i 1793 1618 1618 1876 2754 2928 2929

43 276 566 568 496i 901 1445 1452 2027 3321 3325 3325

894 1610 1610 1827 2752 2928 2928

606 1396 1396

3161 3161

Relative energyf/kcal/mol 0.0

0.0

17.04

17.25

8.05

8.05

a

Distances in Å and angles in degrees. This study. c MP2/cc-pVTZ data and G2 results from Refs. 19 and 20, respectively. d References 17 and 18. e The CC stretch and the CH3 umbrella vibration correspond to the reaction coordinate and to the wagging motion, respectively, as the reaction proceeds. f Including ZPE. b

Ab initio calculations19 predict that the CCO angle decreases by about 15° from the reactant to the transition state, resulting in the tilting of CO toward CH3. For this reason, we have added an interaction term between the CC stretching and the CCO bending internal coordinates. This interaction term reads

librium bond angles! were taken from ab initio structures (CH3CO and CO! optimized at the MP2/cc-pVTZ level of theory.19 Initial values for force constants, curvature parameters, and potential well depths were first taken from the literature, but then refined here to reproduce more accurately the vibrational frequencies and relative energies. Parameters C v and C i were fitted to model the energetics, the geometry V ~ R CC , u CCO! 5k CC,CCO~ R CC2R CC,eq!~ u CCO2 u CCO,eq! , of the transition state, and the ab initio frequencies given in ~15! Table I. Finally, the torsional parameters were fitted to rewhere the force constant k CC,CCO was attenuated by produce the torsional vibrational frequency. r As seen in Table I, our equilibrium geometries, energies, ~16! k CC,CCO5k CC,CCOSW~ 5 ! . and vibrational frequencies compare reasonably well with For dissociation to take place, the CC distance must become the corresponding ab initio and experimental values. For gelarge. Since k CC,CCO was chosen to be positive, a decrease in ometries, it can be seen that the structural parameters of the the potential energy is achieved if the CCO bond angle direactant conformation (CH3CO) calculated by the analytical minishes. PES are in nearly exact accord with the ab initio values. For Table I lists geometrical features, frequencies, and relathe transition state, a difference of 0.04 Å was found for the tive energies taken from the literature3,4,17–20 that were used CC bond distance and of 2.1° for the HCH bond angles. On in this work either for the parametrization of our PES or for the other hand, the relative energies computed by our PES the purposes of comparison. Also included in Table I are the agree very well with the literature values. The barrier height corresponding values obtained by the GENDYN code21 and for the reaction, including the zero-point energies ~ZPE!, our PES with the final parameters shown in Table II. Most was computed to be 17.04 kcal/mol by our PES, a value in geometrical parameters ~equilibrium bond lengths and equiaccord with that predicted by G2 calculations ~17.25 20 at: http://scitation.aip.org/termsconditions. Downloaded to IP: This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms and that estimated experimentally by North kcal/mol! 193.144.81.13 On: Fri, 21 Nov 2014 07:50:10

11326


J. Chem. Phys., Vol. 110, No. 23, 15 June 1999 TABLE II. Parameters for the analytical PES. Stretching parameters Parameter r R CO p R CO

R CC,eq R CH,eq

Value/Å

Parameter

Value/Å21

Parameter

Value/~kcal/mol!

1.187 1.139 1.508 1.088

a CO a CC a CH

2.46 4.26 1.76

D CO D CC D CH

160.0 13.90 105.0

Parameter

Value/degrees

u CCO,eq u CCH,eq r u HCH p u HCH

128.1 109.5 109.5 120.0

Bending parameters Parameter

Value/~mdyn Å/rad2!

r K CCO r K CCH K HCH

1.00 0.60 0.55

Wagging parameter Parameter

Value/~mdyn Å/rad2!

Kv

0.05 Torsional parameters

Coordinate

a r0 /eV

a r1 /eV

a r2 /eV

a r3 /eV

a r4 /eV

a r5 /eV

HCCO

0.0031

0

0

20.0031

0

0

Interaction parameter Parameter

Value/~mdyn/rad!

r K CC,CCO

0.4 Switching functions parameters

Parameter

Value/Å22

Parameter

Value/Å22

C1 C2 C3

5.50 0.90 0.10

C4 C5 Cv

5.00 1.00 0.15

et al.3 ~1761 kcal/mol!. The energy difference between reactant and products computed by our PES is in exact accord with the G2 calculation ~8.05 kcal/mol!.20 Finally, the vibrational frequencies computed numerically by the analytical potential are in reasonable agreement with the ab initio and experimental frequencies. Notable discrepancies, however, appear for the CC stretch and methyl umbrella vibrational modes. III. DYNAMICAL STUDY A. Computational methods

1. Trajectory calculations

Microcanonical rate coefficients for the CH3CO dissociation reaction at total energies ranging from 28.4 to 73.4 kcal/mol above the zero-point energy (ZPE526.6 kcal/mol) were evaluated by classical trajectory calculations employing the GENDYN code.21 The following standard options in GENDYN were selected for the initial conditions: ~a!

Efficient microcanonical sampling ~EMS!22 with the angular momentum restricted to zero. With this option a Markov chain was initiated from the equilibrium geometry taking an initial incubation period of 53105 steps to ensure that the system was relaxed away from the initial configuration. Then a trajectory was integrated until the CC bond distance reached 5 Å, or until

~b!

the simulation time completed 5 ps. A sequence of 10 000 Markov moves was taken from the starting point of the previous trajectory, and the integration/ Markov walk pattern was repeated until an ensemble of 1000 trajectories was generated. At each step of the Markov chain all atoms were moved with a step size of 0.06 Å, leading to an acceptance ratio of about 0.5 for the Markov chain. To confine the sampling to reactant configuration space, maximum bond extensions allowed during the random walk were 2 Å for all the bonds. The angular momentum was restricted in the initial condition selections of these ensembles, making these calculations comparable to those in which normal modes are excited. Normal-mode excitations. With this second method, the initial state of the system is assumed to be such that 65.9 kcal/mol exists into one of the 12 normal modes. An ensemble of phase-averaged classical trajectories was integrated until the CC bond distance reached 5 Å, or until the time of simulation completed 5 ps.

For both types of initial conditions, the lifetime of a trajectory was taken to be the integration time up to the last inner turning point of the CC vibration. Trajectories were integrated by using a fourth-order Runge–Kutta–Gill routine with a fixed step size of 0.08 fs, resulting in an energy conservation of at least five significant figures.




2. Statistical calculations

TABLE III. Markov walk parameters employed in the EMS–TST calculations.

Microcanonical rate coefficients were also evaluated for our PES by means of the efficient microcanonical sampling– transition state theory ~EMS–TST!.16 Following is a brief summary of this method. The microcanonical rate coefficient of a unimolecular reaction can be expressed as an average over the microcanonical ensemble23 k~ E !5

1 * dG d @ H ~ G ! 2E # d ~ q RC2q C ! u q RCu , 2 * dG d @ H ~ G ! 2E #

~17!

where G is the complete set of position and momentum coordinates $q, p%, H(G) is the Hamiltonian of the system excluding the center of mass motion, q RC5q RC(q) is the reaction coordinate, q C is the critical value of q RC required for reaction, and q˙ RC5q˙ RC(q) is the velocity through the critical surface. The integrals over G in Eq. ~17! are understood to be over the reactant part of phase space. The above equation can be rewritten in a form appropriate to be evaluated by EMS algorithms,22 k

1 * dq W ~ q! d ~ q RC2q C ! ^ u q RCu & , ~ E !5 2 * dq W ~ q!

EMS–TST

~18!

where W(q) is the efficient microcanonical sampling weight factor. The appropriate (E,J50) statistical weight reads W ~ q! 5 ~ I a I b I c ! 21/2@ E2V ~ q!# ~ 3N28 ! /2,

~19!

where V(q) is the potential energy at the current configuration, N is the total number of atoms, and I a I b I c is the product of the principal moments of inertia at the configuration q. ^ u q˙ RCu & is the average absolute velocity along the reaction coordinate, given by

^ u q RCu & 5

* dp d @ T ~ p! 2K # u q RCu , * dq d @ T ~ p! 2K #

~20!

where K5E2V(q) is the kinetic energy at configuration q. In order to increase the convergence rate of Eq. ~18!, it is usual to introduce a function with the same form as the EMS weight: the importance sampling function I(q), a

I ~ q! 5 @ E2V ~ q!# ,

~21!

where a is an adjustable parameter. This yields an effective weight factor of W eff~ q! 5 ~ I a I b I c ! 21/2@ E2V ~ q!# @~ 3N28 ! /2# 2 a .

11327

~22!

Good convergence of both numerator and denominator is often obtained16,24 with a value of a such that the exponent in Eq. ~22! amounts between 1 and 1.5. Accordingly, in this work the value of a was chosen to be 4, which for CH3CO (N56 atoms! makes the exponent in Eq. ~22! equal to 1. Hereafter the rate coefficient evaluated from Eq. ~18! will be referred to as the statistical rate coefficient, k stat(E). To complete the EMS–TST procedure, the value of k stat(E) should be minimized with respect to the location of the transition-state dividing surface.25 Details of the minimization routine are given elsewhere.16 In this work, the final value of k stat was obtained by an average over six different values computed with different initial random number seeds.

Markov walk parameter Number of warm-up Markov steps Number of steps in Markov walk Number of atoms moved per Markov step Markov step size/Å Important sampling exponent, a Number of transition state dividing surfaces Spacing between dividing surfaces/Å Width of each dividing surface/Å Location of innermost dividing surface/Å

Value 50 000 10 000 000 6 0.06 4.0 8 0.1 0.1 1.8

In the present case ~a simple bond scission!, the average absolute velocity across the transition-state dividing surface can be evaluated analytically16 by

^ u q RCu & 5

S D 2K pm

1/2

@~ 3N25 ! /2! ]! , @~ 3N24 ! /2! ]!

~23!

where m is the reduced mass for motion along the reaction coordinate. In our study, the reaction coordinate q RC was the length of the dissociating bond and so the dividing surfaces were spheres of radii equal to q C . The reduced mass m was calculated from the masses of the two carbon atoms that form the dissociating bond. Table III gives the Markov walk parameters used in the EMS–TST calculations. B. Results and discussion

1. EMS initial conditions

Decay curves for EMS initial conditions at energies ranging from 28.4 to 73.4 kcal/mol ~above ZPE! are depicted in Fig. 1. Assuming single-exponential behavior, the trajectory rate coefficients for EMS initialization can be evaluated by linear least-squares fits of the trajectory decay curves to the equation ln P52k trajt,

~24!

where P is the fraction of unreacted trajectories at time t. The rate coefficients thus obtained are listed in Table IV, and in what follows they will be identified as k traj. As seen in Fig. 1, nonlinear behavior arises as the energy increases; as a consequence, the single-exponential decay fit is not adequate except for the lowest energies ~up to '41 kcal/mol!. Therefore, k traj is only a rough, time-averaged rate coefficient. For the high-energy ensembles, the trajectory lifetime distribution is nonrandom or, in other words, the initial microcanonical distribution is not maintained during the dissociation. This situation arises when the rate of IVR is lower than the rate of the reaction, that is, when the process of IVR is the rate-determining step. In such cases, the initial dissociation rate, which is the statistical rate coefficient, is faster than the dynamical rate due to long-lived trajectories trapped by phase–space bottlenecks. This kind of nonstatistical behavior was called by Bunker and Hase intrinsic non-RRKM behavior.2 Shalashilin and Thompson26 proposed an approach, called intramolecular dynamics diffusion theory ~IDDT!, to quantify the extent of intrinsic non-RRKM be-


11328



TABLE IV. Trajectory and statistical rate coefficients for the dissociation of the acetyl radical. Ea/~kcal/mol!

(k traj/ps21) b

(k stat/ps21) c

28.4 30.9 33.4 35.9 38.4 40.9 43.4 45.9 48.4 50.9 53.4 55.9 58.4 60.9 63.4 65.9 68.4 70.9 73.4

0.29860.031 0.38660.031 0.47660.020 0.55460.045 0.64060.031 0.68260.077 0.81760.089 0.90560.055 0.93460.126 1.02260.155 1.00560.195 1.24860.202 1.28760.249 1.32560.257 1.46160.243 1.37260.235 1.63460.255 1.63860.190 1.71960.265

0.27060.018 0.35260.016 0.43760.019 0.52260.022 0.61860.050 0.73360.034 0.85560.059 0.95460.084 1.06860.069 1.19960.082 1.32160.052 1.49760.089 1.60360.080 1.74860.140 1.84060.150 2.03960.120 2.12060.115 2.25760.136 2.37360.136

a

Above ZPE. Errors are the standard deviations of the fits. c The values of k stat is the average of six different values, and the errors are the standard deviations ~see the text!. b

FIG. 1. Decay curves computed under EMS initial conditions at ~a! kcal/mol, ~b! 30.9 kcal/mol, ~c! 33.4 kcal/mol ~d! 35.9 kcal/mol, ~e! kcal/mol, ~f! 40.9 kcal/mol, ~g! 43.4 kcal/mol, ~h! 45.9 kcal/mol, ~i! kcal/mol, ~j! 50.9 kcal/mol, ~k! 53.4 kcal/mol, ~l! 55.9 kcal/mol, ~m! kcal/mol, ~n! 60.9 kcal/mol, ~o! 63.4 kcal/mol, ~p! 65.9 kcal/mol, ~q! kcal/mol, ~r! 70.9 kcal/mol, and ~s! 73.4 kcal/mol.

28.4 38.4 48.4 58.4 68.4

efficients as a function of the location of the transition state dividing surface q C is displayed graphically in Fig. 2. For all the energies selected in this study, the minimum in the k versus q C profile was found at 2.1 Å, which is the computed value of the CC distance at the transition state by both ab initio theory and our model PES. In Table IV, the microcanonical rate coefficients obtained from Eq. ~18! are compared with those calculated by trajectory simulations. As can be seen, the EMS–TST rates match the trajectory rates, within statistical uncertainties, at the lowest energies of this study. Therefore, at these energies the system behaves statistically. However, when the energy is increased, the dissociation of the acetyl radical becomes nonexponential and the

havior. The high-energy nonstatistical rate was referred26 to as the dynamical rate coefficient k dyn(E,t). The rate coefficient obtained in this work from Eq. ~24! is related to the dynamical rate coefficient by a time average, that is, k traj(E)5 ^ k dyn(E,t) & . On the other hand, if a microcanonical ensemble is sampled, the statistical rate k stat(E) is the rate coefficient for the initial decomposition of the ensemble @ k stat(E)5k dyn(E,t50) # . Shalashilin and Thompson26 used short-time classical trajectories and Monte Carlo transition state theory ~MCTST! calculations to evaluate the nonstatistical rate k dyn(E,t). The initial statistical rates are computed by using MCTST, whereas the long-time ~IVR-limited! rates are calculated by the IDDT method.26 The variation of the EMS–TST microcanonical rate co-

FIG. 2. EMS–TST results for the dissociation rate coefficient as a function of the transition state-dividing surface location, at the energies selected in this study. This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP:

193.144.81.13 On: Fri, 21 Nov 2014 07:50:10



FIG. 3. Trajectory ~circles!, EMS–TST ~triangles!, and theoretical RRK ~solid line! rate coefficients displayed as a graph of log(k(E)) vs log(1 2E0 /E). The RRK values were calculated with the theoretical value s 512.

statistical rates no longer match the trajectory rates. In particular, for the highest energy selected in this study ~73.4 kcal/mol!, the ratio k stat/k traj is 1.54. Figure 3 shows a comparison between trajectory, statistical ~EMS–TST!, and RRK rate coefficients. They are displayed as a graph of log@k(E)# against log(12E0 /E), where E 0 is the potential energy difference between the transition state and the reactant, neglecting zero-point energy. According to RRK theory, the rate coefficient varies with energy as log~ k ~ E !! 5log n 1 ~ s21 ! log~ 12E 0 /E ! ,

~25!

where s is the number of active oscillators in the molecule ~12 for the acetyl radical! and n is a frequency factor, which will be defined later. The solid line was obtained by a fit of the three lowest trajectory rate coefficients to the RRK expression with the constraint of s512. The broken line was obtained by a fit of the EMS–TST rate coefficients to the RRK expression with the same constraint. From these two fits we obtained the following values for the frequency factors: 1013.44 s21 from the trajectory data and 1013.39 s21 from the EMS–TST results. As seen in Fig. 3, the statistical rate coefficients agree with the RRK behavior, whereas the trajectory rates agree at the lowest energies but differ markedly at the highest energies. Fits of trajectory data to the RRK expression frequently yield s values less than the theoretical 3N26. Trajectory calculations are often undertaken at very high energies to make the computations feasible. At these energies, anharmonicity, which is not considered in the RRK theory, has a great influence on the rates. In the present case, anharmonicity seems to have negligible effects at the energies selected in this study because the k stat values, calculated with an anharmonic potential, fit Eq. ~25! with a value of s 512. The discrepancy between the trajectory and RRK behaviors may be due to the breakdown of other RRK assumption: fast IVR on the time scale of reaction. The ultrafast deep UV mass-resolved photoionization spectroscopy investigation undertaken by Owrutsky and

11329

Baronavski12 reveals that, under their experimental conditions, the acetyl dissociation occurs with a rate coefficient of 0.28–0.38 ps21, which corresponds to an acetyl internal energy, based on RRKM calculations, of about 25 kcal/mol. These authors also suggested that the overall dissociation dynamics of acetone may be described in terms of a fully statistical dissociation mechanism. Our dynamical investigation suggests that at low internal energies the dissociation of the acetyl radical is statistical, corroborating the experimental observations of Owrutsky and Baronavski.12 We have used the RRK theory, with the frequency factor determined from the trajectory data, to calculate the rate coefficient at 25 kcal/mol ~51.6 kcal/mol by inclusion of ZPE!. The theoretical RRK value thus obtained is 0.21 ps21, being in close agreement with the experimental value. At energies close to the threshold ~17–23 kcal/mol!, Shibata et al.14 obtained experimental rate coefficients in the range 0.01–0.10 ps21, which are an order of magnitude smaller than those predicted by RRKM theory. Again, using the RRK extrapolation we obtained 0.06 and 0.16 ps21 for 17 and 23 kcal/mol ~above ZPE!, respectively. These values are somewhat higher than those observed experimentally. Note that at energies near the threshold the ZPE problem inherent to classical calculations may aphysically affect the dynamics. However, as Shibata et al.14 concluded, the disagreement between the experimental and theoretical values may arise from nonstatistical effects. In particular, the initial internal energy may be partitioned to vibrational modes that are poorly coupled with the reaction coordinate, leading to a decrease in the rate coefficient. Our calculations also suggest that apparent non-RRKM behavior may be exhibited in their experiment.14

2. Normal mode excitations

In this section, we explore the possibility of modespecific effects in the dissociation of the acetyl radical. To this end, we evaluated the trajectory rates resulting from initial excitation of each normal mode at 65.9 kcal/mol. In what follows, each ensemble will be distinguished by the frequency ~in cm21! of the normal mode at the reactant geometry that is initially excited. In addition, normal modes will be distinguished by its frequency at the reactant geometry, too, even though the frequencies can change substantially in the course of the reaction. Most of the decay curves obtained under normal mode excitations exhibited nonlinear behavior. For this reason, the unreacted trajectories at time t, N(t), were fitted to N ~ t ! 5N ~ 0 ! fast exp~ 2k fastt ! 1N ~ 0 ! slow exp~ 2k slowt ! .

~26!

Here we have assumed that the trajectories behave according to a double-exponential decay, being N(0) fast and N(0) slow the total number of trajectories that dissociates with a fast (k fast) or a slow (k slow) rate coefficient, respectively. For the sake of example, the fit concerning mode-excited 507 ~CCO bend! is shown graphically in Fig. 4. Table V lists the two sets of rates obtained from these fits. It is clearly seen that ensembles 507, 1079, 1504, and 1939 show strong mode-


11330



FIG. 4. Number of unreacted trajectories as a function of time for ensemble 507. The solid line represents the biexponential fit of the trajectory results.

specific effects. For these ensembles, the fast rate coefficient is much higher than those obtained for the remaining ensembles. To explain the above mode-specificity, we have investigated the energy flow in the acetyl radical using the normal mode approach. It is well known that this approach fails as the simulation time increases because the total harmonic energy quickly exceeds the actual energy of the molecule. However, it is a simple way to explore, from a qualitative viewpoint, the pattern of intramolecular energy flow for a system of this complexity. In this analysis, we employed the local mode approach to monitor the energy changes in the CH3 stretches. The local and normal mode energies were averaged over 200 trajectories. For ensemble 116 ~torsional mode!, these energies are shown in Fig. 5, plots ~a! and ~b!, respectively. In a very short time scale modes 1504 (CH3 umbrella vibration! and 1939 ~CO stretch! gain a great amount of energy, but also modes 1477, 1476, 975, and 958

TABLE V. Trajectory rate coefficients for normal-mode excitations. Ensemblea 116 507 958 975 1079 1476 1477 1504 1939 2799 2908 2909 Microcanonical a

(k fast /ps21) b

(k slow /ps21) b

0.960~79! 4.344~82! 1.731~17! 1.644~20! 2.962~87! 1.764~19! 1.776~18! 3.347~85! 14.902~90! 1.867~92! 1.855~92! 1.870~7!

0.959~21! 1.047~18! 1.714~83! 1.640~80! 0.074~13! 1.758~81! 1.769~82! 0.204~15! 0.097~10! 1.863~8! 1.854~8! 1.867~93! 1.372

FIG. 5. CH local mode energies @plot ~a!#, and normal mode energies @plot ~b!#, versus time averaged over 200 trajectories for selective excitation of mode 116.

are excited after 0.025 ps, approximately. In particular, the CO stretching mode may gain energy due to the strong kinetic coupling between the torsion and the terminal bond, which arises from variations in the torsional moment of inertia when the CO distance changes. On the other hand, the CH bonds seem to lose energy in the simulation period of time. The same sort of analysis is shown in Fig. 6 for ensemble 507. In this case the pattern of energy flow is very different. In the first stages of the dynamics, the energy is shared among a small subset of vibrational normal modes. The CCO bending mode energy leaks out to modes 1939, 1079, 1504, and 958, which, except the latter, show specificity. Coupling with modes 958 ~symmetric rock! and 1504 (CH3 umbrella vibration! might be due to 2:1 and 3:1 resonances, respectively. On the other hand, mode 1079 ~closely related to the reaction coordinate! is highly excited; in this case, the switching function SW~2! and the interaction term given by Eq. ~15! may drive the energy transfer between the two modes ~CO stretch and CCO bend!. Methyl group modes ~except 1504 and 958! did not gain energy upon CCO bend excitation, suggesting a poor coupling with the CCO bending motion. All these features explain more clearly the specificity of mode 507. Figure 7 depicts the local and normal mode energies as a function of time for ensemble 958 ~symmetric rock!. Although modes 1079, 1939, and 1504 gain energy, the time scale for this energy flow is significantly longer than that in ensemble 507. Moreover, the energy seems to be redistributed among a large number of normal modes, as opposed to the behavior found for ensemble 507. In fact, ensemble 958 does not present specific chemistry.

Each ensemble contains 65.9 kcal/mol into one of the normal modes of the molecule. The ensemble is designated by the frequency ~in cm21! of the normal mode that is being excited ~see text!. b The percentage of trajectories decaying with a certain rate coefficient is given in parentheses. For the microcanonical ensemble, k traj is tabulated. This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP:

193.144.81.13 On: Fri, 21 Nov 2014 07:50:10



11331



Figures 8 and 9 display the normal mode energies for mode-excited 1079 and 1939, respectively. The patterns of energy flow for these two normal modes are very similar: in a short time scale the energy is redistributed among normal

modes 1504, 1939, and 1476, for ensemble 1079, and among modes 1504, 1079, and 958, for ensemble 1939. Therefore, we suggest a strong coupling among modes 1939, 1079, and 1504. Mode 1504 corresponds to the umbrella vibration at

FIG. 9. CH local mode energies @plot ~a!#, and normal mode energies @plot FIG. 7. CH local mode energies @plot ~a!#, and normal mode energies @plot ~b!#, versus time averaged over 200 trajectories for selective excitation of ~b!#, versus time averaged over 200 trajectories for selective excitation of mode 1939. mode 958. This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP:

193.144.81.13 On: Fri, 21 Nov 2014 07:50:10

11332



TABLE VI. Trajectory rate coefficientsa for the dissociation of CH3CO and CD3CO. E/~kcal/mol! b

28.4 38.4b 65.9b 65.9c

CH3CO

CD3CO

0.298 0.640 1.372 2.970

0.263 0.535 1.427 2.248

In ps21. With EMS initialization. c With the energy deposited initially in the CCO moiety ~see text!. a

b

served apparent non-RRKM behavior in the dissociation of the acetyl radical at 65610 kcal/mol, with the energy preferentially localized in the CCO moiety. However, they suggest a poor coupling between modes involving the CCO moiety and the C–C stretch reaction coordinate, which is contrary to the results obtained in this study. 3. Kinetic isotope effects


the reactant but as the reaction takes place this movement is transformed into the wagging motion of the CH3 product. These two coordinates ~CCH bend and wag! are coupled to the CC stretch by potential coupling terms @Eqs. ~9! and ~12!#. For ensemble 1079, it can be seen how the CH stretches as well as other modes gain energy as the CC stretch relaxes. This ensemble redistributes the energy, within the same time scale, among more degrees of freedom than ensemble 1939 does. This may cause the rate for ensemble 1079 to be lower than that obtained for ensemble 1939, which is the highest rate. In summary, for ensembles 1079 and 1939, the energy is quite restricted to a small part of phase space that is coupled to the reaction coordinate, resulting in rate enhancement. Finally, Fig. 10 shows the normal mode and local mode energies against time for ensemble 2799 (CH3 symmetric stretch!. In plot 10~a! we display the local mode energies of the three CH bonds as a function of time. In the period of time studied ~0.3 ps!, there is an almost complete, irreversible, energy relaxation of the CH stretching vibrations. There is an energy flow into other modes, which occurs after the first 0.15 ps ~except for mode 1504!. Although modes 1939 and 1504 gain sufficient energy, the time scale for the energy redistribution is much longer than that involved in excitations of modes that present high specificity. In summary, the redistribution of energy for those modes that present strong mode specific effects, most of them involving the CCO arrangement, is restricted to a small subset of modes in a short time scale, which are strongly coupled to the reaction coordinate. This dynamical picture leads to an incomplete IVR and an enhancement of the reaction rate. As remarked in the Introduction, Zewail and co-workers13 ob-

To analyze the influence of deuterium substitution on the rates of reaction, we have computed trajectory rate coefficients, under EMS initialization, for the trideuterated CD3CO radical and evaluated the ratio k H /k D ~where k H and k D are the rate coefficients for the undeuterated and the deuterated molecules, respectively! at energies 28.4, 38.4, and 65.9 kcal/mol. The trajectory rate coefficients for the two species are listed in Table VI. As can be seen, a normal isotope effect is found at the two lowest energies, but a slight inverse isotope effect arises at 65.9 kcal/mol. In particular, the trajectory ratios are 1.13, 1.20, and 0.96 at 28.4, 38.4, and 65.9 kcal/mol, respectively, as shown in Table VII. For the sake of comparison, we also calculated the values of k H /k D by the RRKM theory. The classical harmonic RRKM rate coefficient ~RRK! is given by27 k

F G

P si51 n i E2E 0 5 s21 E P j51 n j

RRK

s21

~27!

,

where the n i andn j are the reactant and transition state harmonic vibrational frequencies, respectively. The remaining parameters appearing in this expression were defined above. On the other hand, the quantum harmonic RRKM rate coefficient for a nonrotating molecule with an internal energy E may be expressed as k RRKM5

W ts~ E * ! , hr~ E !

~28!

where r (E) is the density of states of the reactant molecule, W ts(E * ) is the total number of states on the transition state surface ~with energy E * 5E2E 0 ),15~c!,27 and h is Planck TABLE VII. Kinetic isotope effects found in this study. E/~kcal/mol!

(k H /k D) traj

(k H /k D) RRK

(k H /k D) RRKM

28.4 38.4 65.9

1.13 1.20 0.96, 1.32a

0.68 0.68 0.68

1.54 1.27 1.06

For the ensemble in which only the CCO moiety is excited ~see text!. This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 193.144.81.13 On: Fri, 21 Nov 2014 07:50:10 a



11333

IV. CONCLUSIONS

FIG. 11. Ratio between rate coefficients for CH3CO and CD3CO as a function of energy predicted by RRKM calculations.

constant. Both the classical and quantum harmonic RRKM rate coefficients for CH3CO and CD3CO were calculated in this work. For the latter calculations, the methyl torsion was modeled by a free-rotor at the transition state.15~c! This seems to be a good approximation since the frequency of the methyl torsion at the transition state is only 44 cm21. The RRKM calculations reported in this paper were accomplished using the UNIMOL suite of programs.28 Table VII lists the RRK and RRKM predictions, also for k H /k D at the three energies selected in the trajectory calculations. In addition, the RRKM ratio for the energy range 25–130 kcal/mol is shown graphically in Fig. 11. For energies above '43 kcal/mol, we found previously that the system exhibits intrinsic non-RRKM behavior; therefore, the results of these statistical formulations should be taken with caution. The RRK theory predicts a marked inverse isotope effect that contrasts with the trajectory and RRKM results. The ratio predicted by the RRKM theory decreases significantly with energy, leading to an inverse isotope effect at energies higher than 95 kcal/mol. This trend results from the presence of a free-rotor at the transition state, which provides a much larger sum-of-states contribution in CD3CO than in CH3CO. Overall, our results suggest that at high energies the acetyl dissociation might exhibit an inverse isotope effect. Zewail and co-workers13 investigated the influence of deuterium substitution on the rates of the acetyl dissociation. In particular, they found a ratio k H /k D of 1.67 at internal energies of 65610 kcal/mol. This value cannot be directly compared with those presented in Table VII because in the experimental study the energy distribution in the acetyl intermediate is nonrandom. For this reason, we employed the classical trajectory method to calculate k H /k D for an ensemble where the CCO bending and the CC and CO stretching modes share 65.9 kcal/mol. We fitted the decaying trajectories, for both the undeuterated and the trideuterated acetyl, to a single-exponential to obtain a rough estimate of k H /k D . These computations predicted a normal isotope effect with a ratio of 1.32, which agrees, at least qualitatively, with the experimental findings.

Classical trajectory and statistical EMS–TST calculations were used to study the dynamics of the acetyl radical dissociation. The computations were performed with a potential energy surface that reproduces reasonably well the main characteristics of the ab initio ground-state surface. Under random initial conditions, the calculations predict statistical dynamics at the lowest energies selected ~specifically up to '41 kcal/mol!, which is in agreement with the experimental observations of Owrutsky and Baronavski.12 However, as energy increases, the dissociation dynamics becomes intrinsically non-RRKM. We have also carried out classical trajectory calculations under initial normal mode excitations to analyze the possibility of mode specificity. The results suggest that excitations of normal modes 507 ~CCO bend!, 1079 ~CC stretch!, 1504 (CH3 umbrella vibration!, and 1939 ~CO stretch! enhance significantly the rate of reaction; in particular, excitation of the CO stretch gives a rate coefficient an order of magnitude higher than that obtained under EMS initialization. These mode specific effects are explained in terms of a restricted IVR. Our results contrast with those of Zewail and co-workers,13 which indicate weak couplings between modes involving the CCO moiety and the reaction coordinate. Our analysis of the isotopic substitution suggests that under random initial conditions and at high energies the system may show an inverse secondary isotope effect. On the other hand, under initial excitation of modes involving the CCO moiety, the classical trajectories performed at 65.9 kcal/mol predict a normal isotope effect (k H /k D51.32), in agreement with the experimental results obtained by Zewail’s group.13

ACKNOWLEDGMENTS

A.P.-G. and E.M.-N. thank Fundacioń Gil Da´vila and Xunta de Galicia, respectively, for their grants. We thank R. Rodrı´guez-Fernańdez for his help as system manager of our computing facilities.

D. L. Bunker, Theory of Elementary Gas Reaction Rates ~Pergamon, London, 1966!. 2 ~a! D. L. Bunker and W. L. Hase, J. Chem. Phys. 59, 4621 ~1973!; ~b! W. L. Hase, J. Phys. Chem. 90, 365 ~1986!. 3 S. North, D. A. Blanck, and Y. T. Lee, Chem. Phys. Lett. 224, 381 ~1994!. 4 S. W. North, D. A. Blanck, J. D. Gezelter, C. A. Longfelow, and Y. T. Lee, J. Chem. Phys. 102, 4447 ~1995!. 5 S. W. North, Ph.D. thesis, University of California, Berkeley, 1995. 6 P. M. Kroger and S. J. Riley, J. Chem. Phys. 67, 4483 ~1977!. 7 G. E. Hall, D. Bont, and T. J. Sears, J. Chem. Phys. 94, 4182 ~1991!. 8 S. K. Kim, S. Pedersen, and A. H. Zewail, J. Chem. Phys. 102, 477 ~1995!. 9 S. A. Buzza, E. M. Snyder, and A. W. Castleman, Jr., J. Chem. Phys. 104, 5040 ~1996!. 10 S. K. Kim and A. H. Zewail, Chem. Phys. Lett. 250, 279 ~1996!. 11 K. A. Trentleman, S. H. Kable, D. B. Moss, and P. L. Houston, J. Chem. Phys. 91, 7498 ~1989!. 12 J. C. Owrutsky and A. P. Baronavski, J. Chem. Phys. 108, 6652 ~1998!. 13 S. K. Kim, J. Guo, J. S. Baskin, and A. H. Zewail, J. Phys. Chem. 100, 9202 ~1996!. This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 1

193.144.81.13 On: Fri, 21 Nov 2014 07:50:10

11334 14



T. Shibata, H. Li, H. Katayanagi, T. Suzuki, J. Phys. Chem. A 102, 3643 ~1998!. 15 ~a! Advances in Classical Trajectory Methods, edited by W. Hase ~JAI, London, 1992!, Vol. 1; ~b! T. Baer and W. L. Hase, Unimolecular Reaction Dynamics: Theory and Experiments ~Oxford University Press, Oxford, 1996!; ~c! R. G. Gilbert and S. C. Smith, Theory of Unimolecular and Recombination Reactions ~Blackwell Scientific, Oxford, 1990!; ~d! L. M. Raff and D. L. Thompson, in Theory of Chemical Reaction Dynamics, edited by M. Baer ~Chemical Rubber, Boca Raton, 1985!. 16 ~a! H. W. Schranz, L. M. Raff and D. L. Thompson, J. Chem. Phys. 94, 4219 ~1991!; ~b! H. W. Schranz, L. M. Raff, and D. L. Thompson, Chem. Phys. Lett. 68, 171 ~1990!. 17 ~a! N. E. Triggs, M. Zahedi, J. W. Nibler, D. DeBarber, and J. J. Valentini, J. Chem. Phys. 96, 1822 ~1992!; ~b! G. Herzberg, Electronic Spectra and Electronic Structure of Polyatomic Molecules ~D. Van Nostrand, Princeton, 1969!. 18 J. W. G. Watson, J. Mol. Spectrosc. 45, 530 ~1973!.

19

S. Deshmukh, J. D. Myers, S. S. Xantheas, and W. P. Hess, J. Phys. Chem. 98, 12535 ~1994!. 20 D. L. Osborn, H. Choi, D. H. Mourdant, R. T. Bise, D. M. Neumark, and C. M. Rohlfing, J. Chem. Phys. 106, 3049 ~1997!. 21 D. L. Thompson, GENDYN program, Stillwater, 1991. 22 ~a! H. W. Schranz, S. Nordholm, and G. Nyman, J. Chem. Phys. 94, 1487 ~1991!; ~b! G. Nyman, S. Nordholm, and H. W. Schranz, J. Chem. Phys. 93, 6767 ~1990!. 23 J. D. Doll, J. Chem. Phys. 73, 2760 ~1980!; 74, 1074 ~1981!. 24 R. D. Kay, L. M. Raff, J. Phys. Chem. A 101, 1007 ~1997!. 25 P. Pechukas, Annu. Rev. Phys. Chem. 32, 159 ~1981!. 26 D. V. Shalashilin and D. L. Thompson, J. Chem. Phys. 105, 1833 ~1996!; 107, 6204 ~1997!. 27 P. J. Robinson and K. A. Holbrook, Unimolecular Reactions ~Wiley, New York, 1972!. 28 R. G. Gilbert, M. J. T. Jordan, and S. C. Smith, UNIMOL program suite, Sydney, 1992.
