J. Phys. Chem. 1992, 96, 3662-3669
3662
-
Barrler Heights for Hydrogen Atom Transfer Reactlons. Evaluation of ab Inltlo Molecular Orbital Methods for the Degenerate Exchange H O -t H,O H,O 'OH
+
Asiri A. Nanayakkara, Gabriel G. Balint-Kurti,* School of Chemistry, University of Bristol, Bristol, BS8 1 TS, UK
and Ian H. Williams* School of Chemistry, University of Bath, Bath, BA2 7AY, UK (Received: July 31, 1991)
-
Several different ab initio molecular electronic structure techniques are used to calculate the barrier height for the reaction HO' + H 2 0 H,O + 'OH, which is a prototype for hydrogen abstraction processes involving open-shell radicals and the breaking and formation of strong covalent bonds. The objective of the study is to determine what practical methods are currently available for performing such calculations and to assess their reliability. Methods studied include self-consistent field (ROHF and UHF), multiconfiguration self-consistent field (MCSCF), Merller-Plesset perturbation calculations (up to full fourth order), configuration interaction (CISD), and quadratic configuration interaction (QCISD(T)). It is concluded that, for the calculation of transition-state properties of the important class of reactions studied here, electron correlation, size consistency, and basis set effects are all very important in reducing the computed potential energy barrier from the unrealistically large values obtained with modest levels of theory to the rather small values calculated using more expensive and sophisticated methods. The best estimate for AH* is 44 kJ mol-' at the basis set superposition error (BSSE) corrected, correlation energy scaled QCISD(T)/6-31++G(3d,2p) level. Inclusion of a tunneling correction of about -16 kJ mol-I would further reduce this value in the direction of the estimated experimental barrier height of AH* = 22 f 4 kJ mol-'.
Introduction The most important mode of chemical degradation for many organic molecules in the troposphere is hydrogen atom abstraction by hydroxyl radical HO'
HO'
+ HR
-
H2O + R'
- - -
(1)
followed by a sequence of steps
R'
R02'
RO'
carbonyl
(2)
leading to a carbonyl compound. The rationale for the use of hydrochlorofluorocarbons (HCFCs) as replacements for CFCs is that these species should undergo tropospheric oxidation by this means. Conversion from peroxy R02' to alkoxy RO' is usually effected by NO, species, but in the background middle troposphere (where the concentrations of NO, are very low), the dominant process probably involves hydroperoxy radical H02' (eqs 3 and 4) *
R02'
+ HO2'
ROOH R02' RO2'
-
+
ROOH + 0
+ HO' ROOH + R' ROOH + HOZ'
RO'
+ HR
+ HOOH
2
+
(3) (4) (5)
(6)
Hydrogen atom transfer is also important in combustion chemistry: hydroperoxides produced in H-atom abstractions by R02' (eqs 5 and 6) are responsible for the phenomenon of autoignition since their homolysis (eq 4) leads to branching in radical-chain processes. The difficulty (and expense) associated with experimental studies of such reactions, together with their intrinsic interest and fundamental nature, makes them desirable subjects for theoretical study. Our ability to compute reliably, from first principles, structural and energetic properties of molecular species involved in chemical reactions has increased dramatically over recent years (e.g., refs. 1-5). There are now several extremely comprehensive and ( I ) Pople, J. A.; Head-Gordon, M.; Fox, D. J.; Raghavachari, K.; Curtiss, L. A. J. Chem. Phys. 1989, 90, 5622. (2) Smith, B. J.; Swanton, D. J.; Pople, J. A,; Schaefer, H. F.; Radom, L. J . Chem. Phys. 1990, 92, 1240. (3) Werner, H.-J. Adu. Chem. Phys. 1987, 69,1.
generally available computer codes69 capable of carrying out sophisticated electronic structure computations on almost any chosen molecule. A key element in the immense utility of these codes has been the implementation of efficient analytic methods for the calculation of energy gradients with respect to motions of the nuclei.1° In turn this has permitted the development of procedures for the location of minima and saddle points on potential energy surfaces," which have proved invaluable in investigations of chemical reaction pathways. It is our aim to study the energetics of hydrogen atom abstractions by radicals. The composite radical-molecule systems have open-shell electronic, structures, and the reactions proceed by the concerted breaking and making of strong covalent chemical bonds to the transferring atom. In preparation for studying such reactions, we have set out to investigate the reliability and efficiency of the various readily available theoretical methods and computer codes, viz. GAUSSIAN,~ GAMESS,'J' and CADPAC? as applied to this type of problem. We have chosen to consider for this preliminary exploration the symmetrical abstraction by a hydroxy radical of a hydrogen atom from a water molecule. HO'
+ H2O
-*
H20
+ 'OH
(7)
Since there is no thermodynamic component to the energetic barrier for this degenerate reaction, the utility of the various theoretical methods depends entirely upon their ability to account (4) Truong, T. N.; Truhlar, D. G.; Baldridge, K. K.; Gordon, M. S.; Steckler, R. J. Chem. Phys. 1990, 90,7137. (5) Bowen-Jenkins, P.; Pettersson, L. G. M.; Siegbahn, P. E.M.; Almlof, J.; Taylor, P. R. J. Chem. Phys. 1988, 88, 6977. (6) Frisch, M. J.; Binkley, J. S.; Schlegel, H. B.; Raghavachari, K.; Melius, C. F.; Martin, R. L.; Stewart, J. J. P.; Bobrowicz, F. W.; Rohling, C. M.; Kahn, L. R.; Defrees, D. J.; Seeger, R.; Whiteside. R. A,; Fox, D. J.; Fluder, E.M.; Topiol, S.; Pople, J. A. GAUSSIAN 86; Carnegie-Mellon Quantum Chemistry Publishing Unit: Pittsburgh, PA, 1986; also GAUSSIAN 88; 1988. (7) Schmidt, M. W.; Boatz, J. A.; Baldridge, K. K.; Koseki, S.; Gordon, M. S.; Elbert, S. M.; Lam, B. GAMESS; QCPE Bull. 1987, 7 , 115. (8) Guest, M. F.; Sherwood, P. GAMESS Users Guide and Reference Manual; Science and Engineering Research Council: Daresbury Laboratory, Daresbury, Warrington, UK, 1990. (9) Amos, R. D.; Rice, J. E.CADPAC The Cambridge Analytic Deriud u e s Package, issue 4.0; Cambridge, UK, 1987. (10) Pulay, P. Adu. Chem. Phys. 1987, 69,241. (11) Schlegel, H. B. Adu. Chem. Phys. 1987, 67,249; New Theory and Concepts in the Understanding of Organic Reactions. NATO AS1 Ser. C 1989, 267.
0022-365419212096-3662%03.00/0 0 1992 American Chemical Society
Barrier Heights for H Atom Transfer Reactions
The Journal of Physical Chemistry, Vol. 96, No. 9, 1992 3663
satisfactorily for the bond-making and -breaking aspects of the process. There are many factors to be considered in the calculation of energy barriers to chemical reactions, and several possible ways to approach the problem. From a pragmatic point of view, the particular options available within a given suite of publicly available programs may dictate a certain course to follow. A different set of quantum chemistry codes may offer an entirely different choice of methods, and it may not initially be clear which way is best to proceed. Is it better, for example, to estimate correlation energies using variational procedures, based upon configuration interaction (CI) methods (as available in the two versions of GAMESS), or using perturbation theory techniques such as the M~rller-Plesset (MPn) method (as available in the GAUSSIAN programs)? We were unable to find in the literature an objective comparison of these alternative approaches to the computation of energy barriers for the sort of reactions studied here. Our assessment of the utility of various methods has two aspects. First, we enquire as to how reliable a particular method is for calculating the barrier to the "benchmark" reaction (7) per se. Second, we seek to establish the best viable procedure which will also be applicable to much larger molecular systems. In this study we have considered restricted open-shell and unrestricted HartreeFock (ROHF and UHF) methods, complete active space self-consistent-field (CASSCF) methods,12configuration interacti~n'~.'~ with single and double excitations from a single reference function (CISD) and from a multiple reference function (MRCISD), unrestricted Molller-Plesset perturbation theory methods with and without spin pr~jection,'~,'~ and the quadratic configuration interaction (QCISD(T)) method." We have examined the dependence of our results upon basis set effects and geometrical variation and have evaluated the basis set superposition errorI8(BSSE)and the problem of size consi~tency,'~-'~ Methods employing empirical scaling of the correlation energymV2' have also been examined. We have deliberately chosen to use Pople's widely used basis sets and standard exponents available within the GAUSSIAN program and the other computer codes employed. Our results are presented in detail and discussed in the following section.
Results and Discussion Potential Energy Surface Features of HO' + H,O System. The global minimum on the (H302)*surface corresponds to a planar C, symmetrical structure in which the hydroxyl radical donates a hydrogen bond to the water molecule. A species of this nature has been previously studied at the ROHF/3-21+G although it seems likely that Pataki et al. found not the ground state but a 2Al excited state, since its energy is 235 kJ mol-' higher than that of the separated species. Our calculations at the seven electrons in seven orbitals (7 in 7) CASSCF/6-31G level indicate that the 2B2state (out-of-plane unpaired electron on HO' moiety) is slightly lower (- 1.4 W mol-') in energy than the 2Blstate and is -34 kJ mol-I more stable than the separated HO' + H 2 0 fragments. This hydrogen-bondedspecies is not considered further in the present work. All barrier heights are computed with respect to the separated reactants HO' H 2 0 , unless otherwise stated. Leroy et al.23reported a C2, symmetrical transition structure
+
(12) Roos, B. 0. Adv. Chem. Phys. 1987,69, 399. (13) Saunders, V. R.; van Lenthe, J. H. Mol. Phys. 1983, 48,923. (14) Hehre. W. J.; Radom, L.; Schleyer, P. v. R.; Pople, J. A. Ab Initio Molecular Orbital Theory; Wiley: New York, 1986. (15) Schlegel, H. B. J . Chem. Phys. 1986, 84, 4530; 1988, 92, 3057. (16) Knowles, P. J.; Handy, N . C. J . Chem. Phys. 1988, 88, 6991. (17) Pople, J. A.; Head-Gordon, M.; Raghavachari, K. J . Chem. Phys. 1987, 87, 5968. (18) van Lenthe, J. H.; van Duijneveldt-van de Rijdt, J. G. C. M.; van Duijneveldt, F. B. Ado. Chem. Phys. 1987, 69, 567. (19) Pople, J. A.; Head-Gordon, M.; Raghavachari, K. J . Chem. Phys. 1987,87, 5968. (20) Truong, T. N.; Truhlar, D. G. J . Chem. Phys. 1990, 93, 1761. (21) Gordon, M. S.;Truhlar, D. G. J . Am. Chem. Sor. 1986, 108, 5412. (22) Pataki, L.; Mady, A.; Venter, R. D.; Pokier, R. A.; Peterson, M. R.; Csizmadia, I. G. J. Mol. Srrucr. (THEOCHEM) 1986, 135, 189. (23) Leroy, G.; Sana, M.; Tinant, A. Can. J . Chem. 1985, 63, 1447.
TABLE I: Calculated Geometries and Harmonic Frequencies for Transition Structure (HOHOH)' Using Various Theoretical Methods with the 6-31C Basis' C2 c 2 h
coordinate 0.-H 0-H O...H...O H.4-H dihedral O.-H-O-H s 0 - H stretch a 0 - H stretch a H-0-H bend s H . 4 - H bend sO.-H--Obend a O.-H.-O bend s O.-H--O stretch s torsion s O.-H-O stretch
UHF 1.178 0.963 180.0 103.3 180.0 3900 3899 1637 1095 719 390 575 341i 4765i
ROHF 1.148 0.962 180.0 104.4 180.0 3921 3912 1647 1111 704 427 672 436i 10895
CASSCF UHF ROHF 3 i n 3 7 i n 7 1.174 1.150 1.202 1.199 0.960 0.959 0.962 0.992 148.8 146.1 154.0 153.5 106.7 108.6 103.9 102.5 54.0 54.6 54.3 53.8 3932 3963 3935 3538 3931 3905 3916 3472 1487 1479 1519 1483 1553 1644 1488 1463 912 975 789 822 630 475 555 594 396 407 430 467 371 407 298 304 4657i 8015 39251 34831
a Bond lengths in angstroms, angles in degrees, and frequencies expressed as wavenumbers, in reciprocal centimeters.
Figure 1. Transition structure HOHOH'. (HOHOH)' at the UHF/4-31G level, with an imaginary frequency of 47573. cm-' for the asymmetric O*-H-*O stretching mode. Our calculationsfor this structure at the UHF/6-31G level (Table I) reveal this species to be a second-order saddle point, with another imaginary frequency (omitted from ref 23) of 341i an-'for the torsional mode. The true first-order saddle point for the symmetrical hydrogen atom transfer has C2symmetry (Figure 1 and Table I). A qualitatively correct description is obtained using either UHF or CASSCF, but the ROHF method predicts incorrect behavior for the asymmetric O-H-.O stretching mode. The exceptionally high frequency predicted for this mode by the ROHF method (bottom row of Table I) is indicative of the unsuitability of this method for describing such bond-breaking and -forming processes; a more flexible wave function is required. The simplest CASSCF calculation for the transition structure has a 3 in 3 active space which includes three electrons in the three orbitals corresponding to a simple three-center three-orbital problem. The 7 in 7 CASSCF calculation, which also includes the bonding and antibonding orbitals for the "spectator" 0-H bonds, is the smallest which allows for a rational barrier height estimate to be obtained by comparison with the 3 in 3 and 4 in 4 CASSCF calculations for the HO' and H 2 0 fragments, respectively. The ROHF and UHF methods predict very similar geometries (Table I) for the transition structure, with the exception of the partial 0-.H bond lengths. The CASSCF methods tend to overestimate the importance of bond-dissociating configurations and so tend to predict longer bonds.24 Table I1 compares calculated and experimental geometries for the HO' and H 2 0 fragments for some of the methods considered. UHF/6-31G* is economical but somewhat underestimates the bond lengths, whereas full-valence CASSCF (Le., an active space including all the orbitals of a minimal basis description but excluding the oxygen 1s core orbitals) overestimates the bond lengths. (CASSCF geometry optimizations could not be performed with basis sets including polarization functions owing to disk space limitations on the machine we used for these calculations.) The UMP2/6-311G** geometries for HO' and H 2 0 are in excellent agreement with experiment. (24) Siegbahn, P. E. M. Chem. Phys. Lerr. 1985, 119, 515.
3664 The Journal of Physical Chemistry, Vol. 96, No. 9, 1992
Nanayakkara et al.
TABLE II: Optimized Geometries' for Reactants and Transition Structure for HO' + H z 0 Hydrogen Atom Transfer UHF/6-31G* CASSCFI6-31Gb UMP2/6-31G* UMP2/6-311G**' species coordinate 0.996 HO' 0.979 0.967 0-H 0.958 0.978 0.969 0.958 0-H 0.947 H20 H-O-H 105.5 106.9 104.0 102.4 (HOHOH)' 0m-H 1.165 1.197 1.156 1.140 0-H 0.954 0.991 0.977 0.965 O...H...O 150.3 142.5 136.4 139.7 H.4-H 101.7 103.8 103.2 101.8 dihedral 55.9 57.1 61.4 59.5
expt 0.970 0.958 104.5
Bond lengths in angstroms, angles in degrees. Full-valence CASSCF, Le., oxygen 1s orbitals frozen, 7 electrons in 5 orbitals for HO', 8 electrons in 6 orbitals for H,O. 15 electrons in 11 orbitals for (HOHOH)'. 'The 6-31 1G** basis within the CADPAC program uses six Cartesian d functions, instead of five pure d functions, on the oxygen atoms. TABLE III: TbermocbemicalData Relating to Reaction Barriers for Hydrogen Atom Exchanges between Oxygen Species
AH* reaction"
HO' + HOOH H20 + H02' 0 t H20 HO' + HO' 02 + HOOH H02' + H02'
-
+
+
AHr,, -132.2 70.7 161.1
E, 7.5 76.8 178.2
AH
(expt)b (predicted)c 6.9 6 76.2 72 177.6 176
Experimental data from ref 27. Values in kilojoules per mole. AH*= 'Estimated by use of the Marcus relation, eq 8, with AH*In,= 22 kJ mol-'. a
E,
- RT.
Our objective in the present paper is, as far as possible, to examine separately each of the aspects which contributes to the computation of a barrier height for a chemical reaction. This is not strictly possible, owing to the interdependence of the various factors. Nonetheless, we have chosen to keep the geometries of the reactants and transition structure fixed for the purpose of evaluating various basis sets and methods; this is valid provided that geometrical changes caused by variation of these factors affect the reactants and the transition structure in the same way. The sensitivity of energy barriers to different choices of geometry is considered separately. It should be noted that the error involved in the use of nonoptimized geometries may be no worse than those associated with the use of truncated basis sets or of incomplete treatments of electron correlation. The problem is to know how to strike the most effective compromise of approximations among the various available options. Estimation of Experimental Barrier Height. Although no experimental value is available for the activation energy of reaction 7, it is possible to estimate this quantity by use of the Marcus relation25 AH* = u * i n t + mrxn/2 + (AHrxn)*/16AH*int (8) between enthalpies of activation AH*and of reaction AHrxnfor a series of hydrogen atom exchanges between oxygen-centered species. The quantity AH*intin eq 8 is the intrinsic barrier for a thermoneutral reaction, which is the unknown activation enthalpy for the symmetrical hydrogen atom abstraction (7) by HO' from HzO. Table I11 contains experimental data for three reactions, ranging from very exothermic to very endothermic. The AHrxn values are obtained from known heats of formation,z6 and the activation energies E, are as tabulated in ref 27. A value of A F h t = 22 kJ mol-' leads to a root mean square error of less than 3 kJ mol-' between predicted and experimental activation enthalpies, based upon the tabulated values of AHm. A reasonable estimate for the activation energy for reaction 7 is therefore 22 f 4 kJ mol-'. Depenaence of Barrier Height on Theoretical Method: VariationalMethods. The results in Table IV, for each of the four basis sets considered, show that although the reactant energies do not differ much between ROHF and UHF, the calculated transition structure energies differ substantially and cause the ROHF barriers to be much higher. We believe this difference arises because the single-determinant ROHF wave function does not (25) Marcus, R. A. J . Phys. Chem. 1968, 72, 891. (26) JANAF Thermochemical Tables, 3rd ed.;Chase, M. W., Davies, C. A,, Downey, J. R., Frurip, D. J., McDonald, R. A., Syverud, A. N., Eds.; J . Phys. Chem. Ref. Data 1985, 14 (Suppl. No. 1). (27) Westbrook, C. K.; Dryer, F.L. Combust. Sci. Terhnol. 1979,20, 125.
possess sufficient flexibility to describe the covalent bond-breaking and -formation process or the resonance between differently charged ionic structures present in the transition state. This aspect of the problem will be examined in a future publication. Unlike ROHF, the U H F wave function is not restricted to having a well-defined spin symmetry and is not an eigenfunction of (S2). It (the U H F wave function) is able to use the extra flexibility available to it to compensate partially for the inadequacy of its single-determinantal form. Both the ROHF and UHF calculated barriers are extremely high, greatly in excess of the estimated experimental activation energy of -22 kJ mol-'! Clearly these uncorrelated methods do not provide realistic descriptions of the bond-breaking and -making processes. The wave function for the CASSCF method consists of all molecular orbital (MO) configurations that it is possible to create from a specified number of electrons and MOs (the "active space"). Both the CI coefficients and the MOs are then varied together so as to minimize the energy for this wave function. By this means the nondynamical correlation energy associated with near-degeneracy effects among MO configurations is recovered. The CASSCF calculations whose results are reported in Table I1 are of the full-valence CAS variety. This corresponds to a full CI calculation for each species with a minimal basis set in which the electrons in the 1s orbitals of the oxygen atom are held frozen. The active spaces are 7 in 5 for HO*; 8 in 6 for H 2 0 , and 15 in 11 for (HOHOH)', and the initial orbitals are taken from the ROHF wave functions. The CASSCF barrier heights (Table IV) are much lower than those for the uncorrelated ROHF and UHF methods and tend toward a value of -81 kJ mol-' with increasing size of basis set. (Note, however, that here we have not used very large basis sets; see below for discussion of the basis set dependence of the barrier height. Also, the sophistication of the wave function used-Le., the size of the active spacecould always be increased.) The absolute energies for the reactants and the transition structure (Table IV) are lowest for the CISD method and are therefore "best" in a variational sense. However, the CISD barrier heights, calculated as the difference between the energy of the transition structure and thq sum of the energies of the isolated reactant species, are much larger than the CASSCF predictions. The reason for this is that the CISD method is not size consistent;28,29 the CISD energy for the combined system (HO' + H20), in which the HO' and H 2 0fragments are separated by some vary large distance, is not as low as the sum of the CISD energies of the isolated fragments. This problem is examined in greater detail below. Configuration Interaction and the Size Consistency F'roblem. In principle the CI method is the ultimate ab initio method of choice, as it leads finally to an "exact" variational result. (Other methods, notably MCSCF and coupled clusters, are also capable of yielding exact results, along with the Mdler-Plesset perturbation series when it converges.) In practice the length of the CI expansion must be truncated, commonly to include only single and double excitations out of a single reference configuration (CISD). The size consistency problem, mentioned above, does (28) Meyer, W. J . Chem. Phys. 1973, 58, 1017. (29) Pople, J. A.; Seeger, R.;Krishnan, R. Int. J . Quantum Chem. Symp.
1977, 11, 149.
The Journal of Physical Chemistry, Vol. 96, No. 9, 1992 3665
Barrier Heights for H Atom Transfer Reactions TABLE IV: TOW Energies and Barrier Heights for the HO' CASSCF/6-31G0Geometries
method ROHF
UHF
CASSCF"
CISD'
basis set 6-31'3 6-31G* 6-31G** 6-311G** 6-31G 6-31G* 6-31G** 6-31 lG** 6-31G 6-31G* 6-31G** 6-311G** 6-31G 6-31G* 6-31G** 6-311G**
+ H 2 0 Symmetrical Hydrogen Atom Transfer Computed Us*
HO' -75.361 06 -75.376 98 -75.382 78 -75.404 95 -75.36243 -75.38094 -75.386 70 -75.408 94 -75.388 76 -75.403 62 -75.408 64 -75.431 18 -75.459 00 -75.531 10 -75.541 65 -75.576 50
total energies/hartreeb H20 -75.983 49 -76.008 65 -76.020 93 -76.043 72 -75.983 49 -76.008 65 -76.020 93 -76.043 72 -76.040 30 -76.063 99 -76.074 89 -76.098 25 -76.11431 -76.19767 -76.21969 -76.255 24
(HOHOH)' -151.282 27 -1 5 1.32000 -15 1.338 87 -151.383 36 -151.300 87 -151.340 34 -151.35879 -151.40365 -151.401 01 -151.43612 -151.45284 -151.498 44 -151.538 11 -1 51.684 66 -151.717 15 -1 51.786 27
Full-Valence
barrier height/kJ mol-' 164 172 170 172 118 129 128 129 74 83 81 81 93 116 116 119
a Full-valence CASSCF, Le., oxygen 1s orbitals frozen and active spacas of 7 electrons in 5 orbitals for HO', 8 electrons in 6 orbitals for H20, and 15 electrons in 11 orbitals for (HOHOH)'. 1 hartree = 2626 kJ mol-'. 'Oxygen 1s orbitals frozen in single ROHF reference function.
TABLE V Total Energies and Barrier Heights for the HO' Methods
+ H 2 0 Symmetrical Hydrogen Atom Transfer Computed Using Miscellaneous ~
method CISD/6-31G*
geometry
HO'
CASSCF(15 in 11)/6-31G
-75.531 10 -75.536 18 -75.537 78 -75.538 08 -75.531 69 -75.53959 -75.544 92 -75.39431 -75.415 70 -75.398 14
Pople mrrectn Davidson correctn Siegbahn correctn MRCISD/6-3 1G* CISD/6-3 1+G*
CASSCF(15 in 11)/6-31G CASSCF(15 in 11)/6-31G
Pople correctn FOCI(7 in 7)/6-31(3 FOCI(7 in 7)/6-31G** SOCI(7 in 7)/6-31G a
CASSCF(7 in 7)/6-31G CASSCF(7 in 7)/6-31G CASSCF(7 in 7)/6-31G
total energies/hartree supermolecule H20 (H,O-.HO')' -76.19767 -151.71434 -76.205 29 -76.207 07 -76.207 57 -76.202 32 -76.209 15 -151.73349 -151.761 78 -76.21748 -1 5 1.437 72 -76.043 44 -76.079 11 -76.048 31
(HOHOH)' -151.68466 -151.720 61 -151.719 35 -151.723 00 -151.688 79 -151.70338 -151.73927 -15 1.405 81 -1 51.462 83
barrier height/kJ mol-' wrt wrt superseparated molecule reactants 78 116 55 67 59 119 79 119 59 61 83 83 83
Interfragment separation, 1000 A.
however constitute a very real and major difficulty with this method. The sum of the CISD/6-31G1 energies calculated for isolated HO'and H,O (Table IV) is 0.01444 hartree (38 kJ mol-') lower than the CISD/6-31G* energy calculated for a supermolecular system in which the HO' and H20fragments have the same geometries but are separated by a distance R = 1000 A (Table V). This difference in energy is referred to as the size consistency error. One way to sidestep the problem is simply to compute the bamer height as the difference between the energies of the transition structure and the R = 1000 A reactant supermolecule, such that exactly the same MO cofligurations are used in the calculations for both species. By this means the CISD/ 6-3 lG* barrier height corrected for size consistency errors may be estimated as 78 kJ mol-'. A full CI calculation would not suffer from the size consistency problem but is too computationally demanding to be generally practicable. The CASSCF method is also size consistent and may serve conveniently to provide reference functions for the determination of highly accurate wave functions. In an attempt to go some way toward the objective of a full CI treatment, some multireference CISD (MRCISD) calculations were performed using Dr. M. F. Guest's version of GAMESS.) Single and double excitations were considered out of each of the seven configurations having coefficients greater than 0.05 in the full-valence 15 in 11 CASSCF wave function. The natural orbitals from the CASSCF wave function were .employed, excluding the two oxygen core orbitals. The results of these calculations are given in Table V. Although the absolute values of the calculated MRCISD energies are lower than either the CASSCF or the single-reference CISD results, it may be seen that the bamer height differs only slightly from the CISD value and also that this MRCISD treatment is
not size consistent. Since, however, the original CASSCF wave function is actually dominated by a single configuration having a coefficient of 0.96 (correspondingto the ROHF configuration), this is not really a multireference problem at all, and the MRCISD approach was not pursued any further. Some preliminary studies were carried out using the first-order (FOCI) and second-order (SOCI) configuration interaction methods as implemented in the North Dakota State University (NDSU) version of GAMESS.' These calculations (Table V) used geometries optimized at the 7 in 7 CASSCF/6-31G level (Table I). The FOCI wave function included single excitations out of the active space of all the configurations of the 7 in 7 CASSCF wave function, resulting in 63 18 configuration-statefunctions for the transition state with the 6-31G basis. The 6-31G results indicate that the method is size consistent, but the calculated barrier height does not differ significantly from the CASSCF value. Moreover, the same barrier height was obtained with the 6-31G** basis, from a wave function involving 17640 CSFs for the transition state. The SOCI method, involving all single and double excitations out of the active space of all configurations of the CASSCF wave function, generated 42120 CSFs for (HOHOH)', which exceeded the capacity of the machine available to us. Despite the size of these calculations, they were limited by not including excitationsout of CASSCF "core" orbitals, which are in this case chemically important. In view of the earlier conclusion that this is not essentially a multireference problem, it was not considered fruitful to explore this type of wave function further. Several approximate methods have been proposed to correct CISD energies for the effect of higher excitations. The methods due to P ~ p l e ?D~a ~ i d s o nand , ~ ~Siegbahn3' as implemented in M.
3666 The Journal of Physical Chemistry, Vol. 96, No. 9, 1992 TABLE VI: Total Energies m d Barrier Heights for the HO' Quadratic CI Method
run 1
basis 6-31 1G** (6D)
geometry UMP2/6-311G** (6D)
2
6-311G1*(5D)
UMP2/6-311G**(6D)
3
6-311G**(5D)
CASSCF(l5 in 11)/6-31G
4
6-31+G*
CASSCF(15 in 11)/6-31G
5
6-31+G*
UMP2/6-3 1G*
6
6-31+Gt
UMP2/6-311G**(6D)
7
6-311+G**(5D)
UMP2/6-31 1G**(6D)
Nanayakkara et al.
+ HzO Symmetrical Hydrogen Atom Transfer Computed Using MeUer-Plesset pod method"
total energies/hartree HO' H,O (HOHOHY -75.41059 -76.046 31 -151.407 13 -75.413 23 -76.046 3 1 -151.41855 -75.597 40 -151.869 16 -76.288 61 -75.598 62 -1 5 1.876 89 -76.288 61 -75.598 91 -76.288 61 -1 51.876 75 -76.046 22 -75.4 10 49 -1 51.406 89 -76.046 22 -75.413 10 -151.41828 -75.591 40 -76.282 89 -151.857 23 -76.282 89 -75.59290 -15 1.864 74 -75.57290 -76.263 97 -15 1.819 76 -76.263 97 -75.574 39 -151.827 27 -75.408 84 -76.043 96 -151.403 47 -76.043 96 -75.411 51 -151.416 73 -76.262 49 -75.572 12 -151.814 36 -76.26249 -75.573 70 -15 1.823 66 -76.26673 -75.58403 - 15 1.828 40 -75.585 98 -76.269 89 -151.837 99 -76.274 94 -75.588 62 -1 51.849 78 -75.531 17 -76.209 60 - 151.724 98 -75.541 94 -76.213 52 -151.729 97 -76.217 47 -75.544 52 -151.740 86 -76.220 24 -75.546 20 -151.748 89 -76.220 4 1 -75.546 30 -1 51.747 93 -75.531 32 -76.209 70 -1 51.725 57 -75.546 23 -76.220 31 -1 51.748 47 -75.546 33 -76.220 48 -151.747 56 -75.531 12 -76.209 14 -151.724 22 -75.545 95 -76.2 19 69 -151.746 90 -75.54604 -76.21986 -151.745 92 -75.581 36 -76.274 68 -151.841 60 -75.596 38 -76.287 00 -151.86758 -75.596 39 -151.865 60 -76.286 58 ~
UHF PUHF UMP2( full) PUMP2(full)(KH) PUMP2(full)(S) UHF PUHF UMP2(full) PUMP2(full)(S) UMP2(fc) PUMP2(fc)(S) UHF PUHF UMP2 PUMP2 PUMP3 PUMP4(SDQ) PUMP4(SDTQ) PUMP2 PUMP3 PUMP4(SDQ) PUMP4(SDTQ) QCISD(T) PUMP2 PUMP4(SDTQ) QCISD(T) PUMP2 PUMP4(SDTQ) QCISD(T) PUMP2 PUMP4(SDTQ) QCISD(T)
~~~~
barrier ht/kJ mol-' 131 108 44 27 28 132 109 45 29 45 29 130 102 53 33 59 46 36 41 67 55 46 49 41 47 51 42 49 52 38 41 46
a Post-SCF calculations all employ the frozen-core (fc) approximation unless indicated as "full"; spin projection is by the Schlegel method unless otherwise indicated.
F. Guest's version of GAMESS were applied to the separated reactants and the transition structure at the CISD/6-3 1G* level (Table V). The improvement in the energy of (HOHOH)' relative to HO' and H20causes a reduction in the barrier height of between 49 and 61 kJ mol-I. The resulting barriers are significantly lower than that arising from treating the reactants as a supermolecule. Applying the Pople correction to each of the structures at the CISD/6-31+GS level yields similar values for the barrier when it is calculated with respect to either the separated reactants or the R = lo00 A supermolecule. DelBene and Shavitt found the Pople formula to reduce the size consistency error in CISD calculations of binding energies for hydrogen-bonded complexes to about 4 or 5 kJ mol-' and reported good results from the supermolecule approach even without the Pople correction.32 According to Hamilton and S ~ h a e f e rthe , ~ ~discrepancy between barrier heights calculated at the CISD/DZP level either by the supermolecule approach or with use of the Davidson correction is 6 kJ mol-' for the very endoergic hydrogen atom transfer from CH4 to O,, but 15 kJ mol-l for the very exoergic reverse process.
Dependence of Barrier Height on Theoretical Method: Meller-Pleaset Perturbation Theory Methods. An alternative approach to the computation of electron correlation energy is t o use the Mdler-Plesset perturbation theory34as available in the GAuSSIAN~programs (to full fourth order) and in CADPAC' (to third order). This method has the advantage of being size consistent but is nonvariational. The convergence of the total energy with increasing order of the perturbation series is often poor when using a UHF reference wave function. The value of (S2) for the (30) Langhoff, S.R.; Davidson, E.R.Int. J . Quantum Chem. 1974,8,61. (31) Blomberg, M.R. A,; Siegbahn, P. E. M. J . Chem. Phys. 1983, 78, 5682. (32) DelBene, J. E.; Shavitt, I. Int. J. Quantum Chem., Quantum Chem. Symp. 1989, 23,445. (33) Hamilton, T. P.; Schaefer, H. F. J . Chem. Phys. 1989, 90, 6391. (34) M~ller,C.; Plesset, M. S. Phys. Reu. 1934, 46, 618.
UHF wave function of the transition state (HOHOH)' with the various basis sets employed was typically 0.81, which could be reduced to the correct value of 0.75 for this doublet species by spin projection. Use of the projected UHF method (PUHF) has a substantial effect upon the calculated barrier height, lowering it by over 20 kJ mol-' (Table VI, run 1). Inclusion of electron correlation at the level of second-order Merller-Plesset theory has a dramatic effect in that it reduces the computed barrier heights by as much as 87 kJ mol-I, or slightly less after spin projection, with a flexible basis set (runs 1-3). Two different schemes have been proposed for spin-projected MP2 calculations, one by SchlegelI5(referred to as PMP2 in GAUSSIAN,but denoted here as PUMP2(S)) and the other by Knowles and Handyi6 (referred to as PUMP2 in CADPAC, but denoted here as PUMPZ(KH); our calculations yield very similar barrier heights using both methods (see run 1). By default the 6-311G** basis in GAMESS and CADPAC,with six Cartesian d functions on oxygen, is slightly larger than the default 6-311G** basis in the GAUSSIAN programs, with five pure d functions, and gives a slightly lower barrier height at the spin-projectedMP2 level (compare runs 1 and 2). Run 2 also demonstrates that the frozen-core approximation has an entirely neglible effect upon the calculated bamer heights. The third-order PUMP3 term raises the barrier by 26 kJ mol-', but the fourthorder PUMP4(SDQ) correction-which includes the effects of single, double and quadruple excitations-brings the barrier back down by 12-13 kJ mol-' (runs 3 and 4). Inclusion of the triple excitations at fourth order reduces the barrier height by a further 9-10 kJ mol-' (runs 3 and 4) and clearly has a very important effect. The PUMP4(SDTQ) calculationsyield results rather close to (within 5 kJ mol-') the values obtained at the PUMP2(S) level (runs 3-7), although computation of the full fourth-order correction takes 100 times longer than of only the second-order correction. The total energy for (HOHOH)' at the PUMP2(S)/6-31+G1 level (run 4) is significantly lower than that from the CISD
-
The Journal of Physical Chemistry, Vol. 96, No. 9, 1992 3667
Barrier Heights for H Atom Transfer Reactions
TABLE M: Basii Set Superposition Errors and Barrier Heights for the HO' + H1O Symmetrical Hydrogen Atom Trrmfer Computed Us@ sc;rlcd MsUer-Pksret and Ourdrstic CI Methods (GeometriesOptimized Usiog UMP2/6-31lCt* Method) total energies/hartree
basis 6-31G'
method PUHF PUMP2 6-31 1G** PUHF PUMP2(S) QCISD(T) 6-31 lG(Zd,p) PUHF PUMPS@) 6-311G(3d,Zp) PUHF PUMP2(S) 6-31+G* PUHF PUMP2(S) 6-31 l+G** PUHF PUMPZ(S) QCISD(T) 6-31 1++G(3d,2p) PUHF PUMP2(S)
HO' with ghost distorted as in orbitals (HOHOH)' ofH,O -75.38490 -75.38675 1-75.52434 -75.52857)'
-75.41723 -75.581 35 -75.596 37 -75.42054 -75.60200
-75.418 1 1 -75.58546 -75.60092 -75.421 31 -75.604 52
HZO barrier height/kJ mol-' distorted with ghost B ~ ~ E / k J uncorrected BSSE corrected as in orbitals (HOHOH)' ofHO' (HOHOH)' mol-' unscaled scaled unscaled scaled -75.98377 -75.98668 -151.35340 12.3 110 122 -76.17785 -76.18408 -151.70590 27.5 50 23 77 60 -151.41828 108 -1 5 1.827 27 29 14 -151.85087 39 22 -151.421 16 108 -151.859 17 26 12 -1 51.42524 110 -151.87415 27 20 -151.361 80 117 -151.724 22 42 11 -76.02368 -76.02451 -151.42558 4.5 116 120 -76.25284 -76.255 31 -151.841 60 17.3 38 29 55 53 -76.265 42 -76.267 96 -151.865 60 18.6 46 33 64 60 -76.027 52 -76.02806 -151.432 53 3.4 117 120 -76.277 78 -76.279 27 -151.88760 10.5 36 32 47 47
'PUMPZ(fu1l) results obtained using Knowles and Handy's method in CADPAC as opposed to PMP2(fc) with Schlegel's method in GAUSSIANSS.
calculation with the same basis (see Table V), but the MallerPlesset results are not variational. However, the Pople size consistency corrected energies for the transition structure and for the supermolecule (Table V) are closer to the PUMP4(SDQ) results (Table VI, run 4). Quadratic Configuration Interaction. The full fourth-order MP4(SDTQ) method is generally considered to be roughly comparable with the coupled cluster method with single, double and triple excitations (CCSDT), for systems well described by a single reference f ~ n c t i o n . ~Coupled ~ * ~ cluster theories are size consistent (but nonvariational) methods which consider certain classes of excitation to all orders in perturbation theory, whereas MallerPlesset (many-body) perturbation theory considers all classes of excitation to a certain order only. Pople et al." have formulated a methodso-called quadratic CI-which may be viewed as an approximate coupled cluster technique and is available in the GAUSSIAN programs. The QCISD method includes terms to correct for the lack of size consistency of CISD, and the QCISD(T) method adds further terms to approximate the effects of triple substitutions. Bamer heights calculated with use of the QCISD(T) method (Table VI, runs 4-7) appear to be consistently higher than the Maller-Plesset PUMP4(SDTQ) results by about 3-5 kJ mol-'. However, it is worth noting that the computation times for the QCISD(T) energy were 10 times longer than for the complete PUMP4 correction. With the largest basis used, 6-31l+G**, the calculated QCISD(T) potential energy barrier height is 46 kJ mol-'. Dependence of Barrier Height upon Geometry. Comparison of runs 2 and 3, and of runs 4-6, in Table VI shows that geometrical variations between structures optimized at the CASSCF( 15 in 11)/6-31G, UMP2/6-31GS, or UMP2/6-311G** level are of rather little consequence, leading to differences in calculated barrier heights of no greater than 3 kJ mol-'. It therefore seems quite reasonable to determine geometries at a lower level of theory than that used for evaluation of the barrier height, provided that any important effects of electron correlation upon geometries have been included. UMP2 is adequate for this purpose. Depenaence of Barrier Height on Basis Set. The results for the four variational methods listed in Table IV show an increase in the barrier height arising from inclusion of polarization functions on oxygen in 6-31G* (and larger bases) as compared with 6-31G. Further inclusion of polarization functions on hydrogen (6-3 1G**) or the addition of an extra set of s and p functions in the valence
-
(35) Urban, M.; Dicrckscn, G. H. F.; Cemusak, I.; Havlas, Z. Chem. Phys.
shell (6-31 lG**), however, appears to have little effect. Owing to the incomplete nature of these basis sets, the energy of each of the reactant fragments HO' and H 2 0is lowered by the presence of the orbitals of the other fragment in the (HOHOH)' complex. This spurious stabilizationis called the basis set superposition error (BSSE) and may be eliminated using the Boys-Bernardi meth0d.'*9~' The BSSE for each fragment is calculated as the difference of two energies: first, the energy of a fragment with the geometry it has in the transition structure and, second, the energy of this distorted fragment together with the "ghost" basis functions of the other (distorted) fragment as it is in the transition structure but without its nuclei and electrons. This procedure yields a BSSE energy of 11 kJ mol-' at the ROHF/6-3 1G** level. Thus the energy of the transition structure-and the calculated barrier height-is artificially too low by this amount. The variational methods listed in Table IV do not provide treatments of electron correlation effects as good as those of the nonvariational methods listed in Table VI. Comparison of runs 3 and 4, and of runs 6 and 7, in Table VI suggests that increasing the size of the basis set from 6-31+G* to 6-31 1G** or 6-31 1+G** tends to decrease the barrier height as the description of the transition structure is improved. The BSSE energies at the UHF, PMP2, and QCISD(T)/6-311+G** levels are about 4, 17, and 18 kJ mol-', respectively, suggesting the bamer is underestimated by these amounts. Basis set improvements have two opposing effects. On the one hand, the basis for each constituent part of the supermolecule is made more complete, thereby in general (though not always) reducing the BSSE and tending to increase the calculated bamer height. On the other hand, lower calculated energies are obtained for both the supermolecule and for its isolated constituents. On the whole the energy of the supermolecule, the transition structure HOHOH', is lowered more than the sum of the energies of the isolated constituents, HO' and H 2 0 , so tending to decrease the barrier height. These trends are illustrated by the results in Table VII, which contains barrier heights (uncorrected and unscaled) calculated using a variety of basis sets but in each case using the same UMP2/6-3 llG** optimized geometries. At the PUHF level, increasing the number of polarization functions (d on 0,p on H) has little effect upon the barrier, but inclusion of a set of diffuse s and p functions on 0 raises the barrier by 6 or 7 kJ mol-'. At the PUMP2 level, improvement of the basis set first reduces the calculated barrier height (6-31G*, 50 kJ mol-'; 6-311G**, 29 kJ mol-'; 6-31 lG(2d,p), 26 kJ mol-'). When diffuse functions on 0 are introduced, an increase of -9 kJ mol-' in the barrier results,
Lett. 1989, 159, 155.
(36) Stanton,J. F.; Bartlett, R. J.; Magers, D. H.; Lipscomb, W. N. Chem. Phys. Lett. 1989, 163, 333.
(37) Boys, S. F.; Bernardi, F. Mol. Phys. 1970, 19, 553.
3668 The Journal of Physical Chemistry, Vol. 96, No. 9, 1992
Nanayakkara et al.
TABLE VIII: Basis Set Superpasition Errors and Empirical Scaling Factors for Bod Dissociation Energies D(H0-H)Computed Using Scaled M~dler-Pleset and Quadratic CI Methods (Geometries Optimized Using UMP2/6-311G** Method) total enerniedhartree H' HO' D(HO-H)/kJ scaling factor with with ghost ghost mol-' 9 orbitals distorted orbitals BSSE/ BSSE BSSE basis method" isolated of HO' isolated as in H 2 0 of H' H20 kJ mol-' unwrr corr uncorr corr 6-31G* PUHF -0.49823 -0.498 33 -75.38486 -75.38492 -75.386 11 -76.01027 3.4 334.0 337.4 PUMP2(S) -0.498 23 -0.498 33 -75.52090 -75.522 10 -75.525 23 -76.196 52 8.5 466.8 474.3 0.692 0.726 6-31 lG** PUHF -0.49981 -75.413 10 -76.046 22 350.1 PUMP2(S) -0.499 81 -75.574 39 -76.263 97 498.3 0.843 QCISD(T) -0.499 81 -75.589 25 -76.276 33 491.8 0.806 6-31 1G(2d,p) PUHF -0,49981 -75.41472 -76.047 45 349.0 PUMP2(S) -0.499 81 -75.589 52 -76.279 53 499.5 0.851 6-31 1G(3d,2p) PUHF -0.49982 -75.41636 -76.05092 353.8 PUMP2(S) -0.499 82 -75.595 07 -76.289 55 518.7 0.958 6-31+G* PUHF -0.49823 -75.389 42 -76.01709 339.9 PUMP2(S) -0.498 23 -75.531 12 -76.209 14 472.1 0.71 1 6-311+G'* PUHF -0,49981 -0.49981 -75.417 18 -75.41735 -75.41827 -76.05242 2.4 355.6 358.0 PUMP2(S) -0.49981 -0.49981 -75.581 36 -75.581 24 -75.58602 -76.27468 12.5 508.2 520.7 0.896 0.969 QCISD(T) -0.49981 -0.49981 -75.59639 -75.59622 -75.601 45 -76.286 58 13.7 499.9 513.6 0.847 0.927 6-311++G(3d,2p) PUHF -0,49982 -0.49982 -75.42049 -75.42069 -75.421 35 -76.05660 1.7 357.9 359.6 PUMP2(S) -0.49982 -0.49982 -75.60201 -75.601 92 -75.60484 -76.29935 7.7 518.7 526.4 0.957 1.003 a
Experimental bond dissociation energy, 525.9 kJ mol-I.
provided the basis is already of triple-[ quality in the valence region. Diffuse functions on H appear to have no effect, but extra polarization functions lead to a gradual decrease in the barrier height. The BSSE calculated at the PUHF level diminishes (from 12.3 to 4.5 to 3.4 kJ mol-') as the size of the basis in increased from 6-31G* to 6-311+G** and then to 6-311++G(3d,2~);~~ much larger BSSEs are obtained at the PUMP2 level, but similarly decrease (from 27.6 to 17.3 to 10.5 kJ mol-I) with expansion of the basis set in the same manner. The BSSE calculated at the QCISD(T) level with the 6-311+G** basis is slightly larger than that for the PUMP2 method with the same basis, since a larger proportion of the correlation is recovered. We have followed the recommendations of van Lenthe et al.ls and of Szczesniak and S ~ h e i n e that r ~ ~ the full Boys-Bernardi counterpoise procedure should be used at both the Hartree-Fock and correlated levels of theory. Scaled Correlation Energy Estimates of the Barrier Height. Truong and Truhlar2' have employed the MP-SAC2 (scaling all correlation energy in second-order Molller-Plesset theory) methodm to determine the potential energy barrier height, and then rate constants, for the hydrogen atom transfer from CHI to HO', apparently with some success. In this approach the correlation energy is scaled (Le., divided) by an empirical factor 3, which is a ratio of differences in bond dissociation energies for the making and breaking bonds of interest:
3 (DMP~- gUHF)/(pxpt - g U H F ) (9) Following this idea, we have determined scaling factors 3 for both PUMP2 and QCISD(T) correlation energies (relative to the PUHF energy) for the 6-311G** and 6-311+G** bases, and for PUMP2 energies with the other bases already listed in Table VII. We have considered the bond dissociation energy D(HO-H) of the water molecule. These scaling factors (uncorrected) are given in Table VIII; as expected, the value of 3 tends to unity as the size of the basis increases. The barrier heights calculated for reaction 7 are scaled according to u*Scaied
= M'PUHF + ( M ' p u ~ p z- ~
* P u H F ) / (10) ~
= ~ * P U H+ F ( ~ * Q C I S D ( T )- M * P U H F ) /(11) ~ eqs 10 and 11. The scaled (uncorrected) barrier heights in Table VI1 show that the effect of scaling the correlation energies in this manner is to lower the barrier height by an amount ranging from 27 kJ mol-' with the 6-31G* basis to only 4 kJ mol-' with the u*gaied
(38) Frisch, M.J.; Pople, J. A.; Binkley, J. S. J . Chem. Phys. 1984, 80, 3265. (39) Szczesniak, M. M.; Scheiner, S. J . Chem. Phys. 1986, 84, 6328.
TABLE I X Nonpotential Energy Contributions to the Barrier Heights for the HO' + H 2 0 Symmetrical Hydrogen Atom Transfer contrib to barrier (HOHOH)' ht/kJ
CASSCF/6-3 1G
HO' H,O zero-point energy E,,/kJ mol-' thermal energy (298 K)
&/kJ
20.1
53.8
72.0
-1.9
28.8
63.7
84.5
-8.0
mol-'
largest basis employed, 6-31l++G(3d,2p). The scaled PUMP2 barriers include values as low as 11 kJ mol-l and show a surprising degree of variation with improvement of the basis. The scaled QCISD(T) barriers are 4-8 kJ mol-' higher than the scaled PUMP2 values. With the largest bases it is not practicable to perform calculations at the QCISD(T) level. These scaled correlation energy bamer heights still contain the BSSE. It seems desirable, therefore, to modify the procedure outlined above to take account of the BSSE. We have determined corrected scaling factors 3 based upon BSSE-corrected bond dissociation energies D(HO-H) for three of the bases employed, which are given in the final column of Table VIII. It is of interest to note that these values are consistently closer to unity than the corresponding uncorrected values. The BSSE-corrected, scaled correlation energy barrier heights for reaction 7 are presented in the final column of Table VII. These barriers are considerably higher than their BSSE-uncorrected counterparts. However, improvement of the basis from 6-31G* to 6-31 1+G** and then to 6-31 l++G(3d,2p) now causes a steady decrease in barrier height, to a value of 47 kJ mol-' with the largest basis. With the 6-31 1+G** basis, the BSSE-corrected, scaled barrier is 7 kJ mol-' higher at the QCISD(T) level than at the PMP2 level. The best estimate for the potential energy barrier height AE* for reaction 7 is therefore obtained by adding 7 kJ mol-' to the corrected, scaled PMP2/6-311 ++G(3d,2p) value, giving 54 kJ mol-'. The use of empirical data to correct ab initio calculations has a long h i ~ t o r y . ~ ~ The , ~ ' correlation energy scaling methods discussed above have the desirable feature, which they share with the other methods, that in the limit of an exact calculation no correction is applied to the ab initio calculation. Noapoteatial hergy Contributioap to the Barrier Height. Table IX contains zero-point energies and thermal energies (at 298 K) (40) MoMt, W.Proc. R. Soc. London 1951, A210, 245. (41) Balint-Kurti, G. G.; Karplus, M. J. Chem. Phys. 1969, 50, 478; In Orbiial Theories of Molecules and Solids; March, N. H., Ed.; Ciarendon Press: Oxford, UK, 1974.
3669
J. Phys. Chem. 1992,96, 3669-3674 calculated from the CASSCF(7 in 7)/6-31G vibrational frequencies for (HOHOH)' given in Table I and the corresponding 3 in 3 frequencies for HO' and 4 in 4 frequencies for H20. These nonpotential energy terms contribute = -2 kJ mol-' and A,!?',,, = -8 kJ mol-' to the enthalpy of activation. When added to the best estimate for the potential energy barrier, a total value for AH*298 = 44 kJ mol-' is obtained. Truhlar and Garrett42have noted that what they term as "nonsubstantial" contributions to activation enthalpies, that is, those contributions-such as quantum-mechanical tunnelingwhich are properties of the dynamics and of the potential energy surface rather than of particular molecular structures, may be very significant. In the case of the hydrogen atom transfer from H2 to HO', they estimated AH((substantia1) = 20 kJ mol-' and AH'(nonsubstantia1) = -16 kJ mol-'. If tunneling effects are similarly important in the symmetrical exchange of a hydrogen atom between H 2 0 and HO', then the estimate for the activation enthalpy would be further reduced to a value around 28 kJ mol-'.
Utp
Conclusions We have examined the use of the range of variational and nonvariational methods for electronic structure calculationsoffered by the widely available computer packages GAUSSIANEE,GAMESS, and CADPAC in order to determine the barrier height to the symmetrical hydrogen atom transfer from H20to HO'. This reaction is a prototype for a series of more complicated hydrogen atom transfers of interest to us, and it was our intention to assess the reliability of the various methods and their suitability for use with larger molecular systems. The barrier heights of reactions are of course not susceptible to direct measurement and in many cases the rate constants, from whose temperature dependence barrier heights may be estimated, are extremely difficult to obtain experimentally. Indeed, no direct experimental rate data exist for the title reaction of this study, although an estimate may be obtained by consideration of other hydrogen atom exchanges between oxygenic species. The results of this study indicate that electron correlation, size consistency, and basis set effects are all important in reducing the potential energy barrier from 129 kJ mol-' at the UHF/
-
(42) Truhlar, D.G.;Garrett, B. C. J . Am. Chem. SOC.1989,111, 1232.
6-31G* level to an estimated 54 kJ mol-' at the BSSE-corrected, correlation energy scaled QCISD(T)/6-31 ++G(3d,2p) level. Ow best estimate for the activation enthalpy AH* of reaction 7, 44 kJ mol-', is in reasonable agreement with the Marcus relation derived estimate of = 22 f 4 kJ mol-', if a further correction of reasonable magnitude arising from quantum-mechanical tunneling is assumed. While there is a question concerning the reliability of Moller-Plesset calculations based upon UHF reference wave functions, it does seem that the BSSE-corrected, correlation energy scaled PUMP2 method may yield useful results at a small fraction of the cost of more sophisticated methods. (Typically the ratio of computation times for evaluation of the QCISD(T) vs PUMP2 energy is 1OOO:l.) On the basis of our discussion above and on considerations of computational feasibility, we recommend the procedure laid out below for calculating the barrier heights of reactions involving larger or more complex molecular systems. (1) Geometries of the reactants and of the transition structure should be computed using the UMP2 method with a basis set which includes polarization functions (e.g., 6-31G*). (2) Energies should then be computed using the largest possible basis set (e.g., 6-311++G(3d,2p)) at the PUMP2 level. Corrections for the BSSE should be applied using the full BoysBernardi counterpoise procedure. Further corrections for the remaining correlation error should then be applied using the scaled correlation energy method. (3) Calculations should be performed using a smaller basis (say, 6-31+G*) with both the QCI and PUMP2 methods. The difference between the barrier heights computed by these methods should then be applied as a correction to the barrier height determined in (2); this provides an estimate for the result of a QCI computation using the largest basis set without actually expending the computational effort needed to do so.
-
Acknowledgment. We are grateful to Shell Research Ltd. for financial support and to Dr.C. Morley for many useful discussions. We thank the S.E.R.C. for the provision of computational facilities at the Rutherford Appleton Laboratory (RAL), Dr. M. F. Guest for making the GAMESS program available to use on both the CRAY and IBM computers at RAL, and Professor J. S.Francisco for performing the calculation presented as run 3 in Table V for us.
Geometries, Force Fields, and Vibrational Assignments of Dewar Benzene and Dewar Pytldlne Ruifeng Liu, Xuefeng Zbou, and Peter Pulay* Department of Chemistry and Biochemistry, University of Arkansas, Fayetteville, Arkansas 72701 (Received: September 24, 1991)
The geometries of Dewar benzene and Dewar pyridine were fully optimized by restricted Hartree-Fwk theory and the 6-31G* basis set. Reference geometries for force field calculations were obtained by systematic empirical corrections. The ab initio force fields, after empirical scaling, yield frequencies which are in good agreement with experiments. On the basis of the calculations, some uncertainties in the experimental assignments of the fundamental frequencies of Dewar benzene were resolved, and assignments for the observed matrix IR frequencies of Dewar pyridine were given for the first time. This study demonstrates again the value of the scaled quantum mechanical (SQM) force field method.