J. Phys. Chem. A 1998, 102, 10715-10722
10715
Reaction-Path and Dual-Level Dynamics Calculations of the CH3F + OH Reaction J. Espinosa-Garcı´a,*,† E. L. Coitin˜ o,‡ A. Gonza´ lez-Lafont,§ and J. M. Lluch§ Departamento de Quı´mica Fı´sica, UniVersidad de Extremadura, 06071 Badajoz, Spain, Laboratorio de Quı´mica Teo´ rica y Computacional, Facultad de Ciencias, UniVersidad de la Repu´ blica, MonteVideo, 11200, Uruguay, and Departament de Quı´mica, UniVersitat Autonoma de Barcelona, 08193 Bellaterra, Barcelona, Spain ReceiVed: July 29, 1998; In Final Form: October 16, 1998
Using two approaches, reaction-path dynamics (which uses information on electronic structure energy and energy derivatives calculated ab initio along the minimum energy path) and dual-level dynamics (in which the reaction-path is calculated at a low level of theory, and stationary point information from a high level of theory is used to interpolate corrections to energy quantities, vibrational frequencies, and moments of inertia), the reaction-path for the CH3F + OH reaction was traced. Qualitatively, both methods yield similar results, showing the usefulness of the more economical dual-level approach. The H2O product may be vibrationally excited due to the nonadiabatic flow of energy between the reaction coordinate and the O-H bound mode. The rate coefficents were calculated for the temperature range 200-1000 K using the variational transitionstate theory. The tunneling effect was included in the first approach using the small-curvature tunneling method, and in the second using the microcanonical optimized multidimensional tunneling method. Tunneling plays an important role in this reaction and the calculated rate coefficients show a more pronounced curvature in the Arrhenius plot than the experimental data.
1. Introduction The study of the kinetics of halocarbons with the hydroxyl radical is very important for the knowledge of their combustion and of their role in atmospheric processes, particularly in the chemistry of stratospheric ozone. Experimentally, three techniques have been employed to study the rate coefficient of the reaction of fluoromethane (CH3F) with the hydroxyl radical: discharge flow system,1,2 flash-photolysis methods,3,4 and relative rate measurements, which involve measurement of the fractional loss of the reactant compound compared to a reference compound, in the presence of OH.5,6 In general, the flow measurements yield values lower than those obtained by flashphotolysis methods, and the differences are greater in the lower temperature range. Unfortunately, there have been only a few studies of this reaction’s temperature dependence,2,4,5 and then only over a small temperature range (240-480 K). The first study was performed by Jeong and Kaufman2 in the range 292-480 K, finding a small curvature in the Arrhenius plot. In 1985, Atkinson7 compiled the earlier rate coefficients and recommended the expression
k ) 5.51 × 10-18 T2 exp(-1005 ( 168/T) cm3 molecule-1 s-1 with an estimated uncertainty of (30% at 298 K. Schmoltner et al. 4 studied this reaction in the range 243-373 K, and several different OH sources were used: H2O, HONO, H2O2, and HNO3. They found that the values at 298 K and above were essentially independent of the OH source; however, significant * Author to whom correspondence should be addressed. † Universidad de Extremadura. ‡ Universidad de la Repu ´ blica. § Universitat Autonoma de Barcelona.
differences were observed between the results obtained using different sources at lower temperatures. The data were presented using an Arrhenius plot, and the plotted curve exhibited no curvature. Finally, Hsu and DeMore5 studied the relative rates in the range 283-403 K, and again no curvature was found in the Arrhenius plot. DeMore6 fits the results to
k ) 4.4 × 10-12 exp (-1655/T) cm3 molecule-1 s-1 with an estimated uncertainty of 10% at 298 K. The above experimental disagreement concerning the curvature of the Arrhenius plot is probably simply a result of the small temperature range studied (240-480 K). However, this reaction presents a heavy-light-heavy (HLH) mass combination, which is a good candidate to present a large tunneling effect at low temperatures and, therefore, curvature in the Arrhenius plot. Obviously, this effect will become more patent when a larger range of temperatures is analyzed. The large amount of experimental work contrasts with the paucity of theoretical studies, which have been limited to ratecoefficient calculations using conventional transition-state theory (TST). Thus, Jeong and Kaufman8 carried out calculations using the conventional TST, and estimated the tunneling factors by the Wigner and the Eckart methods, which are known to underestimate this effect. Cohen and Benson9 also used TST, but did not include the tunneling effect in the calculation. To perform the dynamic and kinetic study of the title reaction, we decided to first compare two different approaches to calculating the necessary initial information concerning the potential energy surface (PES). In the first approach, we used ab initio electronic structure information only in the region along the reaction-path,10-12 which is usually referred to as a “direct dynamics” method,11 and hereafter we shall call it the “reaction-path approach”. This
10.1021/jp9832138 CCC: $15.00 © 1998 American Chemical Society Published on Web 12/04/1998
10716 J. Phys. Chem. A, Vol. 102, No. 52, 1998 method describes a chemical reaction by using ab initio information (energies, gradients, and Hessians) only in the region of configuration space along the reaction-path and has presented satisfactory results for some hydrogen abstraction reactions.10,11,13-19 However, it has two major limitations: high computational cost, and the fact that the influence of the PESs large curvature on the tunneling effect cannot be taken into account. This last effect is especially important in reactions with a heavy-light-heavy (HLH) mass combination, such as the title reaction. In the second approach, we used the “dual-level” methodology.20-22 A “higher-level” is used at a few selected points, normally just the stationary points: reactants, products, and saddle point; and a “lower-level” is used at a large number of geometries to generate a wider PES representation. This permits the description of the wide reaction-swath12 of geometries accessed by tunneling paths when the reaction-path curvature is large. To reduce the computational cost, we use semiempirical molecular orbital theory as the “lower-level”. The dynamical calculations are carried out at this “lower-level”, but the energies and frequencies are corrected by interpolated corrections (IC)20 in order to obtain “higher-level” accuracy. To obtain kinetic information we perform variational transition-state theory (VTST) calculations with the inclusion of multidimensional tunneling effects. 2. Methods 2.1. Dual-Level Approach. Lower-LeVel Electronic Structure Calculations. The semiempirical molecular orbital theory has been traditionally used as an economical alternative for theoretical electronic structure calculations. In this paper, we used the PM3 method.23 This economical method permits the description of the wide reaction-swath region12 of geometries accessed by tunneling paths when the reaction-path curvature is large, as in this case. To achieve a better representation of critical regions of the PES, in particular of the geometries and relative energies of the stationary points, an alternative is to develop a set of specific reaction parameters (SRP) for the semiempirical Hamiltonian.24-26 In this reaction, however, the original parameters in the PM3 method already gave an adequate initial description of these stationary points (geometry and energy). Successive trials did not lead to any improvement. Thus, in this particular case, possibly due to a fortunate error cancellation, the original PM3 provides an adequate initial description as lower-level. Higher-LeVel Electronic Structure Calculations. Geometries, energies, and first and second energy derivatives were calculated using the GAUSSIAN 94 system of programs.27 In the first step, we calculated geometries, energies, and vibrational frequencies at all stationary points in second-order Møller-Plesset perturbation theory28 (MP2) with full electron correlation. For the closedshell systems, the MP2 wave function was based on the restricted Hartree-Fock (RHF) wave function and the 6-31G(d,p) basis set.29 For the open-shell systems, the MP2 calculations were based on both an unrestricted Hartree-Fock (UHF) wave function30,31 (UMP2, and the 6-31G(d,p) basis set) and a restricted open-shell Hartree-Fock (ROHF) wave function32-34 (ROMP2, and the 6-311G(2d,2p) basis set). To avoid or, at least, to reduce the spin contamination of the UMP2 wave function, we used the spin projection operator as implemented in GAUSSIAN 94. The energy after the spin decontamination will be PUMP2 (projected UMP2). The expectation values of S2 were about 0.75, indicating minor spin contamination. In these conditions, one may expect a similarity in behavior between the UMP2 and ROMP2 wave functions.
Espinosa-Garcı´a et al. In a second step, to improve the energy description of the stationary points in the reaction-path, we made a single-point calculation at higher levels. Level 1. We made a single-point calculation using the MP2/ 6-31G(d,p) optimized geometries with the single- and doublecoupled cluster approach including a perturbative estimate of corrected triple excitations,35 CCSD(T), using the 6-311++G(2d,p) basis set. We denote this energy as
CCSD(T)/6-311++G(2d,p)//(R-U)MP2/6-31G(d,p) where the double slash (X//Y) denotes geometry optimization at level Y and energies calculated at level X. Level 2. Continuing with the CCSD(T) method, we increased the basis set using the augmented correlation-consistent polarized valence triple-ζ sets developed by Dunning et al.,36 AUGcc-pVTZ. We denote this energy as
CCSD(T)/AUG-cc-pVTZ//(R-U)MP2/6-31G(d,p) Level 3. Finally, the scaling all correlation energy method (SAC)37-40 was applied to correct the calculated energies of the reactants, products, and saddle point. In general, the SAC approach consists of scaling all the correlation energies, calculated as the difference between the Hartree-Fock (HF) and the post-Hartree-Fock (i.e. correlated) energies, by using a factor F, calculated as
F)
De(correlated) - De(HF) De(exp) - De(HF)
(1)
where De is the dissociation energy and “exp” denotes experimental. The SAC method has been successfully used to compensate the deficiences in the basis sets and in the treatment of electron correlation for several reactions.41 We denote these methods by the name of the correlated level, for example, the SAC method with MP2 estimates of the correlation energy is called MP2-SAC. One advantage of these methods is that they can be used for the whole potential energy surface.40 Thus, the scaled energy is
ESAC ) EHF +
Ecorrelated - EHF F
(2)
According to the standard notation22 for dual-level dynamics, these calculations are labeled as
high-level//(R-U)MP2/6-31G(d,p)///low-level This triple slash notation (X//Y///Z) denotes reaction-path construction at level Z, stationary point optimization and frequency calculation at level Y, and single-point energies calculated at higher-level X. 2.2. Reaction-Path Approach. Geometries, energies, and first and second energy derivatives were calculated using the GAUSSIAN 94 system of programs.27 The stationary point geometries (reactants, products, and saddle point) are optimized using restricted (R) or unrestricted (U) second-order MøllerPlesset perturbation theory, MP2, with full electron correlation, using the 6-31G(d,p) basis set. Starting from the saddle-point geometry and going downhill to both the asymptotic reactant and product channels in mass-weighted Cartesian coordinates, we constructed the “intrinsic reaction coordinate” (IRC), or minimum energy path (MEP) at this level, with a gradient step size of 0.05 bohr amu1/2, following the algorithm of Gonza´lez and Schlegel.42,43 In addition, for 27 points along the MEP we
CH3F + OH Reaction
J. Phys. Chem. A, Vol. 102, No. 52, 1998 10717
computed the gradients and the Hessians at this level. To improve the energy description of the stationary points on the reaction-path, we used the single-point calculations at the higher levels analyzed in the previous section. 2.3. Dynamics. In the dual-level approach, the reaction-path was computed using the Page-McIver method,44 with a gradient step-size of 0.01 bohr and with the Hessian matrix being recalculated every 2 steps. The reaction-path was calculated from -2.0 bohr on the reactants side to +2.0 bohr on the products side. In this dual-level calculation, the interpolated corrections scheme20 was applied to interpolate corrections to the frequencies at the MP2 level and energies to reproduce the higherlevel behavior. For both types of calculation, reaction-path and dual-level, s denotes the signed distance from the saddle point (s ) 0) along the MEP in an isoinertial mass-scaled coordinate system, with s > 0 referring to the product side. An isoinertial coordinate system is any coordinate system in which the kinetic energy can be written in the form
TABLE 1: Reactant Properties for CH3F + OH f CH2F + H2O CH3F
OH
PM3a
MP2b
ROMP2c
1.092 1.351
1.087 1.388
1.084 1.383
108.60
109.25
109.20
3141 3087 3087 1533 1356 1356 1202 1016 1016
3241 3241 3135 1566 1566 1557 1230 1230 1115
24.0
25.6
PM3
MP2
ROMP2
0.937
0.971
0.964
3200 3200 3107 1535 1535 1531 1223 1223 1091
3987
3844
3839
25.3
5.7
5.5
5.5
geometryd RC-H RC-F RO-H
