J. Phys. Chem. 1995,99, 3093-3100
3093
Transferable Scaling Factors for Density Functional Derived Vibrational Force Fields Guntram Rauhut" and Peter Pulay* Department of Chemistry and Biochemistry, University of Arkansas, Fayetteville, Arkansas 72701 Received: October 25, 1994; In Final Form: December 13, 1994@
Density functional theory (DFT) using the 6-31G* basis set and two nonlocal exchange-correlation functionals (Becke-Lee-Yang-Par [B-LYP] and the three-parameter compound function of Becke [B3-LYP]) has been used for the calculation of vibrational force fields of a set of 31 organic molecules including a wide range of functional groups. The calculated force constants have been scaled to experimental vibrational frequencies by using (a) an overall scaling constant and (b) a set of 11 factors paying respect to the different kinds of internal coordinates. The comparison of the scaled fundamental frequencies with experiment shows that density functional theory is a reliable tool for the interpretation of IR spectra. The uncorrected DFT frequencies and force constants approximate the experimental ones in a much more uniform fashion than does Hartree-Fock theory. Nevertheless, the use of multiple scale factors leads to further significant improvement. The scaled B3-LYP results are superior to the B-LYP ones, even though the unscaled B-LYP frequencies are, through error cancellation, slightly better than the B3-LYP ones. The reliability of scaled force fields is demonstrated by comparing the calculated and experimental vibrational spectra of aniline.
Introduction In the past two decades, ab initio MO theory has become an important tool for the interpretation of vibrational spectra.'-4 In a number of cases, assignments of fundamentals could be clarified only by using ab initio calculations. Numerous unusual molecules have been identified by comparing calculated vibrational spectra with observed ones. These developments were made possible by the analytical calculation of the frst5 and later the second6 derivatives of the molecular potential energy. Most calculations have been performed at the ab initio SCF level. At this level, the calculated harmonic force constants and frequencies are usually higher than the corresponding experimental quantities, due to a combination of electron correlation effects and basis set deficiencies. The overestimation of the frequencies becomes more severe if the calculated harmonic frequencies are compared with the observed fundamentals,as anharmonicity usually lowers the frequencies. It would be preferable to compare calculated harmonic frequencies with experimental harmonic frequencies, or calculated anharmonic fundamentals with experimental fundamentals. However, neither the calculation of anharmonic vibrational energy levels, nor the experimental determination of harmonic frequencies is routinely practicable for polyatomic molecules at present. As the errors in the calculations (including a significant part of anharmonicity corrections) are largely systematic, they can be corrected by empirical scaling factors. The first suggestion of this nature7 recommended scaling double-5 stretching force constants by 0.9 and bending ones by 0.8. A more systematic scaling procedure has been developed by Blom and Altona,8 also using multiple scale factors in internal coordinate representation. These authors, however, had no access yet to analytical derivative programs and were therefore restricted to relatively small molecules. The procedure of Blom and Altona served as the basis for the scaled quantum mechanical (SQM) force field pro~edure.~ In this model, the molecular force field is expressed in a set of standardized valence internal coordinates.'0." Usually, a set of molecules containing similar motifs are treated together. On the basis of chemical intuition, the @
Abstract published in Advance ACS Abstracts, February 15, 1995.
internal coordinates of all the molecules are sorted into groups sharing a common scaling factor, and the factors for each group are determined by a least-squares fitting to experimental vibrational frequencies. A more systematic procedure for the grouping of coordinates has not yet been attempted because of the large number of combinations; optimization by a neural network is a possible technique which has proved useful for similar problems.'* The quadratic force constants are scaled by the equation
where Fv is the original force constant, Fv is the scaled one, and Ai, 2, are the scale factors. This scaling procedure is invariant against mixing coordinates with the same scaling factor and is thus somewhat preferable to earlier suggestion^.^-^.'^ Experience has shown that a moderate number of separate factors (about 10) suffice to correct the Hartree-Fock force field of most organic molecules to produce good agreement with experimental spectra. In some strongly correlated systems, e.g. benzene, further empirical corrections are necessary at the Hartree-Fock 1 e ~ e l . l ~Obviously, the goal is to obtain an accurate fit to the experimental frequencies with only a few empirical constants. A disadvantage of all selective scaling procedures is the necessity to transform the force constants to nonredundant valence coordinates. The definition of the latter by hand is a tedious and error-prone procedure for larger molecules. Transformation is not needed if a uniform scaling factor is applied, as this can be applied to Cartesian force constants, or simply to the final frequencies. In spite of the fact that the physical origin of the errors in different types of force constants (stretchings, deformations, and torsions) is different, and so are the optimized scale factors, overall scaling is able to produce reasonably accurate force constants by a simple p r ~ c e d u r e . ' ~For ~ ' ~most popular basis sets, the best scale factor for frequencies is closeI7 to 0.9 (e.g. it is 0.89 for the 6-31G* basisI7) and therefore for force constants is 0.92 = 0.81. The accuracy obtained by an overall scaling is naturally lower than selective scaling. Improvements in the quantum chemical method (e.g. introduction of correlation effects via second-order Moller-Plesset theory,
0022-3654/95/2099-3093$09.00/0 0 1995 American Chemical Society
3094 J. Phys. Chem., Vol. 99, No. IO, 1995
Rauhut and Pulay
improvement resulting from multiple scale factors versus a MP2) yield scale factors closer to unity, although it is worth single factor versus no scaling. pointing out that scaling factors fitted to observed (anharmonic) frequencies will deviate from unity even for exact calculations. Calculations For the calculation of zero-point energies and thermodynamic parameters, single scale factors are adequate. However, sigAll calculations have been carried out by using the 1993 nificantly better reproduction of the frequencies can be achieved release (G92-DFT) of the GAUSSIAN suite of programs.@The by using multiple scale factors (for a recent example, see e.g. 6-3lG* basis see' has been used throughout as a compromise ref 18). Due to the better reproduction of the forms of the between accuracy and applicability to large molecules. As normal modes, multiple scaling also improves infrared intensipointed out in the introduction, smaller basis sets suffice in DFT ties. theory because the basis functions do not have to describe the Methods for automatically generating natural internal coorcorrelating orbitals. All calculations were carried out twice, dinateslOT'lhave been developed recently in our group," as well using (a) Becke's exchange functional24 in combination with as in other group^.'^,^^ They remove the difficulties associated the Lee, Yang, and Parr (LYP) correlation functional25 (from with selective scaling procedures for the prediction of vibrational now on B-LYP) and (b) Becke's three-parameter exchange spectra. Using one of these programs, the calculation and functional42in combination with the LYP correlation functional scaling of quantum mechanical force fields can be largely (from now on B3-LYP) as implemented in G92-DlT. The exact automated. The only tasks requiring chemical intuition and definition of the B3-LYP exchange-correlation functional is spectroscopic knowledge are the grouping of the scaling factors and the assignment of the observed fundamentals, based on the E,, = .ET (1 - a ) E Y bE,Bg8 combined spectroscopic and calculational evidence. cEcLYp (1 - c)E,VWN Although much of the error of Hartree-Fock theory for vibrational force fields can be successfully corrected by using where42a = 0.8, b = 0.72, c = 0.81, and the superscripts DS, empirical scaling, this method fails for strongly correlated HF, B88, LYP, and VWN refer to the Dirac-Slater, Hartreesystems (e.g. many organometallics), and its accuracy is limited Fock, B e ~ k eL, e~e ~- Y a n g - P a ~ ~ ,and ~ ~ Vo~ko-Wilk-Nusair~~ in many other cases. Traditional correlation methods are too exchange-correlation functionals. Force constants have been expensive for the large molecules encountered in vibrational calculated at fully optimized geometries using the analytical spectroscopy. Density functional theory (DFT)in the Kohnsecond derivative capability of G92-DFT. We have automated Sham formulation (see e.g. refs 21 and 22) has emerged in the the vibrational calculations in the following way: we extract past few years as a successful alternative to traditional methods the geometries and Cartesian force constants from G92 archive using configuration expansion. The recent popularity of DFT files and use them as input in our vibrational program package. is due to several factors: the introduction of new, accurate The latter generates the natural internal coordinates, transforms exchange-correlation potentials, e.g. refs 23-26, which yield the force constants, calculates the G matrix,43and prepares the substantially improved predictions; the introduction of analytical data for the scaling program, a modified'"' version of the first derivative^;*^-^^ the recent f o r m u l a t i ~ n ~and ~ - ~imple~ SCALE2 code.45 mentation33 of analytical second derivatives; and finally the (a) The Training Set. A set of 20 small molecules development of p r o g r a n - 1 ~ ~which ~ 3 ~ ~approach the numerical ( b e n ~ e n e ' (l), ~ . ~ethene47,48 ~ (2), methanol49(3), formaldehyde9 precision of traditional quantum chemistry methods. An (4), glyoxal9 (5), trans-butadiene9 (6), acroleing (7), dimethyl important advantage of DFT techniques is that the convergence ether50 (8), methylamine5' (9), c y ~ l o p r o p a n e(lo), ~ ~ propene53 with respect to improvements in the basis set is much faster (ll),pyridine54 (12), hydrazine55 (13), acetic acid56 (14), than in traditional correlation techniques since the basis set a ~ e t i d i n (15), e ~ ~ a ~ e t o n i t r i l e(16), ~ ~ pyrrole59(17), furanm (18), serves only to describe the strongly occupied molecular orbitals, formic acid6' (19), formic acid methyl esteS6 (20)), including not the highly oscillatory correlation orbitals. a variety of functional groups and limited to the atoms H, C, Recent evidence shows that modem and n ~ n l o c a l ~ ~ ~N, ~ *and 0, was used to determine scaling factors for the most exchange-correlation potentials produce surprisingly accurate common organic structural motifs. A representative selection molecular force fields and vibrational frequencies. In particular, of geometries of molecules of this set is shown in Table 1. The the recent results of the Pople and Handy g r o ~ p sshow ~ ~ ,that ~~ calculated structures are obviously in good agreement with the nonlocal DFT for molecular force fields is much superior to experimentally obtained values. As already discussed by Hartree-Fock theory. Therefore, we have used two nonlocal Johnson et the B-LYP functional leads to bond lengths functionals in our study here, despite the fact that these methods which are systematically too long, particularly for the C-H bond are slightly slower than the local ones. Nevertheless, the length. The B3-LYP method leads to geometry parameters calculations exhibit systematic errors and thus benefit from which are much closer to experimental data. It still has a slight scaling. This was also shown recently in a study of formatendency to overestimate the bond lengths, but the errors are mide.39 In the present work we would like to establish a on the average less than 0.005 A, and thus they are of the same minimum number of transferable scaling factors for a variety order as the uncertainties in the experimental equilibrium (re) of organic compounds, unlike the in-depth study of a single structures for most polyatomic molecules. The calculated bond molecule.39 We derive transferable DFT scaling factors for angles are similarly accurate to within a few tenths of a degree organic molecules containilig common structural motifs. These in most cases. factors are based on a pool of 347 experimental fundamentals In a first step of the scaling procedure the uncorrected in 20 simple molecules. We paid particular attention to ensure vibrational frequencies have been calculated and compared with experimental data. The experimental frequencies, based on the that the fundamentals used to optimize the scale factors are original literature and later theoretical analyses, have been securely assigned. We have carefully collected previous carefully screened to eliminate misassignments and uncertainties spectroscopicevidence and in some cases results from high level in the interpretation of these complex spectra. Fundamentals ab initio calculations. The optimized factors are used to predict which are uncertain or not directly observed were not used as the 309 fundamental vibrational frequencies and infrared spectra input data to the scaling procedure. Examples of such fundaof 11 medium-sized test molecules. We try to assess the
+
+
+
+
J. Phys. Chem., Vol. 99, No. IO, 1995 3095
Scaling Factors for Force Fields TABLE 1: Molecular Geometries Obtained by the B-LYP and B3-LYP Functionals (Distances, di; Angles, deg)
2000 h
molecule acetonitrile
propene
dimethyl ether
benzene furan
glyoxal
coordinate B-LYP
error
B3-LYP
error
exptd
1.469 1.173 1.102 110.5 1.344 1.511 1.106 1.103 1.099 1.096 1.094 125.4 1.426 1.111 1.101 11 1.9 112.0 106.9 1.407 1.094 1.381 1.372 1.443 1.086 1.088 106.5 110.5 106.3 115.3 127.1 1.223 1.537 1.119 121.8 123.9
0.01 1 0.016 -0.002 1.00 0.008 0.010 0.008 0.018 0.009 0.005 0.013 1.09 0.0 16 0.011 0.010 0.20 1.19 -0.30 0.011 0.01 1 0.019 0.011 0.012 0.011 0.011 -0.05 -0.24 0.24 -0.60 -0.86 0.021 0.0 10 0.010 0.65 0.55
1.461 1.160 1.094 110.3 1.333 1.502 1.099 1.095 1.091 1.089 1.087 125.3 1.410 1.103 1.093 112.4 111.8 107.2 1.396 1.087 1.364 1.361 1.436 1.079 1.08 1 106.8 110.5 106.1 115.6 127.3 1.209 1.525 1.109 121.5 123.8
0.003 0.003 -0.010 0.82 -0.003 0.001 0.001 0.010 0.001 -0.003 0.006 0.96 0.00 0.003 0.002 0.65 1.02 0.03
1.458 1.157 1.104 109.5 1.336 1.501 1.098 1.085 1.090 1.091 1.081 124.3 1.410 1.100 1.091 111.7 110.8 107.2 1.396 1.083 1.362 1.361 1.431 1.075 1.077 106.6 110.7 106.1 115.9 128.0 1.202 1.527 1.109 121.2 123.4
0.000 0.004 0.002
0.000 0.005 0.004 0.004 0.25 -0.18 0.04 -0.25 -0.70 0.007 0.002 0.000 0.35 0.45
Out of plane. In plane. Central atom. Experimental data taken from the following: Costain, C. C. J . Chem. Phys. 1958, 29, 864 (acetonitrile). Lide, D. R.; Christensen, D. J . Chem. Phys. 1961, 35, 1374 (propene). Kasai, P. H.; Myers, R. J. J . Chem. Phys. 1959, 30, 1096L (dimethyl ether). Cabana, A.; Bachand, J.; GiguCre, J. Can. J . Phys. 1974,52, 1949 (benzene). Nosberger, P., Bauder, A.; Giinthard, Hs. H. Chem. Phys. 1973, I , 418 (furan). Birss, F. W.; Braund, D. B.; Cole, A. R. H.; Engleman, R.; Green, A. A.; Japar, S . M.; Nanes, R.; Om, B. J.; Ramsay, D. A.; Szyszka, J. Can. J . Phys. 1977, 55, 390 (glyoxal).
5
> v
ffl
0
&.
1500-
E U t LL 3
D
+
0 3
5 10000 U
0 V
t
3
500 t 500
'
l
1000
1500 Experimental Fundamentals ( 1 / c m )
2000
Figure 1. Unscaled calculated (B3-LYP) vibrational frequencies in comparison to the experimentally obtained data. All units in cm-'.
TABLE 2: Basis Set Effects on the Harmonic Force Constants of Water, Ammonia, and Methane" basis set
6-31G*
6-31G**
6-31 lG(2df,p)
8.030 -0.134 0.7754 0.256
HzO 8.757 8.327 -0.113 -0.074 0.695 1 0.7083 0.268 0.266
6.908 -0.044 0.5714 0.7327 -0.1489
7.002 -0.040 0.5316 0.7029 -0.1536
5.177 0.036 0.5813 0.5383 0.100
5.171 0.03 1 0.5700 0.525 1 0.100
expt 8.454b -0.101 0.697 0.237
3"
(I
mentals are the out-of-plane CH wagging of acrolein (7)9which has not been observed, only estimated from the product rule by Brand62 to be at 980 cm-I; the symmetric CH2 stretching in trans-butadiene (6)9at 2985 cm-I; or the 972 and 1335 cm-' a" vibrations of methylamine (9).51 Most of these neglected fundamentals are torsions and C-H stretchings which often show significant uncertainties. The correlation between the calculated (B3-LYP) and observed results for the 241 securely assigned fundamentals in the fingerprint region (500-2000 cm-l) is shown in Figure 1. Due to the large volume of data, the detailed tables of calculated and experimental frequencies and infrared intensities are available only in the supplementary material. This figure shows that B3-LYP at the 6-31G* level systematically overestimates the vibrational frequencies by about 5%. This is partly due to anharmonicity and partly to the general tendency of most quantum chemical methods to overestimate force constants at the exact equilibrium geometry. Deformations involving central atoms with lone pairs, e.g. the COH deformation in formic acid (19)61at 1223 cm-' or the N H 2 scissoring modes in hydrazine (13)55near 1600 cm-I, show larger deviations. The main reason for this appears to be that the 6-31G* basis set is not quite adequate for the description of the NH2 or OH group. To assess the basis set effects, Table
r
6.935 -0.017 0.5271 0.6857 -0.1562
7.05 f 0.13' 0.01 f 0.19 0.532 f 0.02 0.665 f 0.01 -0.142 f 0.03
CH4 5.086 0.041 0.5545 0.5133 0.103
5.392 f O.Old 0.014 f 0.01 0.584 f 0.01 0.548 i 0.01 0.110 f 0.01
All calculations use the B3-LYP model. Force constants in aJ/.&*, rad) and aJ/rad2, respectively, for stretching, stretch-bend couplings, and bendings. See the experimental data for the coordinate definitions. Hoy, A. R.; Mills, I. M.; Strey, G. Mol. Phys. 1972, 24, 126. Bunker, P. R. J . Mol. Spectrosc. 1979, 74, 1. Duncan, J. L.; Mills, I. M. Spectrochim. Acta 1964, 20, 523. Gray, D. L.; Robiette, A. G. Mol. Phys. 1979, 37, 1901.
aJ/@
2 compares the experimental harmonic force constants of water, ammonia, and methane with those calculated at the B3-LYP level using different basis sets. It is evident that the B3-LYP model yields excellent harmonic force constants in the basis set limit but the 6-31G* deformational constants are too high in water and ammonia. The 6-31G** basis is superior in this respect to the 6-31G*. The root mean square deviation of all vibrations in the fundamental range (0-4000 cm-I) is 29.6 cm-l for the B-LYP method and 74.1 cm-' for the B3-LYP functional. At first, the better agreement of the B-LYP frequencies is surprising, in view of the fact that B3-LYP reproduces the geometries better. This is, however, due to error cancellation: the B-LYP method overestimates the bond lengths, which in turn leads to a decrease of force constants,63 cancelling the general overestimation in the observed frequencies. This is especially important for CH vibrations for which anharmonicity is large, but so are the errors
3096 J. Phys. Chem., Vol. 99, No. 10, 1995
Rauhut and Pulay
TABLE 3: Root Mean Square Deviations (em-') of the Calculated Fundamentals of the Training Set Using Different Exchange Correlation Functions B-LYP
unscaled one factor SQM scaled
B3-LYP
0-500 cm-l
501-2500 cm-I
2501-4000 cm-l
0-4000 cm-I
0-500 cm-I
501-2500 cm-'
2501-4000 cm-'
0-4000 cm-I
22.0 22.4 16.4
22.1 22.5 13.9
46.5 35.8 30.1
29.6 26.2 19.1
20.3 19.0 12.5
45.8 17.0 10.8
131.7 22.2 17.7
74.1 18.5 12.8
TABLE 4: Scaling Factors for the Different Types of Internal Coordinates. (A) For the Training Set of 20 Molecules and (B) for the Training Set and the Test Set (a Total of 31 Molecules)" BLYP
2000
1 ,
B3-LYP
no.
mode
A
B
A
B
7 8 9 10 11
X-Yb stretching X-H stretching XYZ bending XY -H bending H-X-H bending out-of-plane modes NH2 wagging XO-H, XN-H bending torsions of conjugated systems torsions of single-bonded systems linear deformations
1.007 0.977 1.052 1.005 0.964 1.072 0.834 0.980 0.869 0.990 0.986
1.013 0.976 1.057 1.002 0.963 1.062 0.808 0.985 0.819 1.002 0.973
0.922 0.920 0.990 0.950 0.915 0.976 0.806 0.876 0.831 0.935 0.913
0.923 0.919 0.995 0.951 0.913 0.975 0.788 0.887 0.796 0.942 0.900
Note that set B was not used in any further calculations. X, Y, and Z denote heavy atoms.
in the B-LYP bond lengths. The mean deviations in three spectral ranges, 0-500, 501-2500, and 2501-4000 cm-' are shown in Table 3. In a second step an overall scaling factor has been evaluated by using a least squares optimization of the calculated fiequencies to the experimental ones. The least squares procedure used the default weighting described in ref 9. The scaling factor for this training set is 0.990 for the B-LYP method and 0.928 for the B3-LYP functional; these correspond to frequency scale factors of 0.995 and 0.963. The latter values can be compared with the frequency scaling factors of 0.89 for the Hartree-Focld 6-31G* method and appr~ximately'~ 0.95 for MP2/6-31G*, showing that DFT-derived force constants need much less correction than the SCF ones and even less than the MP2 ones. Due to the error cancellation discussed above, the B-LYP method hardly needs scaling but the mean deviation after scaling, 26.2 cm-I, is worse than for B3-LYP (18.5 cm-I). The deviations in the three spectral ranges are included in Table 3. In a third step a set of 11 scaling factors has been chosen for different types of internal coordinates (Table 4). The grouping of the coordinates is obvious for the most part, and we will make only a few comments. We found it necessary to treat separately the NH2 waggings and the XOH, X" bendings (X = heavy atom). As discussed above, the corresponding force constants have large errors at the B3-LYFV6-31G* level, and these modes, particularly the NH2 wagging, have also large anharmonicity. The CH stretching scale factor is close to the other stretchings. However, we treat it separately in order to be able to predict C-D vibrations, using a slightly higher factor to account for the difference in the effective force constants. With use of the optimized scaling factors in Table 4 (sets A), the root mean square deviation decreased to 19.1 cm-' for the B-LYP functional and to 12.8 cm-' for B3-LYP (see Table 3). The improvement is even more significant in the fingerprint region, with deviations decreasing from 22.1 to 13.9 cm-I (BLYP) and from 45.8 to 10.8 cm-' (B3-LYP). After scaling, the B3-LYP method performs uniformly better than B-LYP. The improvement of the SQM method over uniform scaling is
500 500
1000
1500
2000
Experimental Fundamentals ( 1 /cm)
Figure 2. B3-LYP calculated vibrational frequencies using the scaling factors of the SQM method in comparison with experimental data. All units in cm-I.
also significant, close to a factor of 1.6 in the fingerprint region. For reliable predictions it is perhaps even more important that the number of large deviations is significantly reduced, as evident from Figure 2 which compares the positions of 241 fundamentals of the training set in the fingerprint region with experimental data. It is obvious that the SQM method in combination with nonlocal DFT theory is able to reproduce the fundamentals of small organic molecules reliably. Figure 3 shows the deviations in the fundamental frequencies for the individual molecules in the training set. Figure 3 and Table 2 indicate also larger errors for the XH stretching region between 2500 and 4000 cm-'. This may be mainly caused by uncertainties in the assignments and Fermi resonances in this region which cannot be corrected for by scaling. (b) The Test Set. To check the reliability of the fitting, a test set of 11 medium sized molecules was chosen: aniline64 (21), ethanoP5 (22), oxetane66(23), p~ru-benzoquinone~~ (U), nitrobenzene6*(25), phenoF9 (26), cis-furan-2-aldehyde70(27), benzaldehyde7I (28), b e n z ~ n i t r i l e(29), ~~ (30), and acetone74 (31). The assignments for these molecules are somewhat less certain than for the training set. We have tried to exclude from the comparison fundamentals which are not securely assigned, for instance the a-rocking at 840 cm-I (Kydd's assignment) in oxetane (23)66or the CH stretching at 3084 cm-' in uracil (30).73 Even so, some of the weaker fundamentals are uncertain experimentally. The spectra of these molecules, including 297 fundamentals, have been calculated by using (a) no scaling constants, (b) the overall scaling constant of the training set, and (c) the 11 SQM scaling factors of the training set (set A in Table 4). The results are summarized in Table 5 and Figure 4. Table 5 shows that the accuracy of the calculated results for
Scaling Factors for Force Fields 150.00
100.00
3
J. Phys. Chem., Vol. 99, No. IO, 1995 3097 8-LYP
2500-4000
83-LYP
1/cm
1
-
oeboo Unscaled LY. One
Foclor
SOM Scaled
50.00
-One
Unscolcd Foctor SaM Scaled
E
>
>
v
2
Oeeea
50.00
E L
-
100.00
V
0.00
0.00
L
500-2500
e
l/cm
50.00
50.00
c" 2
VI
[L
CY
3 I
I
0-500
I
,
I
,
,
,
,
,
,
,
,
,
,
,
,
,
,
1
9 4 ,
I
0-500
P
0.00
0.00
,
l/cm
=zG
I
,
,
,
,
5
0
50.00
. 6
,
,
I
,
10 Molecule
,
I
,
(
,
,
a
I
0.00
20
15
l/cm
5
0
15
10 Molecule
20
Figure 3. Root mean square deviations (cm-l) for all molecules of the training set using different scaling methods (benzene (l), ethene (2), methanol (3), formaldehyde (41, glyoxal (9,trans-butadiene (6), acrolein (7), dimethyl ether (8), methylamine (9), cyclopropane (lo), propene ( l l ) , pyridine (12), hydrazine (13), acetic acid (14), azetidine (15), acetonitrile (16). pyrrol (17), furan (l8), formic acid (19), and formic acid methyl ester (20)). The left diagram shows the values for the B-LYP method, the right for the B3-LYP method.
TABLE 5: Root Mean Square Deviations (cm-') of the Calculated Fundamentals of the Test Set Using Different Exchange Correlation Functions B-LYP
0-500 unscaled one factor SQM scaled" a
B3-LYP
501-2500 cm-'
2501-4000 cm-'
0-4000 cm-'
0-500
cm-I
cm-'
501-2500 cm-I
20.9 21.2 14.8
20.2 22.8 13.7
52.7 41.9 35.0
28.6 26.9 19.2
18.7 18.3 9.4
39.1 17.5 12.4
2501-4000 cm-I
0-4000 cm-'
136.3 21.5 13.9
65.9 19.7 12.6
Using the scaling factors of the training set (set A) in Table 4.
-
B-LYP 2500-4000
l/cm
eeew Unscalcd
One Foclor SOM Scaled
80.00
150.00
P 100.00
40.00
1
B3-LYP 2500-4000
l/cm
-
j
Oeoeo
Unscolcd
One Foclor 5PM Scaled
i
5 20.00 W v)
220.00
E
40.00
1
4 0-500
1
m 50.00
l/cm
(L
0.00
A
l/cm
500-2500
I
50.00
0-500
l/cm
20.00
0.00
I
I
21
~
23
I
~
I
25 27 Molecule
J
29
I
31
~
I
0.00 J
I
21
23
25
27 Molecule
29
31
Figure 4. Root mean square deviations (cm-I) for all molecules of the test set using different scaling methods (aniline (21), ethanol (22). oxetane (23), p-benzoquinone (24), nitrobenzene (25), phenol (26), cis-furan-2-aldehyde (27), benzaldehyde (28), benzonitrile (29), uracil (30), and acetone (31)). The left diagram shows the values for the B-LYP method, the right for the B3-LYP method.
the test set is similar to the one of the training set, confirming the reliability of the correction of systematic errors. In a last step we have checked the stability of the scaling constants by increasing the training set by the test set and
optimizing the scaling factors for all 31 molecules (644 fundamentals). The values obtained (see Table 4,sets B) differ for most factors only insignificantly from the values of the training set given in Table 4 (sets A); only the scaling factors
3098 J. Phys. Chem., Vol. 99, No. IO, 1995
Rauhut and Pulay
Wavenumbers
1
Wavenumbers
Figure 5. (Top) Experimental spectrum of the Aldrich FT-IR library for aniline (21). Reprinted with permission from ref 76. Copyright 1989 Aldrich. (Bottom) Simulated spectrum using the SQM method within the framework of B3-LYPDFT. For all spectra frequencies are given in cm-I. Note the change of scale on the frequency axis.
for the out-of-plane modes, the torsions of conjugated systems, and the torsions of single-bonded systems show bigger deviations. This may be mainly an effect of the large number of aromatic molecules of the test set, whereas the training set includes only a few aromatic structures. The root mean square deviation in the range 0-4000 cm-' shows only small changes: the values are 19.1 cm-' for B-LYP and 13.4 cm-' for B3-LYP, compared with 19.1 and 12.8 cm-I. (c) Examples. The accuracy of the method presented will be illustrated by the simulation of the spectrum of aniline (21). Because of its highly anharmonic NH2 wagging vibration, this is a more challenging molecule than a typical rigid structure. EvansM and co-workers measured the spectrum of aniline and its deuterated derivatives in 1960 and assigned the fundamentals on the basis of the gas-phase, solution, and liquid spectra. Niu et ~ 1 calculated . ~ ~the vibrational spectrum of aniline by using the Hartree-Fock 4-21G method in combination with the SQM method and reassigned some fundamentals. Aniline is included in our test set both at the B-LYP and the B3-LYP levels, but only the spectrum obtained by the latter will be discussed. Our simulated spectrum (SQM) and the experimental gasphase spectrum of the Aldrich library of FT-IR spectra76are shown in Figure 5. The good agreement between the two for the majority of the bands is evident; the intensities are also very satisfactory. The CH stretching region and the fingerprint region between 1600 and 700 cm-' are very well reproduced in the simulated spectrum. The most obvious deviation between the calculated and the experimental spectrum is the position of the bands associated with the strongly anharmonic N H 2 wagging vibration, corresponding to the ammonia inversion. Both the ground and the first vibrational levels are split due to the double minimum potential. The strong band near 660 cm-l in the gas phase which appears to be the wagging fundamental is just one component of it, the 1 3 (or 0I-) Obviously, a satisfactory description of this mode requires an
-
-
-
-
anharmonic analysis. The experimental value quoted for this fundamental in Table 6 is the average of the 1 3 and 0 2 transitions. A medium band at about 1100 cm-' in the gasphase spectrum is absent in the calculated spectrum, nor is it present in the liquid phase or in solution. It appears to be one of the higher transitions of the NH2 wagging fundamental; its position agrees well with the estimate of Quack and Stockburger7' given for the 0 4 (O+ 2+) transition, 1086 cm-'. Table 6 gives the assignments of the fundamentals based on the present calculations and on the experimental spectrum of Evans.@ The calculations largely confirm the assignment of Evans. They also agree quite well with those of Niu et ~ 1 . with one major discrepancy (see below). We believe that our present, significantly higher level calculations should be more reliable than those of ref 75. Our calculated value, 261 cm-I, for the NH2 torsional mode is in good agreement with the assignment of Larsen et ~ 1 . (277 7 ~ cm-I) and is substantially higher than the 216 cm-I value inferred by Evans@ from a tentative assignment of the first overtone of this mode at 420 cm-I. The assignment of the two close-lying bands at 808 and 823 cm-' cannot be conclusively decided. One of them is a C-N stretching-ring stretching combination; the other is an out-of-plane CH deformation. Both are weak in the infrared and equally strong in the Raman spectrum. Our values (8 10 and 8 18 cm-') support the original assignment of Evans. There is one major discrepancy between the calculations of Niu et ~ 1 and. ours. ~ ~These authors assign, on the basis of their calculations and at variation with Evans, the very weak combination band at 1190 cm-' to a fundamental CH in-plane bending mode. Our calculations are in almost perfect agreement with Evans' experimental results and thus strongly support the latter. We suspect a simple clerical error in ref 75, as the discrepancy is limited to two CH in-plane deformations only, and the 4-21G results, after scaling, should not differ this much for this vibration. Note that the force constants of ref 75 have been hand calculated by finite
- -
~ 5
Scaling Factors for Force Fields TABLE 6: Calculated (BfLYP), SQM Scaled, and Unscaled Vibrational Frequencies (cm-') and Infrared Intensities (Wmol) of Aniline (21)" no. 1
2 3 4 5 6 7 8
9 10
11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
symunassignment metry expt scaled dev intens ring torsion a' 233 217.7 -15.3 5.5 224.7 NHz torsion 22.7 285.9 a" 277 260.7 -16.3 0.2 383.5 CN rock a" 390 380.5 -9.5 ring torsion 0.3 418.1 a" 415 402.9 -12.1 CN def, op, ring torsion a' 501 492.9 -8.1 84.3 508.9 ring def, ip a' 526 529.8 3.8 9.7 536.9 NH2 wag a' (541) 576.7 35.3 270.9 632.4 ring def, ip 0.3 635.9 a" 619 630.2 11.2 ring torsion 6.7 704.6 a' 689 682.5 -6.5 5.4 70.4 763.9 CH bend, op a' 747 752.4 CN str, ring str 808 809.7 1.7 2.9 835.5 a' 823 817.9 -5.1 CH bend, op a" 0.1 827.3 9.6 878.5 874 869.7 -4.3 CH bend, op a' CH bend, op a" 0.0 947.4 957 938.9 -18.1 968 968.3 CH bend, op a' 0.3 0.3 975.0 2.4 1012.0 996 994.4 -1.6 ring def, ip a' ring breath 1.1 1057.5 1028 1020.7 -7.3 a' 2.4 1079.8 1054 1044.9 -9.1 NH2 rock, ring str a" 2.2 2.8 1150.8 NHz rock, CH bend, ip a" 1115 1117.2 CH bend 1152 1156.5 a" 4.5 2.2 1189.9 a' 3.5 7.9 1211.7 CH bend, ip 1173 1176.5 a' CN st^ 1278 1271.1 -6.9 60.1 1316.4 a" ring str 7.5 1368.9 1308 1318.6 10.6 CH bend, ip a" 1340 1343.3 3.3 0.2 1377.8 ring str, CH bend a" 1468 1474.9 6.9 2.0 1519.7 0.2 52.2 1552.2 CH bend, ring str a' 1503 1503.2 ring str a" 1590 1593.5 3.5 4.5 1647.1 5.5 10.5 1668.4 ring str, NHz sciss a' 1603 1608.5 9.1 136.2 1695.4 NHz sciss, ring str a' 1618 1627.1 5.7 3170.1 CH str a" 3025 3039.9 14.9 CH str a' 3037 3038.7 1.7 20.1 3169.0 2.2 4.1 3186.1 CH str a' 3053 3055.2 CH str a' 3072 3060.7 -1 1.3 19.3 3209.5 CH str a" 3088 3077.6 -10.4 47.3 3191.8 NH2 str a' 3418 3399.8 -18.2 9.5 3545.6 NH2 str 7.2 3643.6 a" 3500 3493.8 -6.2
"Experimental values (see text) are based on ref 64. Solution frequencies are listed, as the gas-phase values are not available for all vibrations. The difference between the solution and gas-phase values are only a few cm-' in most cases. The abbreviations ip and op denote vibrations in the ring plane and perpendicular to it, respectively.
differentiation of the gradients. The remaining minor uncertainties concern the very weak fundamentals 14 and 24. It appears that the latter is identical with the weak band at 1340 cm-I in the cluster around 1330 cm-I. not with the 1324 cm-' one. In order to conserve space, the calculated spectra of the other 10 molecules in the test set are given in the supplementary material. Conclusions The scaled quantum mechanical method, applied to density functional force fields with a moderate (6-31G*) basis set, successfully reproduces the vibrational spectra of 3 1 simple organic molecules with 11 scaling factors. The B3-LYP hybrid functional yields more accurate molecular geometries and, in combination with the SQM method, also better vibrational frequencies than the B-LYP functional. The total mean deviation of the B3-LYP frequenices after scaling is only 10.8 cm-' in the fingerprint region (500-2500 cm-l). In addition to multiple scaling, a single uniform scaling factor has been determined for both functionals. This is close to 1.0 for the B-LYP functional, and thus uniform scaling does not improve the B-LYP vibrational frequencies. As an example, the calculation of the IR spectrum of aniline led to a reassignment of some fundamentals. The calculated spectrum is in good agreement with the experimental data.
J. Phys. Chem., Vol. 99, No. 10, 1995 3099 Further work along these lines should concentrate on three subjects: (1) the inclusion of the most important anharmonic corrections in the calculated spectra (this would aid the reliable assignment of the weak fundamentals and greatly increase the predictive power of the method); (2) the enlargement of the present modest database; and (3) the inclusion of further structural motifs and heteroatoms, such as halogens. Acknowledgment. G.R. thanks the Deutsche Forschungsgemeinschaft for financial support. This work was also supported by the U.S. Air Force Office for Scientific Research under Grant No. F49620-94-1 and by the U.S. National Science Foundation (Grant No. CHE-9319929). We thank IBM Co. for the donation of IBM RS6000 workstations and Aldrich Chemical Co. for permission to reproduce a spectrum. Supplementary Material Available: Tables of calculated and experimental frequencies and infrared intensities (23 pages). Ordering information is given on any current masthead page. References and Notes (1) Pulay, P. In Applications of Electronic Structure Theory; Schaefer, H. F., 111, Ed.; Plenum: New York, 1977; p 153. (2) Fogarasi, G.; Pulay, P. In Molecular Spectra and Structure, Vol. 14; Dung, J. R. Ed.; Elsevier: Amsterdam, 1985; p 125. (3) Fogarasi, G.; Pulay, P. Annu. Rev. Phys. Chem. 1984, 35, 191. (4) Carsky, P.; Zahradnik, R. Chem. Rev. 1986, 86, 709. (5) Pulay, P. Mol. Phys. 1969, 17, 197. (6) Pople, J. A.; Krishnan, R.; Schlegel, H. B.; Binkley, J. S. Int. J . Quantum Chem., Quantum Chem. Symp. 1979, 13, 225. (7) Pulay, P.; Meyer, W. Mol. Phys. 1974, 27, 473. (8) Blom, C. E.; Altona, C. Mol. Phys. 1976, 31, 1377. (9) Pulay, P.; Fogarasi, G.; Pongor, G.; Boggs, J. E.; Vargha, A. J . Am. Chem. Soc. 1983, 105, 7073. (10) Pulay, P.; Fogarasi, G.; Pang, F.; Boggs, J. E. J . Am. Chem. Soc. 1979, 101, 2550. (11) Fogarasi, G.; Zhou, X.; Taylor, P. W.; Pulay, P. J. Am. Chem. Soc. 1992, 114, 8191. (12) Peretto, P. An Introduction to the Modeling of Neural Networks; Cambridge University Press: Cambridge, 1992; p 325. (13) Botschwina, P. Chem. Phys. Lett. 1974, 29, 98, 580. (14) Pulay, P.; Fogarasi, G.; Boggs, J. E. J . Chem. Phys. 1981, 74,3999. (15) Hehre, W. J.; Radom, L.; Schleyer, P. v. R.; Pople, J. A. Ab Initio Molecular Orbital Theory; Wiley: New York, 1986. (16) Yamaguchi, Y.; Schaefer, H. F., I11 J. Chem. Phys. 1980, 73, 2310. (17) Person, W. B.; Szczepaniak, K.; Szczeszniak, M.; Del Bene, J. In Recent Experimental and Computational Advances in Molecular Spectroscopy; NATO AS1 Series C, Vol. 406; Fausto, R., Ed.; Kluwer: Dordrecht, The Netherlands, 1993; p 141. (18) Zhou, X.; Pulay, P.; Hargitai, R.; Stirling, A,; Mink, J. Spectrochim. Acra 1993, 49A, 257. (19) Ahlrichs, R.; B&, M.; Ehrig, M.; Haser, M.; Horn, H.; Kolmel, C. TURBOMOLE; University of Karlsruhe: Karlsruhe, Germany, 1991, (20) Baker, J. Biosym Industries, San Diego, CA, private communication. (21) Parr, R. G.; Yang, W. Density Functional Theory of Atoms and Molecules; Oxford: New York, 1989. (22) Labanowski, J. K.; Andzelm, J. W., Eds. Densify Functional Methods in Chemistry; Springer: New York, 1991. (23) Vosko, S. H.; Wilk, L.; Nusair, M. Can. J . Phys. 1980, 58, 1200. (24) Becke, A. D. Phys. Rev. 1988, A38, 3098. (25) Lee, C.; Yang, W.; Parr, R. G. Phys. Rev. 1988, B41, 785. (26) Perdew, J. P.; Chevary, J. A,; Vosko, S. H.; Jackson, K. A,; Pederson, M. R.; Singh, D. J.; Fiolhais, C. Phys. Rev. 1992, B46, 6671. (27) Fan, L.; Verluis, L.; Ziegler, T.; Baerends, E. J.; Ravenek, W. Int. J . Quantum Chem., Quantum Chem. Symp. 1988, 22, 173. (28) Papai, I.; St-Amant, A.; Ushio, J.; Salahub, D. Int. J. Quantum Chem., Quantum Chem. Symp. 1990, 24, 29. (29) Fournier, R.; Andzelm, J.; Salahub, D. J. Chem. Phys. 1989, 90, 6371. (30) Fournier, R. J. Chem. Phys. 1990, 92, 5422. (31) Dunlap, B. I.; Andzelm, J. Phys. Rev. 1992, A45, 81. (32) Komornicki, A,; Fitzgerald, G. J . Chem. Phys. 1993, 98, 1398. (33) Johnson, B. G.; Frisch, M. J. Chem. Phys. Lett., in press. (34) Gill, P. M. W.; Johnson, B. G.; Pople, J. A.; Frisch, M. J. Int. J . Quantum Chem., Quantum Chem. Symp. 1992, 26, 319. (35) Handy, N. C.; Maslen, P. E.; Amos, R. D.; Andrews, J. S.; Murray, C. W.; Laming, G. J. Chem. Phys. Lett. 1992, 197, 506. (36) Berces, A.; Ziegler, T. J . Chem. Phys. 1993, 98, 4793.
3100 J. Phys. Chern., Vol. 99, No. 10, 1995 (37) Johnson, B. G.; Gill, P. M. W.; Pople, J. A. J . Chem. Phys. 1993, 98, 5612. (38) Handy, N. C.; Murray, C. W.; Amos, R. D. J. Phys. Chem. 1993, 97, 4392. (39) Florian, J.; Johnson, B. G. J. Phys. Chem. 1994, 98, 3681. (40) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Gill, P. M. W.; Johnson, B. G.; Wong, M. W.; Foresman, J. B.; Robb, M. A.; Head-Gordon, M.; Replogle, E. S.; Gomperts, R.; Andres, J. L.; Raghavachari, K.; Binkley, J. S.; Gonzales, C.; Martin, P. L.; Fox, D. J.; Defrees, D. J.; Baker, J.; Stewart, J. J. P.; Pople, J. A. GAUSSIAN G92/DFT; Gaussian Inc.: Pittsburgh, 1993. (41) Hariharan, P. C.; Pople, J. A. Theor. Chim. Acta 1973, 28, 213. (42) Becke, A. D. J . Chem. Phys. 1993, 98, 5648. (43) Wilson, J. E.; Decius, J. C.; Cross, P. C. Molecular Vibrations; Dover Publications: New York, 1955. (44) Zhou, Z. Ph.D. Dissertation, University of Arkansas, 1992. (45) Pongor, G. SCALE2, Eotvos, L. University, Budapest, 1978. (46) Thakur, S. N.; Goodman, L.; Ozkabak, A. G. J. Chem. Phys. 1986, 484, 6642. (47) Duncan, J. L.; McKean, D. C.; Mallinson, P. D. J . Mol. Spectrosc. 1973, 45, 221. (48) Fogarasi, G.; Pulay, P. Acta Chim. Acad. Sci. Hung. 1981, 108,
55. (49) (a) Serralach, A.; Meyer, R.; Giinthard, Hs. H. J . Mol. Spectrosc. 1974,52,94. (b) Blom, C. E.; Otto, L. P.; Altona, C. Mol. Phys. 1976,32, 1137. (50) Blom, C. E.; Altona, C.; Oskam, A. Mol. Phys. 1977, 34, 557. (51) (a) Pulay, P.; Torok, F. J. Mol. Struct. 1975, 29, 239. (b) Hamada, Y.; Tanaka, N.; Sugawara, Y.; Hirakawa, A. Y.; Tsuboi, M.; Kato, S.; Morokuma, K. J. Mol. Spectrosc. 1982, 96, 313. (52) Blom, C. E.; Altona, C. Mol. Phys. 1976, 31, 1377. (53) Blom, C. E.; Altona, C. Mol. Phys. 1977, 33, 875. (54) Wiberg, K. B. J. Mol. Struct. 1990, 244, 61. ( 5 5 ) Tanaka, N.; Hamada, Y.; Sugawara, Y.; Tsuboi, M.; Kato, S.; Morokuma, K. J . Mol. Spectrosc. 1983, 99, 245. (56) Pople, J. A.; Schlegel, H. B.; Krishnan, R.; Defrees, D. J.; Binkley, J. S.; Frisch, M. J.; Whiteside, R. A,; Hout, R. F.; Hehre, W. J. Int. J . Quantum Chem., Quantum Chem. Symp. 1981, 15, 269.
Rauhut and Pulay (57) Dutler, R.; Rauk, A.; Shaw, R. A. J . Phys. Chem. 1990, 94, 118. (58) Givan, A.; Loewenschuss, A. J. Mol. Struct. 1983, 98, 231. (59) Xie, Y.; Fan, K.; Boggs, J. E. Mol. Phys. 1986, 58, 401. (60) Scott, D. W. J. Mol. Spectrosc. 1971, 37, 77. (61) Redington, R. L. J. Mol. Spectrosc. 1977, 65, 171. (62) Brand, J. C. D.; Williamson, D. G. Discuss. Faraday SOC. 1963, 35, 184. (63) Pulay, P.; Lee, J.-G.; Boggs, J. E. J. Chem. Phys. 1983, 79, 3382. (64) Evans, J. C. Spectrochim. Acta 1%0, 16, 428. (65) Shaw, R. A.; Wieser, H.; Dutler, R.; Rauk, A. J. Am. Chem. SOC. 1990, 112, 5401. (66) Banhegyi, G.;Pulay, P.; Forgarasi, G. Spectrochirn. Acta 1983,39A, 761. (67) Liu, R.; Zhou, X.; Pulay, P. J . Phys. Chem. 1992, 96, 4225. (68) Kuwae, A.; Machida, K. Spectrochim. Acta 1979, 35A, 27. (69) Bist, H. D.; Brand, J. C. D.; Williams, D. R. J. Mol. Spectrosc. 1967, 24, 402. (70) (a) A d h e k , P.; Volka, K.; Ksandr, Z.; Stibor, I. J. Mol. Spectrosc. 1973, 47, 252. (b) Banki, J.; Billes, F.; Gal, M.; Grofcsik, A.; Jalsovszky, G.; Sztraka, L. J. Mol. Struct. 1986, 142, 351. (71) Green, J. H. S.; Harrison, D. J. Spectrochim.Acta 1976,32A, 1265. (72) Green, J. H. S.; Harrison, D. J. Spectrochim.Acta 1976,32A, 1279. (73) (a) Susi, H.; Ard, J. S. Spectrochim. Acta 1971, 27A, 1549. (b) Szczesniak, M.; Nowak, M. J.; Rostkowska, H.; Szczepaniak, K.; Person, W. B.; Shugar, D. J . Am. Chem. SOC. 1983, 105, 5969. (c) Cshszir, P.; Harshnyi, G.; Boggs, J. E. Int. J . Quant. Chem. 1988, 33, 1. (74) Dellepiane, G.; Overend, J. Spectrochim. Acta 1966, 22, 593. (75) Niu, Z.; Dunn, K. M.; Boggs, J. E. Mol. Phys. 1985, 55, 421. (76) Pouchert, C. J. The Aldrich Librav of FT-IR Spectra, 1st ed., Vol. 3; Aldrich: Milwaukee, WI, 1989. (77) Quack, M.; Stockburger, M. J . Mol. Spectrosc. 1972, 43, 87. (78) Larsen, N. W.; Hansen, E. L.; Nicolaisen, F. M. Chem. Phys. Lett. 1976, 43, 584.
JP9428808
