4418
J. Phys. Chem. A 2010, 114, 4418–4426
Reaction-Path Dynamics Calculations of the Cl + NH3 Hydrogen Abstraction Reaction: The Role of the Intermediate Complexes M. Monge-Palacios and J. Espinosa-Garcia* Departamento de Quı´mica Fı´sica, UniVersidad de Extremadura, 06071 Badajoz, Spain ReceiVed: December 9, 2009; ReVised Manuscript ReceiVed: February 19, 2010
Using ab initio information at the CCSD(T)/cc-pVTZ level, the reaction path for the Cl + NH3 hydrogen abstraction reaction was traced, and the coupling terms between the reaction coordinate and the normal modes were analyzed along it. Two intermediate complexes were located in the entry channel and characterized close to the reactants. One of them presents a typical Cl · · · H-N bond, while the second presents a twocenter/three-electron Cl∴N bond. Both complexes are on the reaction path and contribute to the final rate constants. With this information, the rate constants were calculated over the temperature range 200-2000 K, using variational transition state theory with multidimensional tunneling contributions, and were found to reproduce the experimental evidence in the common temperature range. Finally, analysis of the coupling terms showed qualitatively that vibrational excitation of the N-H stretch and the bending and umbrella modes in the reactant NH3 enhances the forward thermal rate constants, and that, in the products, the H-Cl stretch mode and the bending mode in NH2 could appear vibrationally excited, although the randomization of the energy in the well in the exit channel might diminish this excitation. 1. Introduction The chemistry of the reaction of ammonia with chlorine atoms is very complex with many intermediate fast reactions being involved.1 In some cases, neither the products nor the rate constants are clearly known. Although quantitatively less important than the reactions with H, O, and OH, this reaction also plays an important role in the conversion of fuel-nitrogen to the atmospheric contaminant NO.2 Experimentally, there have been very few studies of the title reaction. To the best of our knowledge, only the following two kinetics measurements have been reported. Using the flash photolysis resonance fluorescence technique, Westenberg and DeHass3 in 1977 determined the rate constant at room temperature, k(298 K) ) 1.23 × 10-13 cm3 molecule-1 s-1. Recently, Gao et al.4 using the time-resolved resonance fluoresce technique presented the first determination of the temperature dependence in the range 290-570 K for this endothermic reaction. Their fit to the forward rate constant was k(T) ) (1.08 ( 0.05) × 10-11 exp(-2.74 ( 0.04 kcal mol-1/RT) cm3 molecule-1 s-1, with a value of 1.06 × 10-13 at room temperature, in close agreement with the previous measurement. Theoretically, this reaction has also received little attention. We know of only three theoretical studies. Kondo et al.5 performed ab initio calculations using the Gaussian-2 (G2) theory, focusing their attention on the geometry and energy of the stationary points. They found an enthalpy of activation at 0 K, [∆H‡(0 K)], of 4.75 kcal mol-1. Later, in 2006, Gao et al.4 constructed the reaction path using the MPWB1K density functional (DFT) method,6 finding ∆H‡(0 K) ) 2.1 kcal mol-1, although using different correlation energy levels and basis sets they report values in the range 1.7-3.1 kcal mol-1, which in any case are considerably lower than the previous value. The calculated rate constant is around a factor of 3 too large compared with experiment, reflecting an error in the width and height of the barrier. Finally, Xu and Lin7 performed a * To whom correspondence should be addressed.
computational study of the mechanisms and kinetics of the reaction. The geometry of the stationary points were optimized using the B3LYP DFT method, and their energies were refined with the modified Gaussian-2 (G2M) theory. They obtained ∆H‡(0 K) ) 4.2 kcal mol-1, in contrast with the preceding values. The predicted rate constants were basically in agreement with available experimental data, although they estimated the tunneling effect using the simple Eckart model, in which it has a negligible influence, 1.2% at 300 K. When they calculated this effect with a more sophisticated model, the small-curvature approach,8 the tunneling factor was 3.99 and 1.62 at 300 and 500 K, respectively, clearly worsening the agreement with experiment. There is therefore notable disagreement between the results of these three theoretical studies with respect to their calculations of the barrier heightsa parameter of crucial importance in the kinetics analysis, and which is clearly strongly dependent on the level of the correlation energy and the basis set. Moreover, in their analysis of the reaction mechanism, these two last groups4,7 found different complexes in the entry and exit channels. In the present work, we perform an exhaustive theoretical study of the title reaction using a higher ab initio level, with the main aim of accurately defining the barrier height and the stability of the complexes. We first locate and characterize the stationary points (reactants, products, saddle point, and complexes in the entry and exit channels) by using high-level electronic structure calculations. Next, we calculate the reaction path with the “direct dynamics method”,9 which describes a chemical reaction by using ab initio information (energies, gradients, and Hessians) only in the region of configuration space along the reaction path. To obtain kinetics information, we perform variational transition-state theory (VTST) calculations with multidimensional tunneling corrections. The theoretical results are compared with the available theoretical and experimental evidence.
10.1021/jp911664t 2010 American Chemical Society Published on Web 03/05/2010
Cl + NH3 Hydrogen Abstraction Reaction
J. Phys. Chem. A, Vol. 114, No. 12, 2010 4419
II. Methods and Computational Details 1. Ab Initio Calculations. Electronic structure calculations were carried out using the GAUSSIAN-03 set of programs.10 The stationary point geometrical parameters were fully optimized and characterized at the single and double coupled cluster theory level with the inclusion of a perturbational estimate for triple excitations,11 CCSD(T), using the cc-pVTZ basis set.12 We denote this the CCSD(T)/cc-pVTZ level. At this level, we constructed the minimum energy path (MEP) over the broad range of the reaction coordinate, s, -1.0 to +1.0 bohr, where s ) 0 corresponds to the location of the saddle point. With these values of s, we are covering the most relevant portion of the reaction path. In the rest of paper, the units of s are bohr, and all calculations are carried out in mass-scaled coordinates with a reduced mass µ equal to 1 amu. Thus, distances through the mass-scaled coordinates in bohr are equivalent to distances through mass-weighted coordinates in bohr amu1/2. 2. Direct Dynamics Calculations. We performed a generalized normal-mode analysis projecting out frequencies at each point along the MEP.13 With this information, as a first step we calculated the ground-state vibrationally adiabatic potential curve G VGa (s) ) VMEP(s) + εint (s)
(1)
where VMEP(s) is the classical energy along the MEP with its zero at the reactants (s ) -∞), and εintG(s) is the zero-point energy at s from the generalized normal-mode vibrations orthogonal to the reaction coordinate. As a second step, the coupling terms, Bk,F(s), were computed. These terms quantify the coupling between the normal mode, k, and the motion along the reaction coordinate, mode F. They control the nonadiabatic flow of energy,14,15 allowing one to calculate accurate semiclassical tunneling factors, that is, dynamical features. The coupling terms are also components of the reaction path curvature, κ(s), defined as
κ(s) ) (Σ[Bk,F(s)] )
2 1/2
∫0∞ k(E)P(E) dE
Qe ) 4 + 2exp(-ε/kBT)
(4)
assuming that the electronic partition function of the transition state is 2, that is, it is assumed that the spin-orbit effect is fully quenched in this zone. Note that the spin-orbit (s-o) effect is more important than in the similar F + NH3 reaction,21 where the excitation energy20 is only 404 cm-1 ≈ 1.1 kcal mol-1. III. Results and Discussion 1. Structures, Vibrational Frequencies, and Energies of the Stationary Points. The reaction of the Cl(2P) atom with ammonia proceeds through a saddle point with complexes in the entry and exit channels. The optimized geometries and harmonic vibrational frequencies at the CCSD(T)/cc-pVTZ level are listed in Table 1 for the reactants and products, Table 2 for the saddle point, and Table 3 for the possible complexes. Figure 1 shows the optimized geometries of all the stationary points at this level of theory. For the reactants and products, the CCSD(T)/cc-pVTZ level reproduces the experimental geometries and vibrational frequencies. Although a direct comparison of the standard enthalpy of reaction (298 K) with experiment is available, this magnitude presents two serious problems which we have to address separately: (i) the accuracy of the experimental enthalpy of reaction, and (ii) s-o coupling. (i) The enthalpy of reaction can be obtained from the reactant and product enthalpies of formation at 298 K, ∆H°, f as 298 ∆H298 ) ∆H298 r f [HCl(g)] + ∆Hf [NH2(g)] -
(2)
Then, the energies, vibrational frequencies and geometries along the MEP were used to estimate the rate constants by using variational transition state theory (VTST). The thermal rates were calculated using the microcanonical variational theory8,16 (µVT) approach, which locates the dividing surface between reactants and products at a point s*, µVT(T) along the reaction path that minimizes the sum of vibrational-rotational states, µVT µVT GT (E), at total energy E. Thus, Nvr (E) ) minsNvr (E,s). The Nvr microcanonical thermal rate constant, k(T), is obtained from the rate constant k(E) at total energy E, taking into account the Boltzmann distribution, P(E),
k(T) ) σ
GAUSSIAN03 set of programs,10 and POLYRATE.18 The rotational partition functions were calculated classically, and the vibrational modes were treated as quantum-mechanical separable harmonic oscillators. In calculating electron partition functions, we included the 2P1/2 excited state of Cl (with an excitation energy20 of ε ) 882 cm-1 ≈ 2.5 kcal mol-1) in the reactant electronic partition function,
(3)
where σ is the reaction degeneracy (the number of equivalent reaction paths, which were taken to be 3 and 1 for the forward and reverse reactions, respectively). As a last step, we considered the tunneling contribution. Since there is only information on the reaction-path, the centrifugal-dominant small-curvature tunneling (SCT)17,18 approximation was used. All kinetics calculations were performed using the GAUSSRATE code,19 which serves as an interface between the
298 ∆H298 f [Cl(g)] - ∆Hf [NH3(g)]
(5)
While the ∆H°f values of NH3, Cl, and HCl are well established,20 the proposed experimental values of the NH2 free radical differ. The JANAF tables20 recommend 45.5 ( 1.5 kcal mol-1, which agrees with a recent experimental determination.22 However, recent theoretical23 and experimental24 studies diminish this value by one unit, 44.6 ( 0.1 and 44.5 ( 0.1 kcal mol-1, respectively, in agreement with our theoretical value,25 43.8 ( 0.6 kcal mol-1, within the uncertainties. That the enthalpy of formation of the NH2 product seems still an open question doubtless strongly influences the standard enthalpy of reaction, giving values between 3.7 and 5.4 kcal mol-1, with a proposed value of 4.4 kcal mol-1 assuming the latest experimental value.24 (ii) With respect to s-o coupling, the chlorine atom has two low-lying fine structure electronic states, 2P1/2 and 2P3/2, with a separation20 of ε ) 882 cm-1 ≈ 2.5 kcal mol-1. A relativistic study is beyond the scope of this work, and to compare our theoretical results with experiment we corrected the experimental value for the s-o effect, reducing it by one-third of the split between the two states, that is, 0.8 kcal mol-1. Thus, the corrected experimental heat of reaction is 3.6 kcal mol-1, in which it has been assumed that s-o effects are completely quenched in molecules. This is the “non-relativistic experimental” value to be compared with the theoretical values in this
4420
J. Phys. Chem. A, Vol. 114, No. 12, 2010
Monge-Palacios and Espinosa-Garcia
TABLE 1: Reactant and Product Propertiesa Cl Theo.
b
NH3 Exp.
c
Theo.
geometry X-H
