J. Phys. Chem. 1995, 99, 14689-14699
14689
A Systematic Theoretical Study of Harmonic Vibrational Frequencies and Deuterium Isotope Fractionation Factors for Small Molecules Nathan J. Harris* Department of Chemistry, Indiana University, Bloomington, Indiana 47405 Received: May 26, 1995; In Final Form: July 27, 1995@
Harmonic force fields were calculated using fourth-order Mgller-Plesset (MP4) theory for C h , NH3, H20, HF, Si&, PH3, HzS, HC1, acetylene, HCN, ethylene, formaldehyde, chloromethane, and fluoromethane. These were compared with the experimental force fields. Most of the diagonal force constants for these molecules are computed reliably (1-2% error) at the MP4/6-311+g(2d,2p) level. Exceptions are force constants for out-of-plane bending in the n-bonded molecules and for stretching of the multiple bonds. For computations using second-order Mgller-Plesset (MP2) theory, the 6-31g* basis is inferior to basis sets having polarization functions on all atoms such as 6-31l+g** and 6-31 l+g(2d,2p). The force fields computed with MP2 theory are as accurate as those computed at the MP4 level except that HX stretching force constants are consistently overestimated by ca. 3% with MP2 theory. Reduced isotopic partition function ratios, (sds~x and fractionation factors (FFs) for hydrogeddeuterium substitution were computed for the same set of molecules. The (sdsl)f values were compared with values computed from experimental harmonic frequencies. The (s2/s1)f values are overestimated by RHF and MP2 theory, but the error is removed when frequencies for the isotopomers are scaled uniformly. The resulting values agree with experiment to within f 2 . 5 % at the MP2 level and f 3 . 3 % at the RHF level.
Introduction* In the harmonic approximation the potential energy for a polyatomic molecule is expressed in the quadratic form 2V = qTq, where q is a vector of Cartesian displacements of the atoms and F is the symmetric matrix of Cartesian force constants. The elements of F are the second derivatives a2Vl aqia% of the potential energy with respect to the atomic coordinates. The harmonic potential can be improved by including higher tems in the Taylor series expansion of the potential energy. It is more convenient to use internal coordinate displacements than Cartesian displacements in dealing with harmonic vibrations. Typically, linear combinations of internal coordinates called symmetry coordinates are used. The vector of Cartesian coordinate displacements q is related to the vector of internal coordinate displacements r by a linear transformation, r = Bq. The matrix B is rectangular with dimensions (3N - 6)x(3N) and can be formed if the molecular geometry is known. Because B is rectangular it has no inverse, but the B matrix can be expanded by applying the Eckart conditions, and the expanded matrix can then be inverted.Ia The matrix of internal coordinate force constants F, can be computed from the Cartesian force constants F from F, = ATA, where A is the inverse of the expanded B matrix. This transformation between coordinate bases is only linear if the vibrational displacements are infinitesmally small, as is assumed in the harmonic approximation.2 The central problem in the theory of harmonic vibrations is setting up and solving the secular equation, GFrL = LA. Here, G is the matrix defined by Wilsonlb as the inverse of the kinetic energy matrix. F, is the force constant matrix corresponding to a chosen set of internal coordinates, L is the orthonormal matrix of eigenvectors of GFr, and A is the diagonal matrix of eigenvalues of GF,. The matrix GF, is not symmetrical, but it can be reduced to a symmetric form and then diag~nalized.~ The matrix, of internal coordinate force constants F, is also called the molecular force field. Determining the force field 'Abstract published in Advance ACS Abstracts, September 1, 1995.
0022-365419512099-14689$09.0010
from observing the vibrational spectrum of a molecule is a difficult task because of the large number of unknown constants compared to the amount of easily accessible data. It is helpful to observe the spectra for several molecules which differ only by isotopic substitution, but still it is generally not possible to determine all elements of F, without additional information. Harmonic Force Fields from ab Initio MO Theory. In the past 20 years it has been common practice to compute some or all of the force field using ab initio molecular orbital theory. With modem computer packages all the elements of the force constant matrix can be computed. A test of theoretical methods is whether they are able to reproduce experimentally determined force fields. Use of restricted Hartree-Fock (RHF) theory with a double 5 plus polarization (DZ+P) quality basis set gives diagonal harmonic force constants which are systematically too large by 10-20%. However, in many cases the off-diagonal elements of the force constant matrix can be determined more reliably from RHF theory than from experimental data.4 The above accuracy in the force constants leads to frequencies which are consistently larger than experimental frequencies by about 10%. Part of the difference comes from the effect of vibrational anharmonicity, which makes observed fundamental frequencies smaller than harmonic frequencies by about 5%.5 Even after observed frequencies are adjusted for anharmonic effects theoretical RHF frequencies are still larger by perhaps 5%. This difference arises from a deficiency in the singledeterminant potential energy surface, which predicts the energy to rise too steeply with increasing distortion from the equilibrium geometry. A related observation is that bond dissociation energies predicted from RHF theory are much too large, because the correlation correction becomes more and more important with increasing internuclear separation. This deficiency in the RHF wave function can be partly corrected by using secondorder M@ller-Plesset perturbation theory (MP2 theory) to account for electron correlation.6 MP2 theory improves the theoretical force field such that calculated frequencies agree with experimental harmonic fre0 1995 American Chemical Society
Harris
14690 J. Phys. Chem., Vol. 99, No. 40, 1995
quencies to within 2% for molecules such as methane and water. To obtain such results from MP2 theory it is necessary to use a high-quality basis set.7 Simandiras et al. have recommended the use of basis sets consisting of triple 5 plus two sets of polarization functions plus one set of f functions (TZf2P+f) for reliable calculations of harmonic frequencies for n-bonded systems with MP2 t h e ~ r y . ~Pulay ~ . ~ et al. investigated the accuracy of the ab initio evaluation of diagonal force constants of HF, HCN, and N H 3 and concluded that the required “chemical accuracy” (about 1% in the frequencies and 2% in the harmonic force constants) can be approached at the level of a TZ+2P basis with single and double substitutions in the CI treatment plus a correction for unlinked clusters.7d However, Sellers and Almlof have argued that even TZf2P-l-f basis sets are inadequate. Their conclusion was based on a study of geometries and harmonic force constants for N2, HF, F2, and 0 2 computed using a multireference configuration interaction theorye9 Harmonic vibrations of polyatomic molecules have also been studied with more sophisticated theories. Examples are studies by Stantonlo using fourth-order Moller-Plesset theory (MP4 theory) and by several different groups using coupled cluster (CC) t h e ~ r y . ~ ~Very J l J ~accurate calculations using CC theory together with large basis sets have been performed for several small molecules,I2 such as methane,12aHCN,Izasbwater,120cand ammonia.12a,d As a method of providing for electron correlation, MP2 theory is fast and efficient and can reproduce experimental data for harmonic frequencies with reasonable accuracy (ca. 2% error). MP2 theory is the most economic method for correcting the deficiencies of RHF theory, but force field calculations including MP2 theory and large basis sets are not feasible for the usual organic molecules containing more than three or four nonhydrogen atoms. A more practical method for computing accurate theoretical frequencies for large molecules is to scale force constants computed with SCF theory to fit experimental frequencies. For most larger molecules only the fundamental frequencies are known, and not the harmonic frequencies, and scaling corrects the theoretical force field for both the effects of anharmonicity and the deficiency of RHF theory. Pupyshev has concluded recently that for calculations with very large basis sets @e., at the RHF limit) the true harmonic force field for a molecule can be reproduced by applying a single scaling factor to the RHF force field.I3 In a recent study Pople found that application of a uniform scaling factor of 0.89 to frequencies computed at the RHF/6-3 lg* level reproduced experimental fundamental frequencies with a root mean square error of 50 cm-I for 122 small molecules containing first- and second-row atoms. Frequencies computed at the MP2/6-31g* level are also improved by ~ca1ing.l~ More accurate methods employ scaling of the force constants rather than the frequencies and use several scaling factors instead of only one. A widely used method was developed by Pulay,I5 in which a scaling factor c, is introduced for each diagonal force constant $, associated with symmetry coordinate Si. The offdiagonal force constants $, are then scaled by the geometric mean (c,c,)1’2 of the corresponding diagonal scaling factors. In Pulay’s scheme, separate scaling factors are introduced for each qualitatively different distortion, whereas identical scaling factors apply for qualitatively similar distortions. For the ethylene molecule, as an example, scale factors would be assigned for HCC bending, CC stretching, CH stretching, CH2 out-of-plane bending, and HCCH torsions. The scaling factors are then optimized using least-squares to give the best fit to the fundamental frequencies of a molecule or to the frequencies of
a group of related molecules. An important use for a scaled ab initio force field is to provide clear assignments of the fundamental vibrational frequencies for a molecule when more than one interpretation of spectral data is possible.16 Equilibrium Isotope Effects from Harmonic Vibrational Frequencies.I7J8 An equilibrium isotope effect (EIE) is defined as the ratio of equilibrium constants for two reactions which differ only by isotopic ~ubstitution.~~ A,
EIE = KJK2 = ~
K, = [B,l/[A,l
~
~
*
1
~
~
1
1
~
~ (1) ~ ~
In eq 1, the subscript 2 denotes substitution by a heavier isotope, such as substitution of hydrogen by deuterium or of I2Cby I3C. The EIE can also be written as the equilibrium constant for a hypothetical isotopic exchange reaction. The expression for the EIE is adjusted to account for symmetry effects by multiplying by a symmetry factor. EIE = ((SA~SB~)/(SA~SB~))([A~I[B~I)/([A~I[B~I)(2) In eq 2, SA^, SA^, SB~,and S B ~are the symmetry numbers of each species. It is convenient to define a quantity (s~/sI’)~A called the reduced isotopic partition function ratio for AI and A2 by
J
where Q2/Ql is the ratio of molecular partition functions for molecules A2 and AI, and the mj refer to masses of the atoms. Bigeleison has shown that the reduced isotopic partition function ratio can be computed in the harmonic approximation using only the normal mode frequencies for AI and Az.I8
In eq 4,ui = hoi/kT, and the product is taken over the 3N - 6 normal mode frequencies oi for A1 and Az. The EIE in eq 2 can be rewritten in terms of ($2/SI’)fAand (SZ/SI)~B as
This is the Bigeleisen equation, which shows that the EIE can be calculated from knowledge of the normal mode frequencies of each species. The [uiduill, [exp((uii - ui2)/2)1, and [(I - exp(-uii))/(l exp(-u,z))] terms in eq 4 are related respectively to mass and rotational inertia, zero-point vibrational energy, and vibrational excitation of A1 and A2. These terms will be abbreviated as R, ZP, and EX. The ZP term is the most important contributor to variations in the (s2/sl)fvalues for different molecules, and the EX term is the least important. Using the convention in which the 2 refers to the heavier isotope, the ZP term is larger than one and the R term is less than one. The ZP term grows smaller with decreasing temperature, approaching unity, while the R term stays constant. The EX term grows closer to unity as the temperature is lowered. At very high temperatures EX tends asymptotically to 1/R, and the product R x ZP x EX tends to unity. Wolfsberg and co-workers have attempted to make corrections for vibrational anharmonicity in computing the reduced
~
1
1
~
A Systematic Theoretical Study of Small Molecules isotopic partition function ratios for several small molec~Ies.~9 The major effect on the molecular partition function from anharmonicity is on the ZP term. If observed fundamental frequencies are used in place of harmonic frequencies in eq 4, (s2/sl)fA is made about 10% smaller due to a change in the ZP term. The R term in eq 4 is also affected by anharmonicity, because the Teller-Redlich product rule is only strictly correct for true harmonic vibrations. Substitution of anharmonic frequencies in place of harmonic frequencies in eq 4 makes the R term larger by about 2% for the molecules we will consider. The EX term is generally close to unity, and the choice of anharmonic or harmonic frequencies in computing (s2Isl)f has a negligible effect on EX. The effects of anharmonicity on the numerator and denominator of eq 5 cancel somewhat, so that the EIE, in contrast to (s2/s1$ is often approximately the same whether or not observed frequencies are adjusted for anharmonic effects, 9b Data for equilibrium isotope effects are sometimes tabulated in the form of fractionation factors (FFs). The FF for molecule A is the reduced isotopic partition function ratio (s~/sI)~A relative to some reference value (s2/sl)F.
Shiner and Neumann have tabulated deuterium FFs for a large number of organic molecules using acetylene as a reference.20 The FF values are useful for comparing organic molecules which differ by a substituent or functional group. In Shiner's method20ba harmonic force field is adjusted by least-squares methods to reproduce fundamental frequencies for a number of isotopic derivatives of a molecule. Frequencies computed using this force field are then used as inputs for the Bigeleisen equation. The FF values from this procedure generally agree with FFs computed from harmonic frequencies?lb and the use of observed fundamental frequencies is preferred because harmonic frequencies are only available for a few small organic molecules. The reduced isotopic partition function ratio can be calculated using harmonic frequencies computed from an ab initio force field.21 A popular Fortran program called QUIVER was developed for this If (sdsl)fis computed using RHF frequencies the value is overestimated, but use of MP2 theory with a reasonably large basis gives improved agreement with values based on experimental frequencies which have been corrected for anharmonicity. Ab initio theory has also been applied to the study of kinetic isotope effects (KIEs), using calculated frequencies of the reactants and the transition state.22 Methods Ab initio and semiempirical calculations were performed on an IBM RS/6000 computer using the Gaussian program.23 Cartesian force constants for 14 small, semirigid molecules were computed from the single-determinant (RHF) wave functions and using Moller-Plesset theory" terminated at second-, third-, or fourth-order (denoted as MP2, MP3, or MP4). RHF and MP2 force constants (Le., the second derivatives of the potential energy with respect to the nuclear coordinates) are computed analytically by the Gaussian program. Gaussian computes MP3 and MP4 force constants by double numeric differentiation of the energy, and the frozen-core approximation is applied. The MP4 calculations include the contribution from triple excitations (MP4(SDTQ)). The calculations were done using the 3-21g* basis,256-31g* basis,26and 6-31 lg** basis2' sets included in Gaussian,28which were supplemented with extra polarization functions and diffuse
J. Phys. Chem., Vol. 99, No. 40, 1995 14691 TABLE 1: Basis Sets A
B C D E
A3-21g* 6-31g* 6-31+g* 6-311g* 6-311g**
F G H I
6-311+g** 6-311+g(2d,2p) 6-311+g(2df,2p) 6-311+g(3d,2p)
J K L M
6-3 11+g(3df,2p) 6-311+g(3d,3p) 6-311++g(3d,3p) 6-311++g(3df,3p)
function^.^^^^^ To demonstrate the notation used for basis sets, the 6-31g* basis set has polarization (*) functions on all heavy atoms, the 6-31 l+g** basis has polarization functions on all atoms (d functions for heavy atoms and p functions for hydrogens) along with diffuse (+) functions of s and p symmetry on heavy atoms, and the 6-31 l+g(2df,2p) basis has two uncontracted sets of d functions and one set of f functions on heavy atoms and two sets of uncontracted p functions on hydrogens. The 3-21g* basis has d functions on second-row atoms and is otherwise identical to the 3-21g basis. In addition, force constants were computed using the semiempirical AM1 Hamiltonian30 (AM1 was replaced by MND03' for molecules containing second-row atoms). The Cartesian force constants were transformed to intemal coordinate force constants using standard methodsIa which are described briefly in the Introduction. For several molecules symmetry coordinate bases were used. Harmonic frequencies were either computed directly by the Gaussian program or were computed from the ab initio valence force field by diagonalizing Wilson's GF matrix.Ib These frequencies were then used in eq 4'* for computing the reduced isotopic partition function ratio. Results and Discussion Harmonic Force Fields. In the first part of this study, the intention is to find systematic errors in theoretical force fields for several molecules and to try to reproduce experimental results by using Moller-Plesset theory along with adequate basis sets. The basis sets used24are listed in Table 1. Throughout this discussion the abbreviated notation for basis sets in Table 1 will be used, RHF/A for RHF/3-21g*, MP3/G for MP3/6311+g(2d,2p), etc. Cartesian force constants for 14 small molecules were computed with the Gaussian program and then transformed to intemal coordinate force constants (with units mdyn/A for stretching, mdyn for stretch-bend, and mdyn*A for bending). Because of the limitations of space, the full harmonic force fields are not presented here, but these data are available as supporting information and include theoretical force fields computed using several different basis sets along with the experimental force fields for comparison. Except for ethylene and hydrogen sulfide, all of the reference force fields and harmonic frequencies were derived from experimental spectral data.32 The reference data for C 2 h comes from the work of C h o ~ g h , ~who ~ p used scaling factors to fit the MP2/6-3 l g force field for ethylene to Duncan's experimentally derived harmonic frequencies.320It should be noted that Duncan has recently readdressed the question of vibrational anharmonicity in C2& and CH2CD2.32q Chough's combined theoretical/ empirical force field is in good agreement with the theoretical results from the present study, except for the CC stretching force constant ( f i l ) , which is underestimated by MP4 theory (see discussion below). Reference harmonic force fields for HIS are taken from two different sources. Mills's force field32dis purely experimental, and Botschwina's was arrived at by making empirical corrections to a force field computed from ab initio calculations using the coupled electron pair (CEPA) approximati~n.~~j Most ab initio studies have focused on comparisons of theoretical and experimental harmonic frequencies, rather than
14692 J. Phys. Chem., Vol. 99, No. 40, 1995
Harris
TABLE 2: Multiple Bond Stretching Force Constants (See Basis Sets in Table 1) exp C2Hz 16.34 HCN 18.71 CzH4 9.48 CH20 12.94
RHFB RHF/F MP2B MP2/G MP3/G MP4/G 20.53 24.71 11.64 16.97
19.85 24.14 11.12 16.55
16.13 17.07 10.02 12.95
15.5 16.68 9.48 12.5
17.44 20.75 9.8 14.21
15.26 16.48 9.12 11.92
comparisons of force constants. Although there are several reasons to favor an analysis based on force the harmonic frequencies computed at several theoretical levels are presented in Table 5 so the present results can be compared with other ab initio studies. Specifically, we wish to compare the MP2/G frequencies for water, methane, ammonia, acetylene, ethylene, and hydrogen cyanide to results from MP2 calculations using large basis sets carried out by Handy and other^.^,^ The greatest failure of the methods used here is unreliable prediction of the multiple bond stretching force constants of hydrogen cyanide, acetylene, formaldehyde, and ethylene. These data are shown in Table 2. For all four n-bonded molecules the M E and MP4 multiple bond stretching constants are too small, and the RHF and MP3 constants are too large.I0 The failure is caused by inadequate treatment of electron correlation and not by an inadequate basis set. In a comparison of the MP4 and coupled cluster (CC) models, Stanton and coworkers found that the CC model is much better than MP4 for giving reliable stretching potentials in n-bonded It should be noted that evaluation of the second derivatives of the experimental rather than the theoretical bond lengths generally gives more accurate stretching force ons st ants.^,'^^,^^ The worst case is the CN stretching force constant for hydrogen cyanide. The RHFF,MP2/G, MP3/G, MPWG, and experimental force constants are 24.14, 16.68,20.75, 16.48, and 18.70 mdyn/A, respectively, and the corresponding CN bond lengths are 112.7, 116.7, 114.6, 116.8, and 115.3 pm. As expected, the force constants vary inversely with the bond lengths. When the MP2/G and R H F F force constants are reevaluated at the experimental bond lengths (by finite differences of the energy gradient in terms of internal coordinates, in order to avoid the use of a nonlinear transformation between internal and Cartesian force constants), the values are improved to 18.3 and 20.6 mdyn/A, respectively. In several cases the basis set was too large to compute force constants at the MP3 and MP4 levels. For ethylene the MP3/G and MP4/G force constants were estimated from the corresponding MP3F and M P 4 F force constants as follows: i,(MP4/G) =xj(MP4/F)
+ (fl,(MP2/G) - f,,(MP2/F))
(7)
2
I
1
I
6
8
10
Experimental Force Constant
L
2
4
6
8
10
Experimental Force Constant
Figure 1. Plots of experimental vs theoretical diagonal H-X stretching force constants ( m d y d a ) : (A, top) AM1 and RHFB and (B, bottom) MP2B and MP2/G. 1
2
8
.I
0.8
3 0.6
-Experiment
0 LL
2 0.4
A RHF/B
0.2
E 01
0
I
I
I
I
0.2
0.4
0.6
0.8
Experimental Force Constant I.
5
-
CI
88
0.6 -
U 0
3 0.4
.-0
E
-
0.2 -
F The MP3/G and MP4/G force constants for formaldehyde were estimated similarly. For the purpose of computing approximate MPYG and MPWG frequencies, it was necessary to adjust the MP3R and MP4/F geometries in a similar way, by linear correction using the difference between the MP2/G and MP2D geometries. Force constants at the MP4M level for HF and HCl were also estimated, from the difference between the MP2M and the MP2/G force constants. The assumption made here is that the change in force constants with the basis set is independent of the level of Moller-Plesset theory (MP2, MP3, or MP4) employed.34 Graphical analyses of the semiempirical, RHF ab initio, and MPn ab initio diagonal force constants are shown in Figures 1-3. The intention is to illustrate the reliability of the various theoretical force fields. The RHF/B diagonal force constants for H-X stretching and for angle bending are consistently overestimated (Figures 1A and 2A). The corresponding AM1
I
4
0'
0
I
I
I
I
0.2
0.4
0.6
0.8
Experimental Force Constant
Figure 2. Plots of experimental vs theoretical diagonal bending force constants (mdyn-a): (A, top) AM1 and RHFB and (B, bottom) MP2B and MP2/G.
force constants are generally in poor agreement with the experimental values, and the errors are not systematic. The MP2B and MP2/G force constants (Figures 1B and 2B) are much closer to the experimental values. The MP2/G force constants are more reliable than the MP2B force constants, particularly for bending (Figure 2B). The error in 39 of the diagonal force constants is analyzed quantitatively in Figure 3. The AM1 force constants (Figure 3A) are randomly distributed around the experimental values, with errors as large as 350%. The RHF/B force constants (Figure 3B) are all larger than the experimental values, with errors falling between 10% and 30%,
A Systematic Theoretical Study of Small Molecules
J. Phys. Chem., Vol. 99, No. 40, I995 14693
7 n
lo
-45
-25
15
-5
35
75
55
% Error 25 20 h
2
Fl
15
IRHFIB
J
!i! 10
U
5 0 -20 -10
0
10
20
30
40
50
% Error
m 0 MP4lG
-12
-8
0 % Error
-4
4
8
Figure 3. Error distribution in 39 diagonal force constants: (A, top) AM1, (B, middle) RHFB and MP2/J3, and (C, bottom) MP2/G and MP4IG.
TABLE 3: Out-of-Plane Bending Force Constants ( m d y d ) in mBonded Molecules exp AM1 RHFlB MP2lB MP2/G MPWG MP4iH CzHz 0.251 0.395 0.365 0.203 0.237 0.224 0.246 HCN C2H4 (wag) CzH4
0.260 0.406 0.274 0.334
0.382 0.342
0.257 0.269
0.251 0.255
0.241 0.247
0.139 0.096
0.165
0.148
0.144
0.138
0.403 0.392
0.510
0.432
0.408
0.396
0.261
(twist) CH2O (wag)
except for the bending force constants for HCN and C2H2, which are overestimated by 45%. The MP2B force constants (Figure 3B) are slightly overestimated on the whole, with errors ranging from -10% to +17%, except for the bending force constant of C2H2, which is 19% too small. Errors in the 39 MP2/G and MP4/G force constants are much smaller, falling in the range from -6% to +4%, with only a few exceptions (Figure 3C). Problems with computing n-bond stretching force constants were already mentioned (Table 2). In addition, the bending force constants in HCN and acetylene and the wagging force constant in ethylene are poorly predicted at the MP2/G or MP4/G levels (Table 3). These errors are caused by an inadequate basis set, which can be improved by adding f-type polarization functions to heavy This is seen by the improvement in the MP4iH force constants in Table 3. The
MP2A-I wagging force constant (f10,lo) for C2H4 is 0.267 mdyd A, which is 4.5% larger than the MP2/G constant and is only 2.6% smaller than the experimental value. Comparison of the MP4A-I and MPWG n-bending frequencies for HCN and acetylene (Table 5 ) also demonstrates the improvement from f functions. Comparison of the MFWH and MP2/G frequencies for ethylene reveals the same effect on the b2g and bl, wagging frequencies, which are respectively increased by 31 and 12 cm-'. The added f functions also cause small increases in the stretching frequencies of all three molecules. The MP2/G frequencies for ethylene (Table 5) are in reasonable agreement with Ahem's frequencies calculated at the MP2nZ+2P level,* except that the w12 (b3") bending frequencies differ by 19 cm-'. Table 4 presents the average ratios of theoretical to experimental force constants for the 14 molecules of this study, along with standard deviations (u) in these average ratios. The intention here is to find systematic errors in force constants computed at various levels. The CJ values are in the order (MP2/ F, MP2/G, MP3/G, MPLUG) < (MP2/D, MP2/B) < (RHFB, RHFF) RHF/A AM1. (See basis sets in Table 1.) This gives the approximate order of reliability of these methods for prediction of force constants after linear correction of systematic errors. It should be emphasized that the selection of molecules here is biased and that larger basis sets than the G basis are required in some cases for reliable calculation of quadratic force constants? The H-X stretching force constants are overestimated at the RHF/B level by 13%, but the error is quite systematic (u = 3~3%). The MP4/G level reproduces the experimental H-X stretching constants and the 14 angle bending force constants (excluding three It-bending constants) with no systematic error (ratios = 1.00). The MP2/G force constants are at least as accurate as the MP4/G force constants, after correction is made for systematic overestimation of stretching constants by 3% and overestimation of bending constants by 1%. It appears that these force constants can be computed reliably from MP calculations with the G basis set, except for stretching and bending constants in n-bonded
molecule^.^ Comparison of the results in Table 4 for the popular A and B basis sets shows that polarization functions on heavy atoms are important for computing force constants at the RHF level. The u values give an indication of the reliability of these methods after linear correction of systematic errors. The errors for H-X stretching are u = f0.06 (A basis) and u = f 0 . 0 2 (B basis), and the errors for angle bending are u = f 0 . 1 8 (A basis) and u = fO.10 (B basis). Change at the MP2 level from the B or D basis sets to the E basis set significantly improves the computed force constants. This is the result of having polarization functions on all atoms (E basis), and not just on heavy atoms (B or D basis). The addition of p functions on hydrogens to form the E basis from the D basis causes significant changes ('2%) at the MP2 level in diagonal force constants for HzO, HzS, S i b , PH3, CH20, HCN, HF, HCl, NH3, CZHZ,CH4, and C 2 b . It is clear that polarization functions are needed on hydrogen atoms as well as heavy atoms for predicting force fields accurately with MP theory. Some instances in which the change from D to E basis causes especially large changes in force constants will be noted. In H20 the bending constant decreases 7.7%; in H2S the stretching force constant increases 8.2%, and the bending force constant decreases 8.3%; in PH3 the stretching force constants (fll and f33) increase 5.5%, and the bending force constants (f22 and&) decrease by 5.6% and 3.8%; in HCN the bending force constant increases 6.4%; in HF and HCl the stretching force constants
14694 J. Phys. Chem., Vol. 99, No. 40,1995
Hamis
TABLE 4: Average Ratios of TheoreticaVExperimental Force Constants (Errors are fl Standard Deviation. See Basis Sets in Table 1) I I1 I11 IV
I I1 I11 IV
All diagonal force constants (n = 39) Diagonal H-X stretching force constants (n = 16) All diagonal bending force constants (n = 17) Diagonal bending force constants excluding acetylene, HCN, and ethylene (flo,lo) bending constants ( n = 14)
AM1
RHFIA
RHFB
1.12f0.24 1.07f0.11 1.061t0.25 0.97 f 0.17
1.19f0.15 1.11fO.06 1.28f0.18 1.21 f 0.06
1.19f0.09 1.13f0.02 1.22f0.10 1.19 f 0.06
RHFF
MP2B
1.153~0.09 1.04fO.06 1.11f0.03 1.04f0.03 1.18f0.11 1.05f0.08 1.14 f 0.06 1.08 f 0.05
increase by 7.9% and 9.6%; in NH3 the bending force constant decreases 8.6%; in C2H2 the bending force constant increases 22.8%; in Si& the stretching force constants increase 3.5% and the f 4 bending force constant increases 2.9%. It is notable that in geometrically similar cases adding polarization functions for hydrogens to the basis causes a qualitatively similar effect. To mention some examples, the bending force constants for acetylene and HCN bending both increase, the bending force constants in PH3 and NH3 both decrease, the bending force constants in H20 and H2S both decrease, and the stretching force constants in HF and HC1 both increase after changing from the D to the E basis. The inclusion of p functions on hydrogens also causes a marked increase in stretching force constants for all the second-row hydrides. The conclusion here is that the basis is improved significantly by adding polarization functions to all atoms rather than just heavy atoms.35 Addition of diffuse f u n c t i o n ~to~ ~ the, ~6-31g* ~ or 6-311g** basis sets produces a significant (>2%) change in diagonal force constants computed at the MP2 level for several molecules. The C-F stretching constant (f33) for CH3F decreases 9% (B to C) or 5% (E to F); the stretching constant in HF decreases 5% (B to C) and 2.5% (E to F); bending constants for both H20 and NH3 decrease 5% (B to C or E to F). In addition, added diffuse functions cause significant (> 2%) changes in the the out-ofplane bending constants of ethylene ( ~ I O , I O and fi2,12) and formaldehyde (fa)and in the bending constant of hydrogen cyanide (fa,). In each of these cases the force constant is decreased by adding diffuse functions. Improvement of the polarization from one set of d,p functions (F basis) to two sets of uncontracted d,p functions (G basis) is quite important for computing H-X stretching force constants using MP theory when X is a second-row atom (Si, P, S , and Cl). Thus, the H-X stretching force constants for S i b , PH3, H2S, and HC1 decrease by 2.2%, 4.4%, 4.1%, and 5.1%, respectively, on changing from the M P 2 F to the MP2/G level. In addition, the C-C1 stretching force constant in CH3C1 decreases by 7.4%. In several other cases the change from F to G basis causes significant effects. The C-F stretching constant in CH3F decreases 2.9%, the bending force constants in C2H2 and HCN decrease by 3.7% and 4.9%, the f 6 6 bending constant in CH3Cl decreases 4.7%, the wagging force constants in C2H4 and CH20 are respectively increased 3.2% and decreased 3.3%, thef22 bending force constant in PH3 decreases 7.2%, and the bending force constant in H2O increases 3.9%. However, many diagonal and off-diagonal force constants are changed less than 1%, and so the E and F basis sets seem to be adequate for computing MP2 force constahts in many molecules. The G basis is not adequate for computations on either water or ammonia, which is seen by comparing the MP2/G frequencies in Table 5 with MP2 frequencies computed by Handy and coworkers7 with a (8s6p3dlf/6s3p) basis set. The MP2/G frequencies for the e and a1 bending modes of ammonia are
MP2D
MP2E
MP2F
MP2IG
MP3lG
MP4lG
1.02f0.06 1.02f0.02 1.04f0.08 1.06 f 0.05
1.03f0.04 1.05f0.02 1.03f0.04 1.03 f 0.03
1.02f0.04 1.05f0.02 1.01fO.04 1.02 It 0.03
1.01f0.03 1.03f0.01 1.OOfO.03 1.01 0.01
1.03f0.03 1.03f0.02 1.01f0.03 1.01 f 0.03
0.98f0.04 l.OOfO.O1 0.99f0.04 1.00 f 0.02
*
1057 and 1691 cm-I, compared to Handy’s values of 1035 and 1675 cm-I. When the basis for NH3 is improved from G to J (see Table 1) the e and al frequencies are improved to 1032 and 1672 cm-I. The MP4IJ umbrella (al) frequency is 1060 cm-’, which agrees with Loutellier’s recent experimental value32cand with the results of recent large basis set calculations using coupled cluster theory.12d Change from the G to the J basis causes a 2.5% decrease in the NH3 bending force constant at both the MP2 and MP4 levels. The resulting MP4/J bending force constant of 0.622 mdyn-A is in near perfect agreement with Loutellier’s value of 0.620 but is 1.7% smaller than Duncan’s value.32b Improvement of the G basis was also explored for the hydrogen halides, in MP2 calculations with the I, K, L, and M basis sets (see basis sets in Table 1 and frequencies in Table 5). For hydrogen fluoride the force constant changes less than 0.5% going from the G to the M basis. However, the force constant fluctuates by about 1% as the basis is improved, showing that the G basis is improved by adding extra d,p polarization functions and f second polarization functions. In the case of hydrogen chloride the force constant is increased by 2.5% going from the G to the M basis. The estimated MP4M force constants for HF and HC1 differ from the respective experimental values32kby only 0.1% and 0.4%. Computation for methane does not have the same problems with basis sets as for hydrides of the more electronegative firstrow elements (ammonia, water, and hydrogen fluoride). The MP2/G frequencies (Table 5) for C& are in very good agreement (8 cm-I) with MP2 frequencies calculated using a (8s6p3dlf/6s3p) basis.7b Basis sets lacking any polarization functions were not explored here. Karplus previously compared the performance of the 6-31g, 6-31g*, and 6-31g** basis sets, along with the corresponding 6-311g, 6-311g*, and 6-311g** basis sets, in MP2 computations of geometries and harmonic vibrational frequencies for several small molecules.35 It seems obvious from the present study and Karplus’s study that polarization functions on all atoms (Le., a D Z f P basis set) is a minimum requirement for accurate computation of harmonic frequencies and molecular geometries. The results here and results from several previous studies show that the D Z f P basis is only a starting point, and that force fields can be improved in most cases by expanding the s,p basis set and adding extra polarization function^.^ It seems worthwhile to mention a few cases where large differences in off-diagonal force constants occur between theory and experiment. (The following data are taken from the force fields provided in the supporting information.) Thefd andfrainteraction constants for NH3 agree better with Loutellier’s force field than with Duncan’s earlier s t ~ d y . ~ * ~ ~ ~ As already mentioned, the diagonal bending force constant (faa) also agrees better with Loutellier’s value. The theoretical constant of 0.04 is in poor agreement with both Loutellier’s value (0.004) and Duncan’s value (0.145). frat
A Systematic Theoretical Study of Small Molecules
J. Phys. Chem., Vol. 99, No. 40, 1995 14695
The theoretical interaction constants 5 2 , fi 3, andf46 for CH3deuterated compounds were used in eq 4. To find the reduced isotopic partition function ratio per deuterium, the resulting value C1 are respectively 0.13 mdyn, 0.08 mdyn/A, and 0.03 mdyn, of (sdsl)f was then raised to the power lln, where n is the all in poor agreement with Duncan's values of 0.03, 0.17, and 0.06.32Y It should be noted that Duncan has recently readdressed number of deuteriums. For example, in the case of silane, (s$ the anharmonicity of CH3C1,32xbut only anharmonic constants si)fdl FZ [(s2/sl)fd4]"4. This method assumes that the partition function ratio is constant (aside from symmetry effects) for each and not harmonic force constants were derived from that study. successive substitution of deuterium, (e.g., Q s ~ D ~ / Q s ~ D ~%H In the case of CHsF, all the theoretical interaction constants are substantially different from Duncan's experimental values.32w Q S I D ~ H / Q S I D ~ H% ~ Q S I D ~ H ~ / Q S I D H% ~ QsIDH,/QsIH~). The poorest agreement is for the 5 2 and f46 stretch-bend The last row of Table 6 gives the ratios exp[2]/exp[l]. These ratios range from 0.86 to 0.94. The variation of 8% in the expconstants, for which the theoretical values are 0.10 and 0.05 mdyn, and the experimental values are -0.22 and 0.18 mdyn. [2]/exp[11 ratios shows that equilibrium isotope effects (EIEs) computed using eq 5 can be expected to vary by as much as The experimental force field for ~ i l a n e ~has ~ gvalues for the 8% for the selection of molecules here, depending on the choice diagonal fi 1 stretching force constant and the f34 stretch-bend of anharmonic or harmonic frequencies.Igb interaction constant of 3.35 mdyn/A and -0.41 mdyn, which The theoretical values of (sz/sl)f in Table 6 are all overestiare in poor agreement with the MPWG values of 3.06 and -0.10 mated, especially values computed from RHF frequencies. Force mdyn. It seems reasonable to assume that the error is in Levin's constants are overestimated by 10-30% at the RHF level, and force field and not in the MPWG force field, since the MP4/G HX stretching force constants are overestimated by 3% at the theory performs well for the other second-row hydrides. The MP2 level (see'Table 4), resulting in overestimation of the (s2/ MP4/G force field yields harmonic frequencies of 2272 and 2212 s1)f values. As mentioned already, agreement between theoreticm-' for the Ai and T2 stretching modes of S i b and 1607 and cal harmonic frequencies and experimental frequencies is 1641 cm-' for the respective modes in SiD4. These values agree improved by uniform scaling of theoretical frequencies l4 or by with Kattenberg's fundamental frequencies32hof 2186 and 2189 nonuniform scaling of the theoretical force field.'5a Uniform cm-' for S i b and 1563 and 1598 cm-' for SiD4, after allowance scaling was used by Williams for computing isotope effects for anharmonicity in the experimental data. more These differences between the present theoretical results and As noted already, changes in the zero-point energy (ZP) term earlier experimentally derived force fields cannot be dismissed account for most of the variations in (sz/s~)f for different as errors in the theory caused by incomplete basis sets or molecules. Uniform scaling of frequencies has no effect on incomplete treatment of electron correlation. Pulay has repeatthe masdrotational term and has very little effect on the edly mentioned that theoretical values are often more reliable term, which is very close to unity. The effect on the than experimental values for off-diagonal force c o n s t a n t ~ . " , ~ ~ ~ excitation .~~~ ZP term is to raise ZP to the power a, where a is the uniform Revisions to the force fields for CH3C1, CH3F, and S i b seem scaling constant. necessary. In the case of NH3 the present results agree better with the more recent experimental force field of Loutellier than (8) ZPscaled = (ZPunscaled)u = n[exp(a(ull - u2i))1 with Mills's force field. I Reduced Isotopic Partition Functions and Fractionation Factors. In the second part of this study attempts are made to Because the scaling constant enters eq 8 as an exponent, small accurately calculate the reduced isotopic partition function ratios deviations of a from unity produce a large effect on the ZP (s2/sl)fand and fractionation factors (FFs) for the 14 molecules value. For a typical value of ZP = 20, a 1.0% overestimation in part 1. The (s2Isl)f values presented in Table 6 were of harmonic frequencies leads to a 3.0% overestimation of ZP, computed from eq 4, using frequencies for the nondeuterated and a 5% overestimation of harmonic frequencies leads to a and monodeuterated isotopomers. Two reference values were 16% overestimation of ZP. Small systematic errors in theoreticomputed for each molecule, one based on experimental cal frequencies should be corrected by scaling to avoid large harmonic frequencies and one based on experimental anharerrors in the (s2/s1)f value. However, errors in theoretical EIE monic frequencies (fundamental frequencies). These values, values resulting from systematic errors in frequencies will be labeled exp[l] and exp[2], are the same quantities which smaller, because EIE values computed using eq 5 contain ratios Goodson and Wolfsbergigbcalled [HARM] and [ANHARM]. of ZP terms, which are close to unity. For example, if ZPB/ The exp[ 11 values in Table 6 were computed from harmonic ZPA= 1.50, then a 5.0% overestimation of frequencies for both force fields taken from the l i t e r a t ~ r e .The ~ ~ exp[l] value for A and B only leads to an error of 2.0% in the theoretical EIE. silane is based on the MPWG force field, because the experiUniform scaling factors were calculated according to eq 9, mental force field3*gcontains serious errors which have already so as to best reproduce experimental zero-point energy differbeen mentioned. The exp[ 11 values for HF, HCl, CH3F, HCN, ences for the 14 molecules of part 1. CH20, C2H2, C b , H20, and NH3 are in excellent agreement with Goodson's [HARM] valuesigbcomputed at 300 K, after a (9) a = c((w -l W2ij)expJ(Wly ij - W2ij)theory) temperature correction is made. The corrected exp[ 11 values IJ for C2H4 and CH3Cl are about 4% smaller than Goodson's values. In eq 9 w l i and wzij are the harmonic frequencies for the light and heavy isotopomers, and the summation runs over the Most of the exp[2] values in Table 6 are Goodson's 3N - 6 normal modes for each of the 14 molecules. The [ANHARM] values, again corrected to 25.0 "C. These values resulting scaling constant is the average ratio of experimental/ are derived from harmonic force fields adjusted to reproduce theoretical zero-point energy differences. The uniform scaling anharmonic fundamental frequencies for several isotopic derivafactors are given in Table 7 for the AM1 RHF/A, RHFB,RHF/ tives of each molecule. The exp[2] values for HzS, PH3, S i b , F, MP2/B, MP2/F, and MP2/G methods. The scaling factors and cyclopropane were computed from fundamental frequencies for the RHF methods are ca. 0.93, and those for the MP2 of the deuterated and nondeuterated isotopomers, without the methods and the semiempirical AM1 method are ca. 0.98. The usual force field refinement.Igb.*OFor SiH4, PH3, and cycloscaling factors for RHF/B and MP2B are about 2% larger than propane, frequencies of the nondeuterated and completely
14696 J. Phys. Chem., Vol. 99, No. 40, 1995
Harris
TABLE 5: Harmonic Frequencies, cm-' (See Basis Sets in Table 1) C2H2 % 816 583 625 624
RHFF MP2/G MP2" MP3lG
nu
a*+
873 737 763 759
2205 1960 1989 2073
%+ 3559 3439 3453 3445
%+ 3673 3525 3542 3544
MP4/G MP4M expb
4+
ng
XU
552 593 624
722 748 747
%+ 3402 3415 3415
1944 1960 2010
a,+ 3489 3504 3496
NH3
e RHFF MP2/G MP2/J MP2' MP3/G
al 1790 1691 1672 1675 1708
1090 1057 1032 1035 1078
e 3697 3533 353 1 3520 3559
e 1082 1077 1060 1022 1057
at
MP4lG MP4D MPWJ expd exp'
3827 3672 3676 3670 3678
e
al 1691 1681 1674 1691 1666
at 3615 3602 3619 3577 3590
3489 3472 3487 3506 3457
CH3Cl RHFF MP2/G MP3/G
al
e
a1
e
al
e
772 753 749
1120 1049 1050
1511 1413 1421
1599 1512 1515
3225 3129 3136
3327 3238 3236
MP4/G ex$ expg
al
e
al
e
a1
e
734 751 740
1036 1037 1038
1402 1396 1383
1500 1496 1482
3089 3088 3074
3190 3183 3166
1502 1490
1516 1498
3058 3031
3150 3131
CH3F RHFF MP2/G MP3/G
1157 1057 1095
1296 1213 1224
1611 1518 1531
1615 1533 1535
3194 3104 3097
3277 3205 3186
MPWG exp
1029 1059
1199 1206
HCN RHFF MP2/G MP2' MP3/G
n
U+
U+
878 715 642 752
2406 2014 2030 2237
362 1 3463 346 1 3490
MPWG MP4A-I exp'
n
U+
a+
698 725 727
2000 2019 2129
3419 3435 3442
czH4 10, b2u RHFF MP2/G MP2/H MP2J MP3IGk MP4IGk expf
8, b2g 1098 934 965 942 945 912 959
890 835 83 1 826 837 826 843
7, blu 1081 974 986 967 978 956 968
4, a, 1135 1070 1074 1058 1064 1045 1047
6, blg 1337 1246 1249 1235 1254 1238 1245
3, a, 1470 1382 1385 1378 1393 1365 1371
12, b3" 1583 1492 1485 1473 1504 1486 1473
2, a,
11,hu
1, a,
5,blg
9,b2u
1813 1677 1681 1677 1706 1659 1654
3270 3172 3186 3172 3178 3137 3138
3292 3191 3204 3192 3197 3157 3139
3347 3263 3275 3263 3259 3221 3211
3375 3289 3301 3292 3286 3248 3234
CH2O RHFF MP2/G MP3/Gk
1333 1200 1212
1362 1283 1299
1649 1550 1575
1996 1747 1850
3097 2971 2985
3168 3047 3060
MP4/Gk exp"
1176 1190
1263 1287
1530 1562
1717 1763
2921 2944
2985 3008
Hydrogen Halides HF
w
w
RHFF MP2/G MP2D
4491 4165 4166
HCl
MP2IM MP3lG MPUG
4156 4259 4146
F w
0
RHFlF MP2/G MP2D
3144 3003 3027
MP2IM MP3/G MP4/G
3043 2991 2959
w
MP2/K MP2L
4138 4138
HC1 w
MP4/Mk 4136 exp" 4138
w
0
MP2/K MP2L
3050 3046
MP41Mk 2999 exp" 2990
CH4 RHFF MP2/G MP2"
t2
e
al
t2
1452 1370 1362
1666 1596 1594
3 149 3083 308 1
3252 321 1 3218
al
e
e
al
1111 1017 1017
1244 1166 1159
2550 2465 2436
2555 246 1 2435
t2
e
a]
t2
MP3/G MP4/G exP
1376 1365 1367
1593 1581 1582
3077 3042 3025
3192 3158 3156
al
e
e
al
MP4/G expq
1010 1012
1149 1141
2408 241 1
2408 2405
PH3 RHFF MP2lG MP3/G
A Systematic Theoretical Study of Small Molecules
J. Phys. Chem., Vol. 99, No. 40, 1995 14697
TABLE 5 (Continued) Si& RHFlF MP2lG MP3lG
t2
e
al
tz
1018 966 95 8
1051 1018 1010
2351 2304 2287
2337 2304 2287
MP4lG exp' exps,'
12
e
at
t2
952 945 913
1005 975 972
2272 2377 2185
2272 2319 2189
H2S a1 2869 2757 2735
al
RHFlF MP2lG MP3lG
1319 1229 1234
bz 2878 2775 2750
MP4lG exp"
a1 1225 1212
a1 2706 273 1
b2 272 1 2745
MP3lG MP4lG exp"
1693 1673 1648
3940 3836 3832
4035 3946 3942
HzO RHFF MP2IG MP2"
1726 1660 1629
4142 3861 3839
4244 3981 3966
a C2Hz MP2 calculations with TZ2PSf basis set from ref 7a. C2H2 ref 321. NH3 MP2 calculations with (8s6p3dlfl6s3p) basis set from ref 7c. NH3 ref 32b. e NH3 ref 32c. fCH3Cl ref 32x. g CH3F ref 32w. * HCN MP2 calculations with (8s6p3dI6s2p) basis set from ref 7c. ' HCN refs 32r-t. C2H4 MP2 calculations with TZ 2P basis from ref 8 (see also ref 7a). Estimated frequencies (see text). C2& ref 320. CH2O ref 3211. HCl and HF ref 32k. C& MP2 calculations with (8s6p3dlfI6s3p) basis set from ref 7b. p CH4 ref 32a. * PH3 ref 32i. Si& ref 32g. SiH4 fundamental frequencies from ref 32h. Observed fundamental frequencies. HzS ref 32j. H20 MP2 calculations with (8s6p4dlfI6s3p) basis set from ref 7b. H2O ref 32e.
+
J
TABLE 6: Reduced Isotopic Partition Function Ratios (.vJ.v~)f at 25.0 "C (See Basis Sets in Table 1) AMla RHF/A RHFB RHFlF MP2lB MP2lF MP2IG exp[ 11 exp[2Ic ratiod
C2H2
NH3
CH3Cl
CH3F
HCN
C2&
CH20
HCl
HF
C&
PH3
Si&
10.28 12.43 11.85 11.51 9.02 9.38 9.19 9.18 8.45 0.92
13.95 16.25 17.56 16.85 15.02 14.41 14.55 13.74 11.82 0.86
14.45 16.91 16.78 15.93 14.59 13.78 13.86 12.81 11.60 0.91
12.77 17.49 17.52 16.65 15.01 14.47 14.61 13.53 12.10 0.89
10.45 13.10 11.93 11.34 9.47 9.35 9.13 9.05 8.25 0.91
11.46 14.86 14.59 13.92 12.54 11.92 12.03 11.59 10.39 0.90
10.82 13.53 13.35 12.62 11.28 10.83 10.81 10.54 9.12 0.87
5.68 6.16 6.31 6.12 5.74 5.90 5.56 5.52 5.17 0.94
13.98 10.72 13.07 14.29 10.58 11.76 11.51 11.34 10.15 0.90
11.34 14.34 14.11 13.32 12.82 12.07 12.22 11.73 10.29 0.88
7.09 6.98 6.92 6.58 6.25 6.16 5.91 5.61 5.09 0.91
5.44 5.87 5.81 5.63 5.38 5.43 5.32 5.16' 4.64 0.90
H2S
HzO
c-C3H6
6.96 7.03 6.69 6.32 6.31 6.05 5.91 5.15 0.87
11.56 14.05 16.65 16.86 13.48 13.91 13.86 13.48 11.71 0.87
12.97 16.73 16.25 15.45 14.20 13.60 11.06
Li MNDO calculations for molecules containing second-row atoms. exp[ 11: experimental value based on harmonic frequencies for deuterated and nondeuterated isotopomers. C2H2, ref 321; NH3, ref 32c; CH3C1, ref 32y; CH3F. ref 32w; HCN, ref 32t; C2&, ref 320; CH20, ref 32n; HC1 and HF, ref 32k; C&, ref 32a; PH3, ref 32i; Si&, see note e above; HzS, ref 32j; H20, ref 32e. exp[2]: experimental value based on fundamental (anharmonic) frequencies, taken from ref 19b except for SiH4, ref 32h; PH3, ref 32i; H A ref 32j; c - C ~ H ref ~ . 322. Ratio of exp[2]/exp[11. See ref 19b. e Value based on MP416-31 l+g(2d,2p) frequencies for SiHdSiH3D.
TABLE 7: Uniform Scaling Factors (See Basis Sets in Table 1) AMla RHFIA RHFIB RHFF MP2IB MP2F MP2lG scale factor 0.982 0.936 0.930 %errort' 11.0 10.1 3.8
0.940 3.7
0.982 5.5
0.983 2.5
TABLE 8: Average Signed Error (e) and Root Mean Square Error (0)in Reduced Isotopic Partition Function Ratios for 14 Molecules (See Basis Sets in Table 1) (A) No scaling of frequencies (B) Uniform scaling of frequencies
0.988 2.2
a MNDO calculation for molecules containing second-row atoms. Root mean square error in ratio ZPobJZPcaIc after scaling frequencies. % error = 1OO[E,(1- ZP,~,.,/ZPca~c,,)2)113]1'2 (ZP = n,[exp(-uj2)/exp(-uti)]; ut = UJkT).
those found by Pople and c o - w o r k e r ~because ,~~ Pople's scaling factors are based on true zero-point energies instead of harmonic zero-point energies. (The energy levels for an anharmonic diatomic oscillatorId are given approximately by E, = ( v f '/2)hwe - ( V f ' / 2 ) 2 ~ e ~SO e , the correction for a diatomic is approximately ZPhm - ZPanh- = 1/4xew,.) The accuracy of the theoretical ZP values after uniform scaling is in the order AM1 RHF/A MP2B (RI-FB, RHF/F) (MP2/F, MP2/G), as indicated by the root mean square errors in Table 7. The large error (11.4%) for AM1 together with the near unit scaling factor (0.982) shows that the AM1 zero-point energies cannot be improved much by scaling. The error in theoretical (~2/sl)fvaluesis analyzed in Table 8. The average signed errors (6) after scaling are close to zero, showing that the systematic overestimation of (s2/sl)f
B
A AM1 RHFIA RHFlB RHFlF MP2lB MP2ff MP2lG exp[ 11
€
C7
0.46 2.17 2.45 2.08 0.59 0.46 0.39
1.24 2.63 2.69 2.30 0.90 0.52 0.52 (0.00)
(0.00)
% erroP
-0.11 -0.05 -0.04
-0.00 0.04 -0.07 0.03
1.15 0.99 0.31 0.32 0.56 0.22 0.24
10.73 9.01 3.19 3.27 5.02 2.50 2.13
Fractional error in (s2sl)f after scaling frequencies. % error = 100[x-l((fobs.r- fcalc.i)/fob~.i)~I a
values is corrected. The root mean square errors in Table 8 show the order of reliability of the various methods in computing (s2Isl)fafter scaling is again AM1 < RHF/A M P D -= (RHF/ B, RHF/F) (MP2/F, MP2/G). The deuterium fractionation factors (FFs) in Table 9 were computed from the ratio [ ( ~ Z / S I ) ~ / [ ( S ~ / where S I ) ~ *(s2/sl)f* ], = 9.18 is the reference value computed from the harmonic frequencies321 for acetylene and acetylene-dl . The (S~/SI )f values
14698 J. Phys. Chem., Vol. 99, No. 40, 1995
Harris
TABLE 9: Hrdrogen/Deuterium Fractionation Factors vs Acetylene at 25.0 "Ca(See Basis Sets in Table 1) CH3Cl CH3F HCN HCl HF CHI S i b H2S C2H2 C2H4 CH20 NH3 PH3 AMlb RHF/A RHF/B RHFF MP2/B MP2E MP2/G exp[1I' exp[2Id
exp[3Ie
1.06 1.10 1.03 1.04 0.93 0.97 0.97 (1 .OO) (1 .OO) (1 .OO)
1.43 1.43 1.50 1.50 1.54 1.48 1.53 1.50 1.40 1.48
1.48 1.46 1.42 1.40 1.50 1.42 1.45 1.39 1.37 1.40
1.31 1.51 1.48 1.46 1.54 1.48 1.53 1.47 1.43 1.46
1.08 1.15 1.04 1.02 0.98 0.97 0.96 0.99 0.98
1.18 1.30 1.25 1.24 1.29 1.23 1.26 1.26 1.23 1.26
1.11 1.19 1.16 1.13 1.16 1.12 1.14 1.15 1.08 1.11
0.60 0.59 0.59 0.59 0.60 0.62 0.59
0.60 0.61
1.45 0.99 1.16 1.30 1.10 1.22 1.21 1.24 1.20
1.17 1.26 1.22 1.19 1.32 1.24 1.28 1.28 1.22 1.25
0.74 0.65 0.63 0.62 0.65
0.64 0.63 0.61 0.60 0.65
0.57 0.55 0.53 0.53 0.56 0.57 0.56 0.56' 0.55
1.19 0.65 0.65 0.63 0.66 0.66 0.64 0.64 0.61 0.65
H2O 1.33 1.26 1.44 1.51 1.39 1.44 1.46 1.47 1.39 1.54
c-C~H~ 1.45 1.38 1.36 1.45 1.40 1.31 1.40
a FF = [(s2/slM/9.18. Values of (sdsl)fwere computed from uniformly scaled frequencies. MNDO calculations were performed for molecules containing second-row atoms. exp[l]: experimental FFs based on harmonic frequencies. exp[2]: experimental FFs based on anharmonic frequencies computed from [(s2/sJfl/8.45. e exp[3]: FFs from Shiner, V. J.; Neumann, T. E. Z. Narurforsch. 1989, 44a, 337-354. fValue based on MP4/63 1 l j g ( 2 d . 2 ~ )frequencies for SiaSiDd.
TABLE 10: Reduced Isotopic Partition Function Ratios for Various Molecules Calculated from Uniformly Scaled MP2/6-311+G(2d,2p) Frequencies R ZP EX (s2/sI)f errop 0.445 C2H2 HCN 0.457 HF 0.725 0.717 HCl 0.562 NH3 0.544 PH3 0.524 CH4 0.505 Si& 0.458 CH3F 0.449 CH3Cl 0.454 C2H4 CHzO 0.502 0.637 H20 0.625 HIS C - C ~ H ~ 0.403 ~ a
Error = (
18.350 18.037 15.366 7.580 24.840 10.450 22.379 10.045 30.35 1 29.171 24.939 20.634 20.968 9.402 30.650
1.088 1.070 1.ooo 1.000 1.004 1.010 1.004 1.021 1.009 1.016 1.024 1.007 1.001 1.003 1.039
8.88 8.82 11.14 5.43 14.02 5.74 11.77 5.18 14.03 13.31 11.59 10.43 13.37 5.89 12.84
s ~ / s ~-& (s~/sl)ferpenmenl. ~ ~ ~ ~ ~ 6-31 l+g**
-0.30 -0.22 -0.21 -0.08 0.28 0.13 0.05 0.02 0.49 0.5 1 0.0 1 -0.1 1 -0.11 -0.01
basis.
in Table 6 were first corrected by uniform scaling before computing the FFs. The last row of Table 9 gives the FFs computed by Shiner and co-workers20 based on anharmonic frequencies at 25.0 "C. As mentioned already, Shiner's FFs are based on a harmonic force field, which is fitted to reproduce fundamental frequencies for several isotopic derivatives of a molecule. The exp[ 11 values (based on harmonic frequencies) in Table 9 agree better with Shiner's values than the exp[2] values (based on anharmonic frequencies), the largest differences being for water (7%) and formaldehyde (4%). Anharmonic effects explain the difference for formaldehyde, but the poor agreement between the exp[ 11, exp[2], and Shiner FFs for water is confusing. Finally, the ab initio and Shiner results are in agreement that the FF = 1.40 for cyclopropane is a little smaller than that for a typical secondary sp3 hybridized carbon.20a The "best" theoretical values for (sz/s~)fpresented in Table 10 were computed from MP2/G frequencies which were scaled uniformly by 0.987. Cyclopropane was too large to calculate with the G basis, so the (s2Isl)f value was computed at the MP2E level with scaling by 0.983. The (s2/~1Ifvalue can be factored into R, ZP, and EX terms, as was discussed in the introduction. The R values in Table 10 range from 0.40 to 0.72. The R values are similar for geometrically similar molecules, such as HC1 and HF, NH3 and PH3, etc. Stretching and bending potentials are decreased in the second-row hydrides compared to their first-row counterparts, which results in relatively small ZP values for the secondrow hydrides. For the first-row hydrides, the ZP values are in the order NH3 > C b > H20 > HF. The same order is observed for the second-row hydrides. The EX terms are all
smaller than 1.1 at 25 "C. The largest EX values are for acetylene and HCN, which possess low-frequency bending modes. The errors for CH3F and CH3C1 in Table 10 are the only ones which are larger than 2 a (a = 0.24; see Table 8). It seems likely that the experimental (s2Is1)fvalues for these two molecules are too small. It was already mentioned that the experimental force fields for both molecules are in poor agreement with the present theoretical force fields, largely because of errors the experimental force f i e l d ~ . ~ ~ The " J (s2/ s1)fvalues calculated at the MP4/G level are 13.30 for CH3C1 and 13.92 for CH3F, in excellent agreement with the scaled MP2/G values. Summary and Conclusions The force fields for the molecules considered here are computed reliably at the MPWG level, except that the G basis is inadequate for computing several bending force constants (see Table 3) and MP theory is inadequate for stretching force constants of multiple bonds (see Table 2). Harmonic force fields computed at the MP2 and RHF levels are improved by correcting a few systematic errors. Specifically, the H-X stretching force constants are overestimated at the MP2 level by about 3% and at the RHF level by about 12%. Bending force constants are overestimated by 1-3% at the MP2 level and by 14-19% at the RHF level (see Table 4). The E and F basis sets, with polarization functions on all atoms, are sufficient for MP calculations in most cases. These basis sets are small enough for MP2 calculations on unsymmetrical molecules containing three or four heavy atoms. The B basis, which has polarization functions only on heavy atoms, does not perform as well in MP2 calculations. The improved polarization functions in the G basis are important for accurately computing stretching force constants for the second-row hydrides. The force constants are improved for ammonia, water, the hydrogen halides, HCN, acetylene, and ethylene from improvements of the G basis set. The MP2/G force field for CH3F would also be improved with a larger basis set. The reduced isotopic partition function ratios are sensitive to the small overestimation of stretching force constants by MP2 theory. Uniform scaling of MP2 frequencies improves the (s2/ sl)f values considerably. Scaling is even more important for computing accurate (s2Isl)f values with RHF theory. Ideally isotope effects should be computed from scaled MP2 frequencies with the E or F basis sets, or a basis set of similar quality. Computation of the isotope effects here with scaled RHFB or W/Ffrequencies is more reliable than computation with scaled MP2B frequencies. Computations at the W / A or AM1 levels are best regarded as semiquantitative. Nonuniform scaling methodsI5 were not used here, but these methods could be used profitably in computing isotope effects with RHF theory, since
A Systematic Theoretical Study of Small Molecules the systematic errors in different types of valence force constants can vary widely. For example, the bending force constants in HCN and acetylene are overestimated by 45% while other force constants are usually overestimated by 10-20% using W/B theory. In cases where MP2 calculations are not practical, the best option is refinement of the ab initio RHF force field by nonuniform scaling to fit experimental vibrational frequencies. 15,19b.20 Acknowledgment. I thank the NSF for support of this work through a grant to J. J. Gajewski and thank Prof. V. J. Shiner for encouragement. Supporting Information Available: Tables of theoretical harmonic force fields for the molecules of this study, along with experimental force fields for comparison (10 pages). Ordering information is given on any current masthead page. References and Notes (1) (a) Califano, S. Vibrational States; Wiley: London, 1976. (b) Wilson, E. B., Jr.; Decius, J. C.; Cross, P. C. Molecular Vibrations; McGraw-Hill: New York, 1955. (c) Gans, P. Vibrating Molecules; Chapman and Hall: London, 1971. (d) Graybeal, J. D. Molecular Spectroscopy; McGraw-Hill: New York, 1988. (2) Hoy, A. R.; Mills, I. M.; Strey, G. Mol. Phys. 1972, 24, 12651290. (3) Miyazawa, T. J. Chem. Phys. 1958, 29, 246. (4) Fogarasi, G.; Pulay, P. Annu. Rev. Phys. Chem. 1984, 35, 191213. (5) (a) Pople, J. A.; et. al. Int. J. Quantum Chem., Symp. 1981, 15, 269-278. (b) Pople, J. A.; Krishnan, R.; Schlegel, H. B.; Binkley, J. S. Int. J. Quant. Chem., Symp. 1979, 13, 225-241. (6) Hout, R. F., Jr.; Levi, B. A,; Hehre, W. J. J. Comput. Chem. 1982, 3, 234-250. (7) (a) Handy, N. C.; Simandiras, E. D.; Rice, J. E.; Lee, T. J.; Amos, R. D. J. Chem. Phys. 1988.88, 3187-3195. (b) Handy, N. C.; Gaw, J. F.; Simandiras, E. D. J. Chem. Soc., Faraday Trans. 2 1987, 83, 1577-1593. (c) Simandiras, E. D.; Handy, N. C.; Amos, R. D. Chem. Phys. Lett. 1987, 133, 324. (d) Pulay, P.; Lee, J.-G.; Boggs, J. E. J. Chem. Phys. 1983, 79, 3382-3391. (8) Ahern, A. M.; Garrell, R. L.; Jordan, K. D. J. Phys. Chem. 1988, 92, 6228-6232. (9) Sellers, H.; Almlof, J. J. Phys. Chem. 1989, 93, 5136-5139. (10) Stanton, J. F.; Watts, J. D.; Bartlett, R. J. J. Chem. Phys. 1991,94, 404-41 3. (11) (a) Besler, B. H.; Scuseria, G. E.; Scheiner, A. C.; Schaefer, H. F., 111. J. Chem. Phys. 1988, 89,360-366. (b) Thomas, J. R.; DeLeeuw, B. J.; Vacek, G.; Schaefer, H. F., 111. J. Chem. Phys. 1993, 98, 1336. (c) Allen, W. D.; Csszh, A. G. J. Chem. Phys. 1993, 98, 2983-3015. (12) (a) Kobayashi, R.; Amos, R. D.; Handy, N. C. J. Chem. Phys. 1994, 100, 1375-1379. (b) Lee, T. J.; Rendell, A. P. Chem. Phys. Lett. 1991, 177, 491-497. (c) Martin, J. M. L.; Francois, J. P.; Gijbels, R. J. Chem. Phys. 1992, 96, 7633-7645. (d) Martin, J. M. L.; Lee, T. J.; Taylor, P. R. J. Chem. Phys. 1992, 97, 8361-8371. (e) Martin, J. M. L.; Taylor, P. R.; Lee, T. J. Chem. Phys. Lett. 1993, 205, 535-542. (13) Pupyshev, V. I.; Panchenko, Y. N.; Bock, C. W.; Pongor, G. J. Chem. Phys. 1991, 94, 1247-1252. (14) Pople, J. A.; Scott, A. P.; Wong, M. W.; Radom, L. lsr. J. Chem. 1993, 33, 345-350. (15 ) (a) Pulay, P.; Fogarasi, G.; Pongor, G.; Boggs, J. E.; Vargha, A. J. Am. Chem. SOC.1983, 105,7037-7047. (b) Allen, W. D.; Csaszar, A. G.; Homer, D. A. J. Am. Chem. SOC. 1992, 114, 6834-6849. (c) Dutler, R.; Rauk, A. J. Am. Chem. SOC. 1989, 111, 6957-6966. (d) Schneider, W.; Thiel, W. J. Chem. Phys. 1987, 86, 923-936. (16) (a) Wiberg, K. B.; Waddell, S. T.; Rosenberg, R. E. J. Am. Chem. SOC. 1990, 112, 2184-2194. (b) Eggimann, T.; Ibrahim, N.; Shaw, R. A.; Wieser, H. Can. J. Chem. 1993, 71, 578. (c) Wiberg, K. B.; Walters, V. A.; Wong, K. N.; Colson, S. D. J. Phys. Chem. 1984, 88, 6067-6075. (17) Wolfsberg, M. Acc. Chem. Res. 1972, 5, 225. (18) Bigeleisen, J.; Mayer, M. G. J. Chem. Phys. 1947, 15, 261. (19) (a) Bron, J.; Chang, C. F.; Wolfsberg, M. Z. Naturforsch. 1972, 28A, 129-136. (b) Goodson, D. Z.; Sarpal, S. K.; Bopp, P.; Wolfsberg, M. J. Phys. Chem. 1982, 86, 659-663. (20) (a) Shiner, V. J., Jr.; Neumann, T. E. Z. Naturforsch. 1989, 44A, 337-354. (b) Hartshorn, S. R.; Shiner, V. J., Jr. J. Am. Chem. SOC.1972, 94, 9002-9012. (21) (a) Gabbay, S.; Rzepa, H. S. J. Chem. Soc., Faraday Trans. 2 1982, 671-677. (b) Hout, R. F., Jr.; Wolfsberg, M.; Hehre, W. J. J. Am. Chem. SOC. 1980, 202, 3296-3298. (c) Williams, I. H. J. Phys. Org. Chem. 1990,
J. Phys. Chem., Vol. 99, No. 40, I995 14699 3, 181-190. (d) Saunders, M.; Laidig, K. E.; Wolfsberg, M. J. Am. Chem. SOC. 1989, 111, 8989-8994. (22) (a) Bigeleisen, J. J. Chem. Phys. 1949,17, 675-678. (b) Saunders, W. H., Jr. In Investigation of Rates and Mechanisms of Reactions, 4th ed.; Bemasconi, C. F., Ed.; Wiley: New York, 1986; Part 1. (c) McLennan, D. J. In Isotopes in Organic Chemistry;Buncel, E., Lee, C. C., Eds.; Elsevier: Amsterdam, 1987; Vol. 7, Chapter 6. (d) Baldwin, J. E.; Reddy, V. P. J. Am. Chem. SOC. 1988, 110, 8555-8556. (e) Halevi, E. A,; Wolfsberg, M. J. Chem. Soc., Perkin Trans. 2 1993, 1493-1496. (f) Glad, S.; Jensen, F. J. Chem. SOC.,Perkin Trans. 2 1994, 871-876. (23) Frisch, M. J.; Trucks, G. W.; Head-Gordon, M.; Gill, P. M. W.; Wong, M. W.; Foresman, J. B.; Johnson, B. G.; Schlegel, H. B.; Robb, M. A.; Replogle, E. S.; Gomperts, R.; Andres, J. L.; Raghavachari, K.; Binkley, J. S.; Gonzalez, C.; Martin, R. L.; Fox, D. J.; Defrees, D. J.; Baker, J.; Stewart, J. J. P.; Pople, J. A. Gaussian 92, Revision C; Gaussian, Inc.: Pittsburgh, PA, 1992. (24) Hehre, W. J.; Radom, L.; Schleyer, P. v. R.; Pople, J. A. Ab Initio Molecular Orbital Theory: Wiley-Interscience: New York, 1986. (25) 3-21g* basis: (a) Binkley, J. S.; Pople, J. A,; Hehre, W. J. J. Am. Chem. SOC.1980, 102,939. (b) Gordon, M. S.; Binkley, J. S.; Pople, J. A.; Pietro, W. J.; Hehre, W. J. J. Am. Chem. SOC. 1982, 104, 2797-2803. (c) Pietro, W. J.; Francl, M. M.; Hehre, W. J.; Defrees, D. J.; Pople, J. A,; Binkley, J. S. J. Am. Chem. SOC. 1982, 104, 5039-5048. (26) 6-31g basis: (a) Hehre, W. J.; Ditchfield, R.; Pople, J. A. J. Chem. Phys. 1972,56,2257. (b) Francl, M. M.; Pietro, W. J.; Hehre, W. J.; Binkley, J. S.; Gordon, M. S.; DeFrees, D. J.; Pople, J. A. J. Chem. Phys. 1982, 77, 3654-3665. (c) Hariharan, P. C.; Pople, J. A. Theor. Chim. Acta 1973.28, 213. (27) 6-311g basis: (a) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. A. J. Chem. Phys. 1980, 72,650-654. (b) McLean, A. D.; Chandler, G. S. J. Chem. Phys.1980, 72, 5639-5648. (28) Frisch, M. Gaussian 90 User's Guide and Programmer's Reference, Revision J.; Gaussian, Inc.: Pittsburgh, PA, 1990 and references therein. (29) (a) Multiple polarization functions for 6-31g and 6-31 l g basis sets: Frisch, M. J.; Pople, J. A.; Binkley, J. S. J. Chem. Phys. 1984, 80, 3265-3269. (b) Diffuse functions: Clark, T.; Chandrasekhar, J.; Spitznagel, G.W.; Schleyer, P. v. R. J. Comput. Chem. 1983,4, 294-301. (c) Diffuse functions: Ohno, K.; Ishida, T. Int. J. Quantum Chem. 1986, 29, 677688. (30) (AM1): Dewar, M. J. S.; Zoebisch, E. G.; Healy, E. F.; Stewart, J. J. P. J. Am. Chem. SOC. 1985, 107, 3902-3909. (31) (MNDO): Dewar, M. J. S.; Thiel, W. J. Am. Chem. SOC. 1977,99, 4907. (32) (a) CH4: Gray, D. L.; Robiette, A. G. Mol. Phys. 1979,37, 19011920. (b) NH3: Duncan, J. L.; Mills, I. M. Spectrochim. Acta 1964, 20, 523-546. (c) N H 3 (frequencies in solid N2 matrix): Loutellier, A.; Perchard, J.-P. J. Mol. Struct. 1989, 198, 51-64. (d) H20 and H2S: Mills, I. M. In Theoretical Chemistry; Dixon, R. M., Ed.; The Chemical Society: London, 1974, Vol. 1, pp 152-153. (e) H20: Hoy, A. R.; Mills, I. M.; Strey, G. Mol. Phys. 1972, 24, 1265-1290. (f) H20: Sverdlov, L. M.; Kovner, M. A.; Krainov, E. P. Vibrational Spectra of Polyatomic Molecules; John Wiley: New York, 1974. (8) SiH4: Levin, I. W.; King, W. T. J. Chem. Phys. 1962, 37, 1375-1381. (h) S i b : Kattenberg, H. W.; Oskam, A. J. Mol. Spectrosc. 1974,49, 52-69. (i) PH3: Duncan, J. L.; McKean, D. C. J. Mol. Spectrosc. 1984, 107, 301-305. (i) H2S: Botschwina, P.; Zilch, A,; Werner, H.-J.; Rosmus, P.; Reinsch, E.-A. J. Chem. Phys. 1986, 85, 5107-5116. (k) HF and HCl: Radzig, A. A,; Smirnov, B. M. Reference Data on Atoms, Molecules, and Ions; Springer Series in Chemical Physics 31; Springer-Verlag: New York, 1985. (1) C2H2: Strey, G.; Mills, I. M. J. Mol. Spectrosc. 1976, 59, 103-115. (m) CzH2: Fast, H.; Welsh, H. L. J. Mol. Spectrosc. 1972,41,203-221. (n) CH2O: Wohar, M. M.; Jagodzinski, P. W. J. Mol. Spectrosc. 1991, 148, 13-19. (0)C2H4: Duncan, J. L.; Hamilton, E. J. Mol. Struct. (THEOCHEM) 1981, 76, 65-80. (p) C2H4: Chough, S. H.; et al. J. Mol. Strucr. 1992,272, 179-186. (4) C2&: Duncan, J. L.; Robertson, G. E. J. Mol. Spectrosc. 1991, 145, 251-261. (r) HCN: Nakagawa, T.; Morino, T. Bull. Chem. SOC.Jpn. 1969,42,2212-2219. (s) HCN: Quapp, W. J. Mol. Spectrosc. 1987,125, 122-127. (t) HCN: Strey, G.; Mills, I. M. Mol. Phys. 1973, 26, 129-138. (u) HCN: Smith, A. M.; Coy, S. L.; Klemperer, W.; Lehmann, K. K. J. Mol. Spectrosc. 1989, 134, 134-153. (v) DCN: Suzuki, I.; Pariseau, M. A.; Overend, J. J. Chem. Phys. 1966, 44, 3561. (w) CHsCl and CH3F: Duncan, J. L.; McKean, D. C.; Speirs, G. K. Mol. Phys. 1972, 24, 553-565. (x) CHsC1: Duncan, J. L.; Law, M. M. J . Mol. Spectrosc. 1990, 140, 13-30. (y) CH3C1: Duncan, J. L.: Allan, A,; McKean, D. C. Mol. Phys. 1970, 18, 289-303. (z) c - C ~ H ~ : Blom, C. E.; Altona, C. Mol. Phys. 1976, 31, 1377-1391. (33) (a) Michalska, D.; Schaad, L. J.; Carsky, P.; Andes Hess, B., Jr.; Ewig, C. S. J. Comput. Chem. 1988, 9, 495-504. (b) Pulay, P.; Fogarasi, G.; Pang, F.; Boggs, J. E. J. Am. Chem. SOC. 1979, 101, 2550. (34) Ignacio, E. W.; Schlegel, H. B. J. Comput. Chem. 1991, 12,751760. (35) Guo, H.; Karplus, M. J. Chem. Phys. 1989, 91, 1722. JP95 1469X
