Theoretical calculations of the thermal rate constants for the gas

in good agreement with each other in temperature regions where they overlap, but ... neling, and attempts have been made to substantiate this claim by...
0 downloads 0 Views 2MB Size
J . Phys. Chem. 1990, 94, 1096-7106

7096

radiative transport is important at small concentrations of acceptor, as observed before by Birks et and influences the decay curves and fluorescence spectrum of the donor. When the radiative processes are taken into account by the expressions presented, reliable values of the nonradiative energy-transfer rate constant are obtained.

Acknowledgment. This work was supported by the I.N.I.C. (Instituto Nacional de Investiga@o Qentifica) and J.N.I.C.T. (Junta Nacional de InvestigaGZo Cientifica e Tecnolijgica). Financial support from the FundacZo Calouste Gulbenkian for the purchase of a Spex Fluorolog 1 12 spectrofluorometer is gratefully acknowledged.

Theoretical Calculations of the Thermal Rate Constants for the Gas-Phase Chemical Reactions H NH, C H, 4- NH, and D 4- ND, C D, 4- ND,

+

Bruce C. Garrett,*,+ Chemical Dynamics Corporation, Upper Marlboro, Maryland 20772

Michael L. Koszykowski, Carl F. Melius, Sandia National Laboratories, Livermore, California 94550

and Michael Page Naval Research Laboratory, Washington, D.C. 20375 (Received: December 12, 1989; In Final Form: April 23, 1990)

Rate constants for the title reactions are computed by using variational transition-state theory with semiclassical ground-state adiabatic transmission coefficients for the temperature range from 200 to 2400 K. The rates are computed from selected information about the potential energy surface along the minimum energy path as parameters of the reaction path Hamiltonian. The potential information is obtained from ab initio electronic structure calculations with an empirical bond additivitycorrection. The accuracy of this semiempirical technique for obtaining the potential information is tested by comparing the results of the underlying ab initio calculations with higher quality multiconfiguration SCF and multireference CI calculations and by using increasingly higher quality ab initio electronic structure calculations before applying the bond additivity correction. For the reactions H + NH3 H2 + NH2, H2 + NH2 H + NH3, and D + ND3 D2 + ND2, the ultimate test is given by comparison with recent experimental results. Although the agreement is good in general, the comparisons of experiment and theory indicate that the computed barrier height is overestimated by about 1 kcal/mol. Quantitative agreement is found for all three reactions listed above by suitably scaling the potential information. The experimental results are extended down to 200 K and up to 2400 K for these three reactions, and predictions are provided for the reaction D2 ND2 D + ND3 for the temperature range from 200 to 2400 K.

-

-

-

+

I. Introduction Hydrogen atom abstraction reactions are ubiquitous in the combustion of hydrogen-containing compounds. As an example, the gas-phase chemical reaction H NH3 H2 + NH2 (R1)

+

+

is important in the thermal decomposition and combustion of ammonia. The rate of reaction R I has recently been measured by using three different experimental techniques;]+ these experiments combine to yield rate data from 500 K-where the slowness of the reaction tests the limits of the sensitivity of the experimental apparatus-up to 1777 K. These experiments are in good agreement with each other in temperature regions where they overlap, but overall they disagree with older experiment^^.^ that indicated the reaction was much slower. Compared to reaction R I , its reverse H2 NH2 H + NH3 (R2)

+

-

has been the subject of fewer direct measurement^,"^^^^ although it is also important in ammonia combustion. The most recent direct measurement of the rate of reaction R2 extends over the temperature range from 673 to 1003 K using a single experimental technique4 and combined with the earlier study7 provide reaction 'Present address: Molecular Science Research Center, Pacific Northwest Laboratory, Battelle Boulevard, Richland, WA 99352. Pacific Northwest Laboratory is operated for the US. Department of Energy by Battelle Memorial Institute under contract DE-AC06-76RLO 1830.

0022-3654/90/2094-1096$02.50/0

-

rates from 400 K up to 1003 K. By also studying reaction R I in the same experimental apparatus, estimates of the equilibrium constant for reactions R I and R2 and the thermochemistry of the N H 2 radical have also been ~ b t a i n e d .In ~ another recent study,9 the equilibrium constant for reactions R1 and R2 has been evaluated over the temperature range from 900 to 1620 K, and combined with earlier determinations of the rate of reaction R1 , I reaction rates for R2 were also obtained. The deuterated analogue of reaction R 1 D

+ ND3

-+

D2

+ ND2

(R3)

has also been studied by using a single experimental technique over the temperature range from 590 to 1220 K.Io The wealth (1) Michael, J. V.; Sutherland, J. W.; Klemm, R. B. Int. J . Chem. Kinet. 1985, 17, 315; J. Phys. Chem. 1986, 90, 497. (2) Marshall, P.; Fontijn, A. J. Chem. Phys. 1986, 85, 2637. (3) KO,T.; Marshall, P.; Fontijn, A. J . Phys. Chem. 1990, 94, 1401. (4) Hack, W.; Rouveirolles, P.; Wagner, H. Gg. J . Phys. Chem. 1986, 90, 2505. ( 5 ) Baulch, D. L.; Drysdale, D. D.; Horne, D. G. Eualuared Kinetic Duto f o r High Temperature Reactions; Butterworths: London, 1973; Vol. 2. (6) Hanson, R. K.; Salimian, S. Combustion Chemistry; Gardiner, W. C., Jr., Ed.; Springer: New York, 1984; Chapter 6. (7) Demissy, M.; Lesclaux, R. J . Am. Chem. SOC.1980, 102, 2897. (8) Holzrichter, K.; Wagner, H. Gg. Symp. (Inr.) Comb. [Proc.] 1980, 18, 769. (9) Sutherland, J. W.; Michael, J. V. J . Chem. Phys. 1988, 88, 830. (IO) Marshall, P.; Fontijn, A. J . Phys. Chem. 1987, 91, 6297.

0 1990 American Chemical Society

The Journal of Physical Chemistry, Vol. 94, No. 18, 1990 7097

Rate Constants for Gas-Phase Chemical Reactions of high quality experimental data for reactions RI-R3 provides an excellent testing ground for theoretical methods for predicting gas-phase reaction rates since all three reactions involve the same potential energy surface. To date, most comparisons for computed and experimental reaction rates are based upon conventional transition-state theory." The experimental rates of reactions RI and R3 have been compared with theoretical ones in refs 1, 2, and 9. These conventional transition-state theory calculations utilize information about the potential energy surface obtained from electronic structure calculationsI2 using a bond additivity correction (the BAC-MP4) method.I3-l7 The agreement of the theory and experiment indicates that the calculated barrier height for reaction R1 of about 18 kcal/mol (or about 16 kcal/mol when corrected for changes in zero-point energy) is reasonableS2 Curvature in the Arrhenius plot of the rate for reaction RI at low temperatures has been attributed to quantum mechanical tunneling, and attempts have been made to substantiate this claim by using simple tunneling models.2 One of the purposes of the present paper is to provide a more rigorous treatment of the dynamics controlling the rates of reactions RI and R3 and to also provide the first theoretical estimates of the rate for reaction R2 and predictions for its deuterated analogue Dz

+ ND2

4

D

+ ND3

(R4)

Previous calculations on reactions R1 and R3 employed conventional transition-state theory in which the dynamical bottleneck is assumed to occur at the saddle point for the reaction. As will be shown later, because of the nature of the semiempirical electronic structure calculations of the potential-energy-surface information utilized in these calculations, the use of conventional transition-state theory leads to inconsistencies in the treatment, but by basing the calculations upon variational transition-state theory,18-26these inconsistencies are no longer a problem. Another advantage of using variational transition-state theory is that it allows for a consistent treatment of quantum tunneling effects, which are important at low temperatures. Other goals of this work are to reexamine the importance of quantum mechanical tunneling in these reactions utilizing consistent semiclassical adiabatic transmission coefficients and to extend the rate data to lower temperatures. Compared to other dynamical methods for obtaining reaction rates, variational transition-state theory has the practical advantage

( I I ) Glasstone, S.;Laidler, K . J.; Eyring, H. The Theory of Rate Processes; McGraw-Hill: New York, 1941.

( I 2) Reported as private communications by C. F. Melius in refs 1,2, and 10. (13) Melius, C. F.; Binkley, J. S. Paper WSS/CI 83-16, 1983 Fall Meeting of the Western States Section of the Combustion Institute, University of California, Los Angeles, Oct. 17-18, 1983. (14) Melius, C. F.; Binkley, J. S. The Chemistry of Combustion Processes (ACS Symp. Ser. No. 249); Sloane, T. M., Ed.; American Chemical Society: Washington, DC, 1984; p 103. (15) Melius, C. F.; Binkley, J . S . Symp. (In!.) Comb. [Proc.] 1985, 20, 575. (16) Perry, R. A.; Melius, C. F. Symp. (fnr.)Comb. [Proc.] 1985,20,639. (17) Ho, P.; Coltrin, M. E.; Binkley, J . S.; Melius, C. F. J . Phys. Chem. 1985, 89, 4647; J . Phys. Chem. 1986, 90, 3399. (18) Isaacson, A. D.; Truhlar, D. G.; Rai, S. N . ; Steckler, R.; Hancock, G. C.; Garrett, B. C.; Redmon, M. J . Computer Phys. Comm. 1987, 47, 91. (19) Truhlar, D. G.; Garrett, B. C. Arc. Chem. Res. 1980, 13, 440. (20) Walker, R. B.;Light, J . C. Ann. Reu. Phys. Chem. 1980, 3 / , 401. (21) Pechukas, P. Ann. Rev. Phys. Chem. 1981, 32, 159. (22) Truhlar, D. G.; Isaacson, A. D.; Skodje, R. T.; Garrett, B. C. J . Phys. Chem. 1982.86, 2252; J . Phys. Chem. 1983, 87, 4554. (23) Truhlar. D. G.; Hase, W. L.; Hynes, J. T. J . Phys. Chem. 1983.87, 2664. (24) Truhlar, D. G.;Garrett, B. C. Ann. Reu. Phys. Chem. 1984,35, 159. (25) Truhlar, D. G.; Isaacson, A. D.; Garrett, B. C. The Theory of Chemical Reaction Rates; Baer, M., Ed.; CRC Press: Boca Rgton, FL, 1985; Vol. 4, p 65. (26) For Version 1.1, see ref 18. For Version 1.5, see: Isaacson, A. D.; Truhlar, D. G.; Rai, S . N . ; Hancock, G.C.; Lauderdale, J. G.; Truong, T. N.; Joseph, T.; Garrett, B. C.; Steckler. R. POLYRATE Program manual, 1988.

of not requiring a global representation of the potential energy surface. From knowledge of the potential in a region around the minimum energy path connecting reactants to products, accurate estimates of the reaction rate can be made. This type of potential information can be directly obtained from modern electronic structure methods as the energy and derivatives of the energy. Previous conventional transition-state theory calculations of the rates for reactions R1 and R3 have utilized potential information obtained by the BAC-MP4 One of the most important purposes of this paper is to critically evaluate the method used to obtain this potential-energy-surface information. Section 2 of this paper presents the theoretical methods employed: variational transition-state theory is reviewed with a focus on the potential-energy-surface information that is required in these calculations, the electronic structure calculations are summarized with a focus on the bond additivity corrections, and the methods used to estimate anharmonic contributions are briefly described. Section 3 presents computational details of the electronic structure and rate constant calculations. Discussions of the calculations and comparisons with experiment (when available) for the title reactions are provided in section 4 and conclusions are presented in section 5. 2. Theory 2.1. Variational Transition-State Theory. Variational transition-state theory, including semiclassical transmission coefficients, have been extensively reviewed in the literature.18-26In the present work we briefly review the method with an emphasis on the potential energy surface information needed for the calculations. Conventional transition-state theory" reduces the calculation of the rate constant to one of quasiequilibrium statistical mechanics; the equilibrium rate is approximated as the equilibrium flux headed toward reactants through a dividing surface located at the saddle point. With this approximation, the rate constant kt( 7') for temperature T takes on the simple familiar form

where u is a symmetry factor, kB is the Boltzmann constant, h is the Planck constant, Qt( 7') is the partition function for the bound 7') is the reactants degrees of freedom at the saddle point, aR( partition function, and V' is the value of the potential at the saddle point. Thus conventional transition-state theory requires information about the potential energy surface only in the saddle-point and reactant regions. In particular, the value of the potential at the saddle point (relative to the reactant value) is required and if the partition functions are computed by using a harmonic approximation, then the matrix of second derivatives (the Hessian matrix) suffices. In variational transition-state theory, the dividing surface is viewed as a tentative dynamical bottleneck to flux in the product direction, and the best bottleneck (the dividing surface allowing the least flux) is located by a variational procedure. The generalized expression for the thermal rate constant is given as a function of the location s of the dynamical bottleneck along the reaction coordinate

where PT(T,s)is the generalized partition function for the bound degrees of freedom orthogonal to the reaction path at s and V,,,(s) is the value of the potential along the reaction path at s. One version of variational transition-state theory, the canonical variational theory (CVT),27results from minimizing eq 2 with respect to s: kCVT(T)= min kGT(T,s) S

(3)

(27) Garrett, B. C.; Truhlar, D. G. J . Phys. Chem. 1979,83, 1079; J . Phys. Chem. 1983, 87, 4553.

7098 The Journal of Physical Chemistry, Vol. 94, No. 18, I990 The improved canonical variational theory ( ICVT)28 also variationally optimizes the location of the transition-state dividing surface for a given temperature but provides an improved treatment of threshold energies by using an ensemble that removes energies below the ground-state adiabatic threshold. More information about the potential energy surface is required than for a conventional transition-state theory calculation to compute the rate constant using either the canonical or improved canonical variational theory; information about the potential in a region around the reaction path is also required. The reaction path is defined as the minimum energy path connecting the saddle point with both the reactant and product regions. The minimum energy path is located by following the path of steepest descents in both directions from the saddle point in a mass-weighted coordinate system such that each degree of freedom has the same effective mass in the kinetic energy expression. The first step from the saddle point is taken along the eigenvector of the negative eigenvalue of the Hessian, and subsequent steps are taken in the direction of the negative of the gradient (the first derivative of the energy). For a harmonic treatment of the partition functions, the matrix of second derivatives along the minimum energy path will suffice. In this case the potential information needed is the energy and its first and second derivatives along the minimum energy path. In the above expressions for the rate constant, the bound degrees of freedom are treated quantum mechanically (i.e., the partition functions are evaluated quantum mechanically) but the motion along the reaction coordinate (the degree of freedom corresponding to the negative eigenvalue of the Hessian at the saddle point) is treated classically. Quantum mechanical effects on the reaction coordinate motion (e.g., quantum mechanical tunneling) are included by a multiplicative factor-the transmission coefficient.27 A consistent route to including quantum mechanical effects on reaction coordinate motion is provided by the vibrationally and rotationally adiabatic theory of reaction^.^^-^' When reaction coordinate motion is treated classically, one form of variational transition-state theory-the microcanonical variational theoryyields an expression for the thermal rate constant which is equivalent to that obtained from the adiabatic theory of react i o n ~ . ~Thus ~ , ~it ~is consistent within variational transition-state theory to treat tunneling as occurring through the one-mathematical-dimensional vibrationally-rotationally adiabatic potential (4) where a is a collective index of the quantum numbers for the bound modes and e,GT(s)is the bound energy level for state a at the generalized transition state located at s along the reaction path. For thermal rate constants the tunneling is approximated by using only the ground-state adiabatic potential curve (a = 0), which is denoted The adiabatic approximation is made in a curvilinear coordinate system, and, although the potential term is very simple, the kinetic energy term is complicated by factors dependent upon the curvature of the reaction path. For systems in which the curvature of the reaction path is not too severe, the small-curvature semiclassical adiabatic ground state (SCSAG) m e t h ~ d ~includes l - ~ ~ the effect of the reaction-path curvature to induce the tunneling path to "cut the cornern and shorten the tunneling length. Combining the SCSAG transmission coefficient with the improved canonical variational theory rate constant yields kICVT/SCSAG(T ) = KSCSAG(T)~ICVT(T ) (5)

c(x).

Garrett et al. The type of potential information required is identical with that needed for the variational transition-state theory calculation to construct the adiabatic potential. For the SCSAG calculation it is also necessary to know the curvature of the reaction path, which can be obtained from the first and second derivatives of the potential along the reaction path. Thus, no new information about the potential energy surface is required to provide a consistent and accurate estimate of the tunneling. The information about the potential energy surface is used in the POLYRATE program26in the form of parameters of the reaction path H a m i l t ~ n i a nat~ ~a set of points along the reaction path. Because of the expense of the electronic structure calculations the information about the potential is often limited, and the parameters of the reaction path Hamiltonian are interpolated from a sparse grid and can be extrapolated out to the asymptotic reactant and product 2.2. Bond Additivity Corrections. Despite the developments in quantum chemistry in the past decade, it remains a challenge to calculate accurate potential energy surfaces useful for determining rates of activated chemical reactions. This is largely due to the fact that a small change in the barrier height translates into a large change in the reaction rate (e.g., a change of only 0.5 kcal/mol at room temperature leads to a change in the rate of over a factor of 2 ) . Currently, ab initio calculations of energetics to within chemical accuracy is affordable only for relatively small systems. Even so, progress has been made in developing methods for determining reaction energetics, but these have relied on empirical corrections to ab initio calculations. The concept of a molecule as a collection of individual bonds between pairs of atoms is central to the chemist's view of molecular structure and its rearrangement during reactions. The order of a bond (e.g., single bonds between monovalent atoms, double bonds between divalent atoms, etc.) has been used in establishing empirical rules relating bond order to bond lengths36 and bond strength^.^' The empirical knowledge that the dissociation energy for a given bond with a given bond order is, in general, roughly the same in different molecules has lead to the chemist's intuition that bond strengths are approximately additive. This concept of bond additivity has been successfully applied in obtaining estimates of the thermochemistry of stable species.38 Within the last several years, this concept has been used in a more quantitative approach (the BAC-MP4 method) in which corrections are made to ab initio calculations.I 3 - I 7 Although based upon simple bond additivity, the BAC-MP4 method goes beyond pair interactions to include three-body terms. The BAC-MP4 method has been successfully applied in predicting the thermochemistry of a wide variety of stable species and transient radical species that occur in combustion environments. Inherent in the BAC-MP4 method is the assumption that geometries and frequencies can be obtained to a good approximation from a Hartree-Fock calculation using a moderately sized basis. Therefore, a BAC-MP4 calculation begins with a spin-unrestricted Hartree-Fock (HF) optimization to the equilibrium geometry xHFs0 (Le., the location at which the gradient vanishes) with a contracted atomic orbital basis. In the present application we use a basis that is denoted 6-31G*.39 At this same critical geometry using the same basis, the Hessian matri&HHF(xHF*O) is also computed at the H F level. The energy levels e a (xHFs0)for the bound degrees of freedom are computed within the harmonic approximation by t,HF =

x ( a a+ '/,)hwFF

(6)

n

where the sum is over the bound degrees of freedom and the (28) Garrett, B. C.; Truhlar, D. G.;Grev, R. S.;Magnuson, A. W. J , Phys. Chem. 1980, 84, 1730; J . Phys. Chem. 1983, 87,4553. (29) Marcus. R . A. J . Chem. Phys. 1966,45,4493,4500 J . Chem. Phys. 1967, 46, 959; J . Chem. Phys. 1968, 49, 2610. (30) Truhlar. D. G . J . Chem. Phys. 1970, 53, 2041. ( 3 1 ) Skodje, R . T.; Truhlar, D. G.; Garrett, B. C. J . Chem. Phys. 1982, 77. 5955. (32) Garrett, B. C.; Truhlar. D. G. J. Phys. Chem. 1979, 83. 1052;J . Phys. Chem. 1983.87, 4553. (33) Skodje, R. T.; Truhlar, D. G.; Garrett, B. C. J . Phys. Chem. 1981, 85. 30 I9

(34) Miller, W. H.; Handy, N . C.; Adams, J. E. J . Chem. Phys. 1980, 72, 99. (35) Truhlar, D. G.; Kilpatrick, N. T.; Garrett, B. C. J . Chem. Phys. 1983, 78, 2438. (36) Pauling, L. J . Am. Chem. Soc. 1947, 69, 542. (37) Johnston, H. S. Gas Phase Reaction Rare Theory; Ronald Press: New York, 1966. (38) Benson, S. W. Thermochemical Kinetics; Wiley: New York, 1976. (39) Krishnan, R.; Binkley, J. S.;Seeger, R.; Pople, J. A. J . Chem. Phys. 1980. 72. 650.

The Journal of Physical Chemistry, Vol. 94, No. 18, 1990 1099

Rate Constants for Gas-Phase Chemical Reactions

+

TABLE I: BAC-MP4 Parameters for the H NHI Reaction System HFIBAC-MP4 MP2/BAC-MP4 Af/foOH (kcal/mol) AfHoON(kcal/mol) EFP4(hartree) E#'' (hartree) FH FN AH, (kcal/mol) ANH (kcal/mol) aHH aNH

Ro (A)

51.63 112.53 -0.498232 -54.473256 0.0 0.2174 18.98 69.66 2.0 2.0 1.4

51.63 112.53 -0.498232 -54.473256

0.0 0.2174 18.92 68.06 2.0 2.0 1.4

harmonic frequencies are obtained by diagonalization of the mass-weighted Hessian matrix [UHF12 = UT[M1/2HHFMl/2]U

(7)

where the matrix M is a diagonal matrix with the appropriate mass for each degree of freedom on the diagonal and U is the orthogonal matrix that diagonalizes the mass-weighted Hessian. For a stable equilibrium geometry there will be 5 or 6 zero eigenvalues for a linear or nonlinear molecule, respectively. The remaining nonzero harmonic frequencies appear in eq 6. The absolute energy for a given geometry is much harder to accurately compute and requires inclusion of a large percentage of the electron correlation energy. As a first step in this direction, at the optimized H F geometry, the energy EMP4(xHF9O) is computed by full fourth-order Mdler-Plesset perturbation theory40 (MP4-SDTQ), including the effects of all single, double, triple, and quadruple excitations in the contracted atomic orbital basis denoted 6-31G** 39 in which the 1s electrons of the nitrogen are frozen. Although the MP4 method does include a large fraction of the electron correlation energy, there are still errors due to basis set and higher order correlation effects, so a correction based upon bond additivity-the BAC-is applied. For stable equilibrium geometries the HF/BAC-MP4 method approximates the heat of formation at 0 K by EHFlBAC-MP4 = EMP~(XHF,O)+ ,$F(xHF.O) + C(A,HO,, i

matrix are computed at the H F level as well as the fact that the absolute energies are computed by the MP4 method. In addition, the values of the parameters depend on the nature of the basis sets used in the ab initio calculations, but we do not study basis set behavior in the present work. The parameters for the HF/ BAC-MP4 method used in this study are presented in Table I. The HF/BAC-MP4 method has also been used to calculate the energetics of transition states for chemical reaction^.'^,^^ Similar to stable equilibrium geometries, the saddle-point geometry xHF.+is obtained at the H F level by locating the geometry at which the H F gradient vanishes. Energy levels for the bound degrees of freedom are obtained as they are for stable equilibrium geometries by using eqs 6 and 7 except that now one of the eigenvalues of eq 7 is negative and corresponds to the unbound motion leading from the saddle point down to reactants on one side and products on the other. The zero-degree heat of formation of the transition state is determined by a MP4 calculation followed by a bond additivity correction and is given by eq 8 but with the geometries xHFs0replaced by xHFvt. The HF/BAC-MP4 method can also provide a practical method for obtaining information about the potential along the minimum energy path. After the saddle point is located at the H F level, the minimum energy path xHF(s) is located by following the path of steepest descents into the reactant and product regions. The reaction coordinate s measures the distance along the minimum energy path from the saddle point and by definition is negative on the reactant side of the saddle point. At each point along the minimum energy path the gradient vector gHF(s)and Hessian matrix HHF(s)are computed at the H F level by using the 6-31G* basis and these are also used in determining the next point along the minimum energy path.42v43 (Functions of s such as gHF(s) mean that they are evaluated at the geometry xHF(s)along the minimum energy path.) Each point along the minimum energy path represents a generalized transition state, and the heats of formation of these generalized transition states are obtained by using the 6-31G** a calculation of the MP4-SDTQ energy EMP4(s) basis followed by the bond additivity correction to give an expression analogous to eq 8 EHF/BAC-MP4(~) = EMP4(s) ttF(s) C(AfHooi - ETP4)-

+

+

n

- C Aij eXp[-~~ijRi,(X~~")]f,k(XHF") (8) ijbonds

k#ij

where tOHF(xHF*O) is obtained from eq 6 with CY = 0, AHo . is the experimental atomic heat of formation of atom i, and E$' is the atomic energy computed at the MP4 level for atom i. The final term is a sum over all bonds in the molecule; specifically for the reactant and product species considered here we include one H H bond in H2. two N H bonds in NH2, and three N H bonds in NH3. The parameters A ' . and aijdetermine the strength and range of the correction, a n J the termJjk damps out the correction for a bond between atoms i and j when another atom k is bonded to either atom. This damping factor takes the form Lp(XHF") = { I - FjFk eXp[-2CUik(Rik(XHF*o)- &)I} x { I - FiFkexp[-2ajk(Rjk(xHF,O)- Ro)]J (9) where the parameters Fi are dependent only on the type of atom i.

The parameters of the HF/BAC-MP4 method are fitted to a database of experimental heats of formation of molecular species in a least-squares sense. By definition, the atomic heats of formation are exactly reproduced. The parameters A,, cyi,, Fi, and Ro are fitted to reproduce diatomic and polyatomic heats of formation. Because the parameters are fitted to heats of formation of molecular species, the energy EHF/BAc-MP4 contains corrections to the zero-point vibrational energy of the molecules as well as to the absolute energies. The values of the parameters are therefore dependent upon the fact that the geometry and Hessian (40) Pople, J . A.; Binkley, J . S.;Seeger, R. I n ( . J . Quantum Chem. Symp. 1976, 10. I . Pople, J . A.; Krishnan, R.; Schlegel, H. B.; Binkley, J. S . Int. J . Quantum Chem. Symp. 1978, 14, 545.

where the zero-point energy along the minimum energy path is computed by projecting out the motion along the gradient vector gHFbefore diagonalizing the mass-weighted Hessian matrix.j4 The final term is a sum over all bonds in the reaction complex; specifically for the reactions considered here we include one H H bond and three N H bonds throughout the reaction. For the purposes of carrying out variational transition-state theory calculations, we need to define the classical potential along the minimum energy path, neglecting zero-point energies. The heat of formation at 0 K of a generalized transition state at s is equal to the ground-state adiabatic potential @(s) defined in eq 4 (to within an additive constant AEo defining the zero of energy) EHF/BAC-MP4(s) = AEo (1 1)

c(s)+

where the zero of energy is chosen to be the classical energy of the reactants so that AEo = EHF/BAC-MPd(s = -a)- c f T ( s = -a) (12) With use of eq 4 for the adiabatic potential, an expression for the HF/BAC-MP4 approximation to the energy along the minimum energy path is given by 14vFAC-MP4(~) = EHF/BAC-MP4(~) - tfT(s) - AEo (1 3) With this definition, the ground-state adiabatic barrier and, (41) Marshall, P.; Fontijn, A.; Melius, C. F. J . Chem. Phys. 1987, 86. 5540. (42) Garrett, B. C.; Redmon, M. J.; Steckler, R.; Truhlar, D. G.; Baldridge, K. K.; Bartol, D.; Schmidt, M. W.; Gordon, M. S.J . Phys. Chem. 1988, 92, 1476. (43) Page, M.; McIver, J . W., Jr. J . Chem. Phys. 1988, 88, 922.

7100

The Journal of Physical Chemistry, Vol. 94, No. 18, 19919

therefore, the heat of reaction are independent of how the zeropoint energy level is computed in the variational transition-state theory calculations. The gradient vector gHFand Hessian matrix HHFused in the HF/BAC-MP4 calculations of saddle point energetics correspond to an actual saddle point (Le., the gradient vanishes and the Hessian has one negative eigenvalue). However, this is a saddle point on the H F potential ener surface and the energy along which is calculated at the minimum energy path G$'c.Mp4(s), the MP4 level and includes the bond additivity correction, is not necessarily at its maximum at this saddle point. This inconsistency of the potential information makes the use of conventional transition-state theory questionable since it is based on a knowledge of the potential energy surface only near the H F saddle point and the BAC-MP4 potential energy can increase appreciably upon moving off of the saddle point along the minimum energy path. Variational transition-state theory utilizes information all along the minimum energy path and the saddle point is not singled out as a special location. The variations of both the potential and the frequencies of the bound modes normal to the reaction path are included in the variational transition-state theory calculation, and the fact that the gradient vector vanishes at one particular location along the minimum path is of no practical consequence in the calculations. Thus, the HF/BAC-MP4 method provides a practical method for obtaining information about the potential, and variational transition-state theory provides for a consistent use of this information in calculating rate constants. A major concern in using the HF/BAC-MP4 method is the accuracy of the computed potential energy surface information. One test of this is the comparison of computed rates with experimental ones. Tests of this nature can sometimes be inconclusive because cancelations of errors can lead to fortuitous agreement. Another test is to improve the quality in the underlying ab initio calculation to decrease the magnitude of the empirical bond additivity correction and test the effect on the computed rate constant. As a first step in this direction, we have reevaluated the potential information using the MP2/BAC-MP4 method. In or this method the geometry of the equilibrium structure xMP2,O transition state xMP2*+ is determined by using a full second-order Mdler-Plesset perturbation theory (MP2) calculation including the effects of all single and double excitations in the contracted atomic orbital basis denoted 6-31G*. The matrix of second derivatives HMP2is also calculated at the MP2 level at the same geometry. The minimum energy path xMP2(s)is located by following the path of steepest descents by using the gradient vectors and Hessian matrices computed at the MP2 level. For each point is computed along the minimum energy, the energy EMP4[xMm(s)] by the MP4-SDTQ method using the 6-31G** basis followed by a bond additivity correction. The MP2/BAC-MP4 approximations to the heats of formation of the equilibrium geometries, saddle point, and generalized transition states are given by equations similar to eq 8 and 10 except that the superscript H F is replaced by MP2 with their obvious meanings. Although both the HF/BAC-MP4 and MP2/BAC-MP4 methods calculate the energy by the MPCSDTQ method, because the geometries, vibrational frequencies, and MP4 energies are different at the H F and MP2 levels, the bond additivity correction is also different. For the MP2/BAC-MP4 method the parameters are refitted to experimental heats of formation as in the HF/BAC-MP4 method, and the parameters used in the MP2/BAC-MP4 calculations are shown in Table I. 2.3. Anharmonicity. In a variational transition-state theory calculation based upon a global potential energy surface, the effects of anharmonicity in the potential energy surface are included, using the independent-normal-mode approximation by including the principal anharmonicities in each mode.25 When the potential energy surface information is included as parameters of the reaction path Hamiltonian, anharmonic effects can still be incorporated by including third and fourth derivatives of the potential along the independent normal modes. The effects of anharmonicity are important and must be included to obtain accurate estimates of the rate constants.

Garrett et al. The bond additivity correction is based upon information about the potential energy surface only through quadratic terms but is designed to fit the experimental heats of formation of bound molecules, which include the effects of anharmonicity. Therefore, the BAC-MP4 methods have anharmonic effects built in to them in an averaged sense, not on a mode-by-mode basis. The harmonic frequencies computed at the Hartree-Fock level are known to be systematically overestimated by an average of about 10-12%, and the BAC method corrects for this discrepancy as well as including the effects of anharmonicity in calculating the zero-degree heat of formation. At finite temperatures, effects of anharmonicity enter the thermal rate constants through partition functions that require calculations of excited-state energy levels for each mode. Since even the harmonic frequencies of stable species are known to be accurate to only within IO-12% when calculated by HF, the effect of anharmonicity in excited-state energy level is neglected. As a test of the importance of excited states in the calculation of partition functions, the ICVT calculations including excited-state energy levels computed directly from the H F frequencies uHFare compared with those in which the H F frequencies are systematically lowered, in this case by a scale factor of 0.893, which is representative of an average of frequency corrections for a variety of compounds containing C, N, 0, and H.44 Anharmonic effects are most important for the lowest frequency modes, and for the reactions studied here the lowest frequency mode of the reaction complex is a bound twist vibrational motion in the interaction region but becomes a hindered rotor in the asymptotic regions. For this mode the harmonic approximation becomes more suspect especially at high temperatures, and a better treatment is offered by a hindered rotor Previous calculations on this system' have found that at the H F saddle point the hindered rotor partition function and harmonic oscillator partition function differ by only 15% at 1500 K. Therefore, the error in the rate constant from the harmonic approximation to this hindered rotor mode will be negligible compared to other errors in the calculations, such as uncertainties in the computed barrier height. In the present study, all modes are treated within the harmonic approximation with anharmonic effects entering through the bond additivity correction. 3. Computational Details 3.1. Quantum Chemistry Calculations. Details of the electronic structure calculations underlying the BAC-MP4 method are given in section 2.2. The Hartree-Fock, MP2, and MP4 calculations are performed by using the GAUSSIAN86 electronic structure code.% As a test of the accuracy of these electronic structure calculations and to assess the reliability of the BAC-MP4 method for computing the reaction barriers, a limited number of multiconfiguration S C F (MCSCF) and multireference configuration interaction (MRCI) calculations are carried out. Details of the MCSCF and MRCI calculations are presented here. Typical of hydrogen abstraction by radical^,^' the significant electronic structure changes accompanying reaction R 1 are predominantly confined to three electrons and three molecular orbitals ( M a ) . For the reactant, these include the hydrogen atom electron and MO as well as the two electrons in the attacked N H bond of NH3 with the associated bonding and antibonding MOs. At the product side, these three electrons and MOs are the nitrogen radical electron and MO of N H 2 as well as the electrons and bonding-antibonding pair of MOs of the hydrogen molecule. We ~~

~~~

(44) Pople, J . A.; Schlegel, H. B.; Krishnan, R.; DeFrees, D. J.; Binkley, J . S.; Frisch, M. J . : Whiteside, R. A. Inr. J . Quantum Chem. Symp. 1981, I S , 269. (45) Pitzer, K. S.; Gwinn, W. D. J . Chem. Phys. 1942, 10, 428. Pitzer, K. S. J . Chem. Phys. 1946, 14, 239. Pitzer, K. S. Quantum Chemistry; Prentice-Hall: New York, 1953. Davidson, N. Staiisrical Mechanics; McGraw-Hill: New York, 1962. (46) Gaussian 86. Frisch, M. J.; Binkley, J . S.;Schlegel, H. B.; Raghavachari, K.; Melius, C. F.; Martin, R. L.; Stewart, J. J. P.; Bobrowitz, S. W.; Rohlfing, C. F.; Kahn, L. R.; Defrees, D. J.; Seeger, R.; Whiteside, R. A,; Fox, D. J.: Fluder, E. M.; Topiol, S.; Pople, J. A. Carnegie Mellon Quantum Chemistry Publishing Unit, Carnegie Mellon University, Pittsburgh, PA. (47) Dunning, T. H . , Jr. J . Phys. Chem. 1984, 88, 2469.

Rate Constants for Gas-Phase Chemical Reactions choose, as a multiconfiguration (MCSCF) reference description for this reaction, a complete active space S C F (CASSCF)@wave function including three electrons in three active orbitals. This eight-configuration MCSCF wave function provides a qualitatively correct description of the orbital changes mentioned above and allows the orbital occupation numbers to smoothly change from their reactant to their product values. The molecular orbitals generated at the reference CASSCF level are used in large-scale MRCI calculations to obtain accurate energetics. The MRCl calculations are obtained in the frozen core approximation; the molecular orbitals that correlate with the Is core orbital of nitrogen and with the associated virtual orbital are excluded in the configuration selection and are thus forced to be doubly occupied and unoccupied respectively. The MRCI are fully second order; all doublet configurations are included that can be generated by single or double substitutions from all of the MCSCF reference configurations. Wave functions of this type have been shown-by extensive comparisons with full-C1 calculations in the same basis set-to provide accurate descriptions of bond breaking.49 In addition to the 6-31G* basis set,39two other one-electron basis sets, denoted pVDZ and pVTZ, are used in the MCSCF and MRCl calculations. Recently introduced by Dunning,so these new basis sets have been optimized for highly correlated wave functions and are referred to as correlation-consistent basis sets to emphasize the well-defined energy-lowering patterns. There are 35 contracted basis functions in the pVDZ basis set including N:(3s,2p,ld) and H:(2s,lp). There are 95 contracted functions in the pVTZ set including N:(4s,3p,2d,lf) and H:(3s,2p,ld). The MRCl calculations include about 50000 configurations with the pVDZ basis set and over 500000 configurations with the pVTZ basis set. The MCSCF and MRCI electronic structure calculations were all performed by using the MESAs' system of programs on the Cray XMP-24 at the Naval Research Laboratory. 3.2. Rate Constant Calculations. Rate constant calculations were carried out by using the POLYRATE program.26 Several algorithms for determining the minimum energy path have been tested42and one of the most efficient methods when numerical accuracy of about 15% is required is the lowest order Euler one-step algorithm. This method uses only the information about the gradient at the current geometry to predict the next geometry along the minimum energy path. In the present study the minimum energy path was located by using a gradient-following algorithm suggested by Page and M ~ I v e rwhich , ~ ~ is based upon a local quadratic approximation to the potential and utilizes both the gradient vector and Hessian matrix in determining the new geometry along the reaction path. The initial step from the saddle point was taken along the eigenvector of the negative eigenvalue of the Hessian matrix. We performed tests of the convergence of the ICVT/SCSAG rate constant with respect to the step size taken in the gradient-following algorithm and the step size at which the energy and derivative information was stored for use in interpolation of the parameters of the reaction path Hamiltonian. Since the Hessian matrix is used in the gradient-following algorithm, the two step sizes are the same. Decreasing the step size from 0.10 to 0.05 a. was found to change the computed rates at 300 K by less than 20% and a step size of 0.05 a. is used in this work. The tunneling correction factor computed by using the SCSAG method is sensitive to the extent of the reaction path. Extending the reaction path to 0.8 a. toward the H N H 3 asymptote and 1.2 a. toward the H 2 + N H 2 asymptote was sufficient to converge the computed rates to within 10%at 300 K and about 25% at 200 K. For the D + ND, reaction and its reverse (R4), the reaction path was extended out to 0.36 a. toward the D + ND, asymptote

+

(48) Roos, B. 0.;Taylor, P. R.; Siegbahn, P. E. M. Chem. Phys. 1980,48, 152. (49) See, e&: Bauchlicher, C. W.; Taylor, P. R. J . Chem. Phys. 1987, 86, 5600 and references therein. (50) Dunning, T. H., Jr. J. Chem. Phys. 1989, 90, 1007. ( 5 1 ) MESA (Molecular Elechonic Structure Applications), Saxe, P.; Lengsfield, B. H.; Martin, R.; Page, M.

The Journal of Physical Chemistry, Vol. 94, No. 18, 1990 7101

-

-

TABLE 11: Thermochemistry for Reactant and Product Species for the Reactions NH3 + H NH2 + H2 and ND3 D ND2 Df (Energies in kcal/mol) NHi H NH, NDi D ND, ArH' ( T = 0 K) -9.36 51.6 46.8 -15.5 51.6 43.3 ArH' ( T = 300 K) A@' ( T = 300

K)

A@' ( T = 600 K) A@' ( T = 1000 K) A@' ( T = 1500 K)

+

-9.3' -1 1.0 -11.0 -4.0 -3.9 3.7 3.8 14.7 14.8 28.9 28.8

51.6 52.1 52.1 48.5 48.6 44.8 44.8 39.5 39.6 32.5 32.6

46.2 46.0 45.5 48.3 47.8 50.8 50.3 54.6 54.1 59.5 59.1

-12.3 -17.2 -14.0 -9.4 -6.2 -1.1 2.1 10.4 13.7 24.9 28.1

+

52.5 52.1 53.0 48.4 49.3 44.6 45.5 39.1 40.1 32.0 33.1

45.0 42.6 44.3 45.2 46.9 48.0 49.8 52.0 53.8 57.1 59.0

'By definition the heats of formation and free energies of formation of H2 and D2 are zero a t all temperatures. bTop entry is from HF/ BAC-MP4 calculations. 'Lower entry is from J A N A F tables (ref 53).

-

TABLE 111: Saddle-Point Geometries for the H + NH3 H2 + NH2 Reaction Computed at Four Levels of Theory' HF MP2 MP4 CASSCF 1.244 1.010 1.010 0.949 166.7 105.4 100.5 100.5

1.305 1.023 1.023 0.869 158.5 103.5 99.1 99.1

1.277 1.026 1.026 0.893 158.3 103.3 99.3 99.3

1.283 1.010 1.010 0.962 170.0 107.6 99.2 99.2

'The hydrogen atom denoted H' is the one furthest from the N that is forming a bond with hydrogen H. Atoms H" and H"' are spectators not directly involved in the abstraction reaction.

+

and 0.88 a. toward the D2 ND2 asymptote. This was sufficient to converge the ICVT/SCSAG rate to better than 15% a t 300 K but an uncertainty of about 50% in the 200 K rate remains. The tunneling correction factors are also sensitive to the calculations of the reaction path curvature components. As computed here, these require evaluation of either derivatives of the eigenvectors of the Hessian with respect to the reaction coordinate or derivatives of the gradient vector with respect to the reaction coordinate. Since the gradient vector and Hessian matrix are available only on a fairly coarse grid, these derivatives are not well convergedes2 Tests indicate that this uncertainty in the reaction path curvature components will cause errors of about 30-40% at 200 K but less than 20% for 400 K and above. 4. Results and Discussion 4.1. Thermochemistry of the Reactants and Products. As a

guideline to the accuracy expected from the HF/BAC-MP4 method in computing thermochemical data for bound molecular species, a comparison is presented of experimentals3 and computed heats of formation and free energy of formation of the reactant and product species of the reactions studied here. Note that by definition the experimental heat of formation and free energy of formation of H2 and D2 are zero at all temperatures. The HF/BAC-MP4 results for these species are 0.03 kcal/mol for H2 and -1.9 kcal/mol for D2 for both the heat of formation and free energy of formation at all temperatures. The thermochemical data for the remaining species are shown in Table 11. The (52) Although not done in this study, we note that the curvature compcnents can be calculated analytically given the Hessian matrix and the gradient at only one point, avoiding numerical difficulties associated with finite difference techniques (see ref 43). The accuracy of these components, however, is still indirectly dependent on the step size for evaluating the derivative information since the curvature components are only correct right on the minimum energy path. (53) Chase, M. W., Jr.; Davies, C. A,; Downey, J. R., Jr.; Frurip, D. J.; McDonald, R. A.; Syverud, A. N. 'JANAF Thermochemical Tables, 1985 Supplement", J. Phys. Chem. Ref. Daia 1985, 14.

7102

Garrett et al.

The Journal of Physical Chemistry, Vol. 94, N o . 18, 1990

-

TABLE IV: Saddle-Point Harmonic Frequencies and Zero-Point Energies (cm-I) for the H NH3 H2 + NH2 Reaction Computed at Five Levels of Theory CASSCF

+

zero-point energy

HF'

MP2"

MP4b

6-31G*

DVDZ

28333 723 732 1258 1485 1541 1731 3647 3750 7434

20211' 718 730 1146 1315 1610 2077 3489 3606 7345

19053 726 901 1186 1347 1629 1921 3497 3553 7381

26041' 682 704 1253 1434 1435 1725 3644 3748 7314

24133 666 674 1197 1366 1535 1678 3629 3719 7232

'Computed with the 6-31G* basis set. bComputed with the 631G** basis set.

BAC-MP4 results are in excellent agreement with the J A N A F tables for all species. The thermochemistry of the deuterated species have the largest uncertainties, but the HF/BAC-MP4 method presents a consistent procedure for obtaining reliable estimates at all temperatures. We note that to obtain the excellent agreement with experiment, the H F frequencies were scaled by 0.893 in the BAC-MP4 calculations of the thermochemistry. We will also employ the frequency scaling in calculating rate constants for comparison with experiment. 4 . 2 . H 4- N H , H 2 4- N H 2 and H 2 -+ N H 2 -c H -+ N H , Reaction Rates. Although most of the electronic structure calculations were carried out at the HF/BAC-MP4 or MP2/BACMP4 level, limited calculations of saddle point geometries and frequencies were also carried out by using MP4 and CASSCF calculations. The CASSCF and CAS-MRCI calculations also tested the dependence of the reaction barrier on basis set. These multireference calculations indicate that the electronic structure of the saddle point is dominated by a single configuration. Therefore, the H F calculations that are the basis of the BAC-MP4 calculations are a good zero-order approximation, giving us more confidence in the BAC-MP4 calculations of the reaction-path Hamiltonian information. I n Table 111 we examine the saddle point geometries at four levels of theory-HF, MP2, MP4-SDTQ (denoted simply MP4), and CASSCF. The HF, MP2, and CASSCF calculations all utilize a 6-31G* basis, whereas the MP4 calculations use a 631G** basis. The geometries change only slightly as the level of the ab initio electronic calculation is improved, indicating that the lowest order calculations, the H F method, give a reasonable estimate of the geometries. The largest changes are for the N H and HH' bond lengths that are most actively involved in the bond rearrangement in the reaction. Note that the geometries do not change monotonically with improvement in the level of treatment from H F to MP2 to MP4, but the H F values for RNH and RHH, are closer to the MP4 values than they are to the MP2 values. The saddle point frequencies are compared to the HF, MP2, MP4, and CASSCF levels for the reaction of H with NH3 in Table IV. Again the H F and MP2 calculations use a 6-31G* basis and the MP4 calculations use a 6-31G** basis. The CASSCF calculations were performed with two basis sets, the 6-31G* basis and a valence double-< (pVDZ) basis. The imaginary frequencies at the top of each column correspond to the unbound mode for motion along the reaction coordinate. Note that of all the frequencies this imaginary one changes the most as the level of theory is improved. The harmonic zero-point energy given at the bottom of each column is one-half the sum of the frequencies for all bound modes. Although the frequencies for some individual modes do change by as much as 206, the changes in zero-point energy are much smaller. A change of 1 I O cm-' in the zero-point energy (the difference between the H F and MP2 values) will shift the adiabatic potential curve by only 0.3 kcal/mol. The energetics of the reaction computed at several levels of theory are given in Table V . The first column is the difference in I'ME~(s) evaluated at products and reactants. When the dif-

-

TABLE V: Reaction Energetics (kcal/mol) for the H + NH2 Reaction Computed at Seven Levels of Theory

H F (6-31G*) H F (6-3lG**) MP2 M P4 CASSCF (pVDZ) CASSCF (pVTZ) CAS-MRCI (pVDZ) CAS-MRCI (pVTZ) HF/BAC-MP4 (w unscaled) HF/BAC-MP4 ( w scaled) MP2/BAC-MP4

+ NH,

reaction energy'

forward barrierb

-1.22 -1.50 8.68 3.28 1.63 1.94 4.28 5.26 8.02 1.63 7.69

24.58 23.70 21.59 17.68 22.69 24.17 15.24 16.56 16.44 16.23 16.70

-

H2

reverse barrief 25.80 25.19 12.92 14.40 2 I .05 22.22 10.96 1 I .30 8.42 8.60 9.01

a Energies of products at their equilibrium geometries relative to the energy of reactants at their equilibrium geometries. bEnergiesof the classical barrier maximum relative to the energy of reactants at their equilibrium geometries. Energies of the classical barrier maximum relative to the energy of products at their equilibrium geometries.

ference in zero-point energies between products and reactants is added to this reaction energy, the heat of reaction at 0 K results. Since the HF/BAC-MP4 gives good agreement with experiment for the heat of reaction, the values for the classical reaction energy using HF/BAC-MP4 results with the frequencies scaled give the best estimate for this quantity. The MP2/BAC-MP4 parameters are also fitted to reproduce the heat of reaction and thus it is also in good agreement. The H F results with the 6-31G* and 6-31G** basis sets are in very poor agreement for the reaction energy. Including electron correlation by the MP2 and MP4 methods improves the results but it is evident that the improvement is not monotonic-the MP2 results give a better estimate of the reaction energy. The CASSCF and CAS-MRCI results have no empirical corrections and it is encouraging that at the highest level of theory the reaction energy is reproduced to within 2.4 kcal/mol. The last two columns of Table V are the values of V M E ~ (atS ) the saddle point relative to reactants and products, respectively. Both the H F results and the MP2 and MP4 results are obtained for the saddle point geometry calculated at the H F level with the 6-31G* basis. The H F calculations give barrier heights that are much too high to be consistent with experiment. The MP2 and MP4 methods reduce the barrier height, but given that the error in the MP4 reaction energy is over 4 kcal/mol, the uncertainty in the barrier height is at least that large. All of the CAS results for the barrier heights are calculated at the geometry of the CASSCF calculation with the pVDZ basis. The comparison of the CAS results with the pVDZ and pVTZ basis set indicate that the basis set errors are about 1.5 kcal/mol. The MRCI calculations significantly lower the calculated barrier height; however, as mentioned earlier, the MRCI wave functions have been shown to agree well with full-CI calculations in the same basis set. Thus we do not expect any further significant changes of the barrier height for higher order calculations. The classical barrier heights for the HF/BAC-MP4 calculations are computed at the H F saddle point even though this does not correspond to the maximum of V M E p ( s ) along the reaction path. The maximum occurs about 0.1 a, on the product side to the saddle point and the value of VMEP(s) is 0.64 kcal/mol higher for unscaled frequencies and 0.67 kcal/mol higher for scaled frequencies. The differences compared to the CAS-MRCI results are consistent with the expected uncertainties in these methods. Potential energy curves, VMEp(s),and the change in the ground-state adiabatic potential curves, AC(s) (relative to the reactant value), for the HF/BAC-MP4 method with and without scaling of the harmonic frequencies and for the MPZ/BAC-MP4 method are given in Figure 1. The scaling of the frequencies in the HF/BAC-MP4 calculations does not change the ground-state adiabatic barrier but does change V M E p ( s ) . The greatest differences are between the HF/BAC-MP4 and MP2/BAC-MP4 results. A comparison of the rates for reaction R1 computed by the ICVT/SCSAG method using the potential information obtained by the HF/BAC-MP4 and MP2/BAC-MP4 methods is

The Journal of Physical Chemistry, Vol. 94, No. 18, 1990 7103

Rate Constants for Gas-Phase Chemical Reactions

5

I

',/

"MEP

>'

1

-0.5

0.0

0.5

1.o

I

s(a o)

Figure 1. Comparison of potential energy VMep(s) along the minimum energy path and change in the ground-state adiabatic potential ALf(s) as a function of reaction coordinate s for the reaction H + NH, H2 + NH2. The zero of energy is the energy of the reactant at their equilibrium geometries and ALf(s) is the difference in the value of the ground-state adiabatic potential at location s along the reaction path and the reactant's value (at s = --). The solid curves are computed by using the HF/BAC-MP4 method with no scaling of the harmonic frequencies, whereas the long-dashed curves are computed by using the HF/BACMP4 method with the harmonic frequencies scaled as described in section 2.3. For the ground-state adiabatic potential, these two curves coincide. The short-dashed curves are computed by using the MP2/BAC-MP4 method.

-

-

TABLE VI: Comparison of ICVT/SCSAG Reaction Rates (em3 molecule-' s-') for the Reaction H + NH3 H2 + NH, Computed by Using Three Different Sets of Information about the Potential Energy Surface H F/ BAC-MP4 scaled MP2/BAC-MP4 T,K unscaled frequencies frequencies 200 9.8 x 10-24 1.1 x 10-2, 5.5 x 10-24 300 1.8 x 10-20 2.1 x 10-20 I .2 x 10-20 400 2.0 x 10-18 2.3 X 1.5 X 600 4.5 x 10-16 5.1 X 3.9 x 10-16 IO00 6.2 x 10-14 7.5 x 1044 6.1 x 10-14 1.3 X 1.1 x 10-12 I500 1 . 1 x 10-12 2400 1.3 X IO-" 1.7 X IO-" 1.3 X IO-" given in Table VI and Figure 2. The potential along the minimum energy path and the adiabatic potential curves show very similar shapes for the two methods, although the MP2/BAC-MP4 results are shifted slightly to the left. The classical barrier heights are very similar, 17.1 kcal/mol for HF/BAC-MP4 (without scaling the frequencies) compared to 16.7 kcal/mol for MP2/BAC-MP4, and the adiabatic barriers are 15.5 and 15.7 kcal/mol above the reactant zero-point energies for these respective methods. The rate constants computed by using this potential information agree well for temperatures above 300 K. Scaling of the harmonic frequencies changes the rates only slightly, by 12% at 200 K, 17% at 300 K, and up to 3 1 % at 2400 K. Larger differences occur between the HF/BAC-MP4 and MP2/BAC-MP4 results. At 300 K and below, tunneling is very important and the rates are sensitive to finer details of the adiabatic potential curves. At 300 K, the rates computed by using the HF/BAC-MP4 and MP2/ BAC-MP4 sets of potential information differ by 75% and at 200 K the difference increases to a factor of 2. These results are very encouraging for obtaining accurate rates for a modest computing effort at all temperatures except those where very deep tunneling is important.

2.5

-

3.0

3.5

4.0

4.5

5.0

loOO/T (K')

Figure 2. Rate constants as a function of temperature for the reaction H + NH3 H2 + NH2 computed by using the improved canonical variational theory with small-curvature semiclassical adiabatic groundstate transmission coefficient (ICVT/SCSAG). Temperature ranges from 417 to 2500 K and from 196 to 435 K are shown in (a) and (b), respectively. The ICVT/SCSAG rate constants are computed from potential information obtained by using the HF/BAC-MP4 method with no scaling of the harmonic frequencies (solid curves), the HF/BAC-MP4 method with the harmonic frequencies scaled as described in section 2.3 (long-dashed curves), and the MPZ/BAC-MP4 method (short-dashed curves). TABLE VII: Comparison of Three Levels of Variational Transition-StateTheory Rate Constants (cm3 molecule-' s-I) for Reaction R1 Computed by Using HF/BAC-MP4 Potential Information with HF Frequencies Scaled by 0.893 T, K TST ICVT ICVT/SCSAG 1.1 x 10-27 1.1 x 10-23 200 1.4 x 10-26 4.2 x 2.1 X 300 2.3 X 400 9.7 X 2.8 X 2.3 X 2.1 x 5.1 x 600 4.7 x IO+ 5.4 x 10-14 7.5 x 10-14 1000 9.2 x 10-14 1.3 X 1500 1.7 X lo-'' 1.2 X 1.7 X IO-" 2400 2.2 X IO-" 1.6 X IO-" In Table VI1 a comparison is made between the different dynamical methods of computing the rate constants for reaction R1 using the potential information obtained by the HF/BAC-MP4

1104

The Journal of Physical Chemistry, Vol. 94, No. 18. 1990

1 .o

0.5

-

Garrett et al.

I

I .5

2.0

lOo01T (K-')

Figure 3. Comparison of computed and experimental rate constants for

the reaction H + NH, H2 + NH2. The experimental results are from Michael, Sutherland, and Klemm' ( 0 ) ;KO, Marshall, and Fontijn) (triangles);Hack, Rouveirolles, and Wagner4 (squares). The solid curve is the result of the improved canonical variational theory with smallcurvature semiclassical adiabatic ground-state transmission coefficient (ICVT/SCSAG), using the potential information obtained from the HF/BAC-MP4 method with HF frequencies scaled by 0.893. The dashed curve is the result of ICVT/SCSAG calculations using the HF/BAC-MP4 (HF frequencies scaled by 0.893) in which I'MEp(!) was scaled so that the potential at the maximum of the ground-state adiabatic barrier was lowered by I kcal/mol.

method with the H F frequencies scaled by 0.893. The results of conventional transition-state theory (denoted TST) and the improved canonical variational theory (ICVT) neglect quantum mechanical tunneling effects, whereas the improved canonical variational theory with a small-curvature semiclassical adiabatic ground-state transmission coefficient (ICVT/SCSAG) includes tunneling within the vibrationally and rotationally adiabatic approximation. The TST calculations are carried out as previo u ~ l y , ' ~using ' ~ ~ ~the ' potential information from the BAC-MP4 calculations obtained at the H F saddle point. As described earlier, the HF/BAC-MP4 energy along the minimum energy path is not a maximum at the H F saddle point; therefore, it is not surprising that the lCVT results are much smaller than the TST ones. At 1000 K the TST results are higher by about 70% and at 300 K they are higher by about a factor of 5.5. Tunneling is very important for this reaction, contributing 98.7, 91.9, 71.7, and 42.9% to the thermal reaction rate at 300, 400, 600, and 1000 K, respectively, by the SCSAG method. A large enhancement in the rate is observed when tunneling is included-the ICVT/ SCSAG rates are higher than the ICVT ones by factors of 48.8, 8.2, 2.5, and 1.4 at 300, 400, 600, and 1000 K, respectively. The reaction rates for reaction R1 computed by the ICVT/ SCSAG method using potential information obtained by the HF/MP4-BAC method with the H F frequencies scaled by 0.893 are compared with the experimental reaction rates'-4 in Figure 3. The overall agreement is reasonable; at the lowest temperature near 500 K the error is about a factor of 2. Arrhenius fits to the experimental results give activation energies of 16.0 f 0.2 kcal/mol for the temperature range from 908 to 1777 K,' 13.8 kcal/mol for the temperature range from 490 to 960 K,3 and 14.6 f 1 .O kcal/mol for the temperature range from 673 to 1003 K.4 These can be compared with computed values of 17.4, 14.3, and 15.2 kcal/mol for the same temperature ranges. The theoretical results systematically underestimate the experimental reaction rate and overestimate the activation energies by 0.5-1.4 kcal/mol. While this is within the uncertainty of f 2 kcal/mol in the calculated potential information, a better fit of theory to experiment may be achieved by scaling VMEP(s)by a factor of 0.94 to lower its

-

Figure 4. Comparison of computed and experimental rate constants for the reaction H, + NH, H + NH,. The experimental results are from Hack, Rouveirolles, and Wagner4 ( 0 ) ;Sutherland and Michael9 (triangles); and the Arrhenius fit to the results of Demissy and Lesclaux' (long-dashed curve). The solid line is the result of the improved canonical variational theory with small-curvature semiclassical adiabatic groundstate transmission coefficient (ICVT/SCSAG), using the potential information obtained from the HF/BAC-MP4 method with HF frequencies scaled by 0.893. The short-dashed curve is the result of ICVT/ SCSAG calculations using the HF/BAC-MP4 (HF frequencies scaled by 0.893) in which was scaled so that the potential at the maximum of the ground-state adiabatic barrier was lowered by 1 kcal/ mol.

maximum value by about 1 kcal/mol. The resulting computed rates shown by the dashed curve in Figure 3 are in very good agreement with the experiment with possibly a slight overestimate at the highest temperatures. The computed activation energies for the same temperature ranges listed above change to 16.5, 13.5, and 14.3 kcal/mol, also in very good agreement with experiment. Figure 4 compares the reaction rates computed by the ICVT/SCSAG method with the HF/MP4-BAC potential information (with the H F frequencies scaled) for reaction R2 with the experimental ones of Hack, Rouveirolles, and Wagner4 and Sutherland and Michael9 and with the Arrhenius fit of Demissy and Lesclaux' given by

kFL(r ) = 2 . 2 x I O-I2 exp(-4300/ r ) cm3 molecule-' s-' (400 K I T 5 520 K) The agreement with the low-temperature results of Demissy and Lesclaux is fairly good (differences of less than 50% over the temperature range from 400 to 520 K) but the computed activation energy of 7.5 kcal/mol for the range 400-520 K is lower than the experimental value of 8.6 f 0.5 kcal/mol. For temperatures from 673 to 1003 K the computed rates are lower than the measurements of Hack, Rouveirolles, and Wagner by as much as a factor of 2.4, and the computed activation energy of 10.0 kcal/mol is higher than the experimental value of 9.1 f 0.7 kcal/mol. For the temperature range from 900 to 1620 K, the computed rates are factors of 1.6 to 2.5 lower than the results of Sutherland and Michael, and the computed activation energy of 12.3 kcal/mol is higher than the experimental value of 1 1.6 f 0.4 kcal/mol. To be consistent with the 1 kcal/mol adjustment of the barrier for the N H 3 H reaction, we have also scaled VMEp(s)for the NH, H, reaction to lower it by 1 kcal/mol. In doing so the computed results are now in excellent agreement with those of Hack et al. but still systematically underestimate the results of Sutherland and Michael by as much as a factor of 1.7 and greatly overestimate the results of Demissy and Lesclaux at lower temperatures. Lowering the classical barrier also lowers

+

+

The Journal of Physical Chemistry, V o l . 94, No. 18, 1990 7105

Rate Constants for Gas-Phase Chemical Reactions TABLE VIII: ICVT/SCSAC Thermal Rate Constants for Four Reactions Computed by Using HF/BAC-MP4 Potential Information with the HF Frequencies Scaled by 0.893 and with the Potential along the MEP Scaled To Reduce the Barrier by 1.0 kcal/mol

T,K 200 300

400 600 1000

1500

2400

+

H NH, 3.6 x 10-23 6.4 X 6.2 X 1.1 x 10-15 1.2 X 1.9 x 10-12 2.1 X IO-"

H 2 + NH2 4.7 x 10-19 1.5 X lo-!' 1.7 X 3.5 x 10-15 6.8 X 4.7 x 10-13 3.3 X

D+ND3

3.7 X 8.1 X 2.9 x 10-16 5.5 X lO-I4 1.1 x 10-12

D2+ 1.7 x 1.9 X 3.8 X 1.2 x 3.3 X 2.8 x

1.4 X IO-"

2.2

2.9 x 10-25

ND2 10-20

1O-I' 10-15 1O-I' 10-13

X

the computed activation energies to 6.9, 9.2, and 11.4 kcal/mol for the temperature ranges 400-520, 673-1003, and 900-1620 K, respectively. This improves the agreement with Hack, Rouveirolles, and Wagner and Sutherland and Michael but worsens the agreement with Demissy and Lesclaux. The rates of reaction R2 are related to those of reaction R1 by detailed balance. Knowledge of the equilibrium constant as a function of temperature allows calculation of rates for one of the reactions from rates for the other; this approach was used in ref 9 and in the present work to derive rates for reaction R2 from those for R1 and the equilibrium constants. The computed rates obey detailed balance; therefore, the excellent agreement of the computed rates with those of Michael, Sutherland, and Klemm' for reaction RI along with the underestimate of the computed rates for reaction R2 by 30% to a factor of 2 compared to Sutherland and Michael9 indicates an uncertainty in the computed equilibrium constant. As seen in Table 11, the thermochemistry of H, NH3, and H2is computed very accurately by the BAC-MP4 method, and the largest uncertainties compared with the JANAF tables are for the N H 2 radical. An error in the equilibrium constant of a factor of 2 in the temperature range from 900 to 1600 K corresponds to an error in the free energy of reaction of 1.2 to 2.2 kcal/mol; however, the heat of formation of N H 2 at 298 K is estimated to be 45.5 kcal/mol, which is only 0.5 kcal/mol lower than the computed value. Decreasing the computed heat of formation of N H 2 will shift the computed rate constant curve up in Figure 4, worsening the agreement with Demissy and Lesclaux even more. With the classical energy lowered by about 1 kcal/mol for both the forward and reverse reactions, the computed reaction rates are in good agreement with the experimental results of KO, Marshall, and Fontijn,' Michael, Sutherland, and Klemm,' and Hack, Rouveirolles, and Wagner4 for reaction RI and are consistent with those of Hack, Rouveirolles, and Wagner4 and Sutherland and Michael9 for reaction R2 but are inconsistent with the experimental results of Demissy and Lesclaux.' Our best estimate of the rates of reaction RI and R2 using the H F / BAC-MP4 potential information with I'MEP(S) scaled are provided in Table VI11 over the temperature range from 200 to 2400 K. 4.3. D ND3 D2 ND2 and D2 ND2 D ND3 Reaction Rates. The ICVT/SCSAG rate constants for reaction R3 computed using the HF/BAC-MP4 potential information (with the H F frequencies scaled by 0.893) are compared with an experimental one of Marshall and Fontijn'O in Figure 5. Once again the agreement is good but the theoretical results tend to underestimate the reaction rates. This underestimate is consistent with our observations for the H + NH, H 2 N H 2 reaction that the barrier for reactions R I and R2 is too high. To be consistent with the energetics for reactions R1 and R2, we have scaled VMEp(s) for reaction R3 so that the potential is also lowered by 1 kcal/mol at its maximum. This provides better agreement over most of the experimental temperature range except at the lowest end where the calculated results tend to overestimate the rates. The computed activation energy of 15.9 kcal/mol is smaller than the experimental value of 17.5 for the same temperature range from 591 to 1215 K. For completion, Table VI11 also presents our best estimates of the rates of reactions R3 and R4. These are computed by the ICVT/SCSAG method using potential information obtaining from HF/BAC-MP4 calculations in which VM&) is scaled to lower the barrier by 1 kcal/mol.

+

-

+

-

+

-

+

+

1CQO/T (K-l) Figure 5. Comparison of computed and experimental rate constants for

-

the reaction D + ND3 D2 + ND2. The experimental results are from Marshall and Fontijn'O ( 0 ) .The solid line is the result of the improved canonical variational theory with small-curvaturesemiclassical adiabatic ground-state transmission coefficient (ICVT/SCSAG), using the potential information obtained from the HF/BAC-MP4 method with HF frequencies scaled by 0.893. The short-dashed curve is the result of ICVT/SCSAG calculations using the HF/BAC-MP4 (HF frequencies scaled by 0.893) in which V,,,(s) was scaled so that the potential at the maximum of the ground-state adiabatic barrier was lowered by 1 kcal/ mol. 5. Conclusions Rate constants have been computed for the gas-phase chemical reactions H N H 3 2 H2 N H 2 and D + ND3 D2 ND2 over the temperature range from 200 K to 2400 K. The rates are computed by variational transition-state theory with semiclassical adiabatic ground-state transmission coefficientslSz6using limited information about the potential energy surface along the reaction path. This type of information about the potential can be obtained directly from ab initio electronic structure calculations of the energy and its first and second derivative with respect to coordinates. In the present calculations, a b initio information is empirically adjusted by the BAC-MP4 m e t h ~ d I ~ -to~ yield ' more reliable predictions of the reaction energetics. The use of variational transition-state theory in these calculations is mandated for several reasons: (1) it includes the factors most important in controlling the rate of chemical reactions and is currently the most cost effective method for obtaining reliable predictions of rates for a variety of gas-phase reactions; (2) it is capable of utilizing limited information about the potential along the minimum energy path without requiring a global potential energy surface; (3) because of the nature of the semiempirical potential information utilized in this work, the dynamical bottleneck for the reaction does not occur at the saddle point for the reaction and a variational procedure to locate the optimum dividing surface is needed; and (4) it provides a consistent method for incorporating quantum mechanical tunneling effects that are crucial for accurate predictions of the rates, especially at temperatures below 600 K. The BAC-MP4 method is a cost effective means of obtaining reliable information about the potential energy surface that is needed as input into the variational transition-state theory calculations. With current computational capabilities, it is impractical to calculate the necessary potential information from a b initio electronic structure calculations, and the use of semiempirical methods is required. The accuracy of the semiempirical potential information has been critically tested in this paper. Improving the quality of the underlying ab initio calculation from a Hartree-Fock (HF) calculation to a full second-order Mailer-Plesset

+

+

* +

J . Phys. Chem. 1990, 94. 7106-7 110

7106

perturbation theory (MP2) calculation decreased the magnitude of the empirical bond additivity correction. Comparison of the rate constants computed by using the potential information from the HF/BAC-MP4 and MP2/BAC-MP4 methods showed little difference except at temperatures below 300 K where tunneling is more sensitive to details of the potential information. The good agreement obtained here is encouraging for the applicability of the more practical HF/BAC-MP4 to generate potential information for a variety of chemical reactions. As a further test of the methods employed here, the computed rates were compared with experiment for the reactions H NH, ;2 H2 N H 2 and D ND, D2 ND,. For the reaction H NH, Hz + NH2 (RI), theory and experiment are in excellent agreement over the entire experimental temperature range from 500 K to 1770 K. The computed reaction rates for the reaction H 2 NH2 H + N H , (R2) are in very good agreement with experimental ones of Hack, Rouveirolles, and Wagner4 and Sutherland and M i ~ h a e lbut , ~ they disagree with those of Demissy and Lesclaux.' Reaction R2 is the reverse of R1 and the rates are related by detailed balance: therefore, the excellent agreement between theory and experiment for the former reaction at the lower temperatures is inconsistent with the larger discrepancies seen for the reverse reaction R2 in the same temperature range. The

-

+

+

+

+

-

+

+

-

-

computed and experimental rates for the deuterated reaction D + ND3 D, ND2 are in good agreement over the entire experimental temperature range from 590 to 1220 K. The agreement at the lower temperatures for this reaction is not quite as good as for reaction R1: the theoretical results tend to overestimate the rates slightly. The overall good agreement of the computed rate constants using different levels of theory for the potential information and the good agreement between the computed and experimental rate constants have given more confidence in the theoretical methods utilized here. Besides the comparisons with these experimental results, extensions of the rates down to 200 K and up to 2400 K are provided for reactions R I , R2, and R3 and predictions are made for reaction R4.

+

Acknowledgment. B.C.G. thanks Dr. J. V. Michael for helpful suggestions concerning the use of scaled HF frequencies to improve the agreement of the computed equilibrium constant with the experimental ones. The work at Chemical Dynamics Corporation was sponsored by SDIO/IST and managed by the Office of Naval Research under contract number N00014-87-C-0746. The work at NRL was supported by the Mechanics Division of the Office of Naval Research under contract number N0014-89-WX-24024.

Picosecond Time-Resolved Absorption and Emission Studies of the Singlet Excited States of Acenaphthyiene Anunay Samanta,*vt Chelladurai Devadoss, and Richard W. Fessenden* Radiation Laboratory and Department of Chemistry and Biochemistry, University of Notre Dame, Notre Dame, Indiana 46556 (Received: December 27, 1989; In Final Form: April 18, 1990)

Radiative and nonradiative processes in acenaphthylene have been examined by a combination of steady-state and picosecond time-resolved emission and transient absorption measurements. Fluorescence from a higher excited state (S,) is found to compete with the nonradiative decay of this state. A lifetime of