J. Phys. Chem. 1995,99, 15837-15845
15837
Evaluation of the Density Functional Approximation on the Computation of Hydrogen Bond Interactions Juan J. Novoat and Carlos Sosa* Cray Research, Inc., 655 E Lone Oak Drive, Eagan, Minnesota 55121 Received: April 11, 1995; In Final Form: July 24, 1995@
The intermolecular interaction energies for hydrogen-bonded complexes have been investigated using density functional theory (DFT). Contributions to the exchange energy were estimated with the Slater-Dirac (S), Becke (B), and Xu functionals. Exchange-correlation effects were estimated with the exchange functional of Becke and the correlation energy functional of Vosko, Wilk, and Nusair (B-VWN), the exchange functional of Becke and the correlation energy functional of Lee, Yang, and Parr (B-LYP), and the Slater-Dirac exchange and Lee Yang and Parr (S-LYP) and Vosko, Wilk, and Nusair correlation functionals (S-VWN). In addition, the hybrid method B3-LYP has been tested. The following complexes were fully optimized with the above functionals: HF-HF, H20-Hz0, C z H r H 2 0 , CK-HzO, and N H 3 - N H 3 . All the optimizations were carried out with Pople's split valence 6-31++G(2d,2p), basis set. For each of these complexes, interaction energies have been calculated using the correlation consistent polarized valence double-zeta (cc-pVDZ), Pople' s split valence (6-3 1++G(2d,2p) and augmented correlation consistent polarized valence triple-zeta (aug-cc-pVTZ) basis sets at the MPn ( n = 2, 3 , and 4) and DFT levels. The importance of the basis set superposition error (BSSE) is addressed for each density functional approximation and basis sets. DFT calculations with small basis sets (cc-pVDZ) are more sensitive to BSSE than MPn methods. The interaction energies of the HFHF, H20-Hz0, C2H2-H20, C&-H20 and NH3-NH3 hydrogen-bonded complexes estimated with B3-LYP/ 6-31++G(2d,2p) are -4.68, -4.82, -2.65, -0.12, and -2.67 kcal/mol, respectively.
Introduction The importance of hydrogen-bonded dimers has become apparent in the last decade. Hydrogen bond interactions are known to play an important role in determining the structure of molecular crystals' and biological systems.2 This has prompted experimentalists3as well as theoritician~~ to focus their attention on these small systems. A wealth of high-level ab initio molecular orbital calculations have been reported for hydrogen-bonded dimerse4 Some of these systems have proven to be highly sensitive to the level of theory and basis sets? They provide rigorous test cases for new methodologies. Density functional theory (DFT)5offers a promising tool that may be applied to large systems. It includes correlation effects in a form in which the cost is only p,where N is the number of basis functions. In practice, the scaling can be This is in contrast to conventional ab initio methods such as MP2 or MP4 that scale as Ns or N ,respectively.6b Also, because of recent developments of analytic gradient techniques? progress in functionals for gradient correction^^^^ to the local spin density (LSD) approximationlo and the availability of versatile software,6a,11the density functional is becoming a popular tool in chemistry. Calculation of thermodynamic data for polyatomic molecules' 1d.12.13 shows that the gradient-corrected functionals provide energetics typically better than the Hartree-Fock (HF) method. The gradient-correctedresults are closer to correlated methods (e.g., MP2) or better.11d In recent papers,I3 we presented systematic studies using density functional methods for bond dissociation, rotational barriers and transition structures. For these systems (single bond breaking), we found that the density functional can be used to study chemical reactivity. In t Permanent address: Department de Quimica Fisica, Facultat de Quimica, Universitat de Barcelona, Av. Diagonal 647, 08028 Barcelona, Spain. Abstract published in Advance ACS Abstracts, October 1, 1995. @
0022-3654/95/2099-15837$09.00/0
particular, gradient corrections are important in the proper description of the dissociation process. In the present study we extend our studies to include the performance of several DFT functionals for various hydrogenbonded complexes selected to cover interaction energies between HF, H20, and NH3 dimers, CZHZ-H20, and CH4-H20. The goal of this work has been (a) to establish how well the local and gradient-corrected functionals reproduce the MP2 or MP4 interaction energy for a given basis set, (b) to determine the basis set dependence of the DFT functionals and how it compares to MP2, (c) to compute the size of the basis set superpositionerror (BSSE)I4for these functionals compared to the MP2 values and how well it is corrected using the counterpoise method (CP),I5and (d) to determine how well the optimum energy structures compute using DFT compare to the predicted structures using MP2 . The complexes selected here are the HF-HF and Hz0-H20, as examples of moderately strong complexes for which there are estimates of the interaction respectively, corrected energy (4.56 f 0.2916 and 5.4 kcaVm01,~~ for the change in zero-point energies), and the NH3-NH3, CzHz-HzO, and C&-H20 complexes, which are good examples of moderate to weak hydrogen-bonded complexes. In each case, the interaction energies have been computed with medium to large basis sets. Also, basis set superposition errors (BSSE) for the interaction energy were estimated using the CP method.I5 There is still a widespread feeling that the BSSE values obtained using the CP method are only approximate and that this method overestimates the size of the BSSE.IS Recent theoretical studiesIg on van der Waals systems and hydrogenbonded complexes20,21suggest that it is possible to predict interaction energies within experimental uncertainty. They used the CP method in conjunction with a large basis set, i.e., a basis set which is flexible enough to provide a good description of the electrons in both the atomic core and the long-range region of the atoms. These tests have shown that by use of adequate 0 1995 American Chemical Society
15838 J. Phys. Chem., Vol. 99, No. 43, 1995
Novoa and Sosa
basis sets, one can remove virtually all the BSSE and obtain corrected interaction energy (AEcp) values very close to the available experimental data,'9-21 while, as was already noticed before,I8 small basis sets tend to overcorrectthe BSSE.'4f-gThree basis sets have been employed in this study to perform a similar test using several DFT functionals and to explore if the same conclusions are valid. Finally, all the geometries were fully optimized for the selected hydrogen-bonded complexes using several DFT functionals. The results were compared to Hartree-Fock and MP2 results and experimental data where available.
6-3 1++G(2d,2p) level. We have selected for the NH3-NH3 complex its optimum HF/6-3 1++G(2d,2p) geometry, while for the CI&-HzO complex is the geometry of the RCHs-HzO groups as found in the orthoamide crystal,38 with the R group replaced by a hydrogen located at a C-H distance of 1.09 8,. This simplifies the analysis since it eliminates differences in the interaction energies that would result from slightly different geometries used with the different methods. We have computed the interaction energy (0 and the counterpoise-correctedinteraction energy (hEcp)39 for the HFHF complex with the cc-pVDZ, 6-31++G(2d,2p), and aug-ccpVTZ basis sets at the MP2, MP3, and MP4 methods. For all Computational Method the other complexes, 6-3 1++G(2d,2p) was used throughout these calculations. The MPn methods are size consistent, and First principles calculations were carried out using the since there is no near degeneracy in any of the systems studied Gaussian92AIFT program22on a Cray C90. The Kohn-Sham here, no multireference formulation seems to be required. We equations were solved as described by Pople and ~o-workers.6~.~~ should expect, according to previous calculations for some The Coulomb and exact exchange integrals were evaluated hydrogen-bonded system^,'^^^^-^^ that the MP2 method gives analytically. All the integrals for the approximate exchange results very similar to those of the MP4 method and close to and correlation contributions to the energy were canied out over the CCSD(T) values.40 a numerical grid.24 In the present study we have tested within the LSD apResults and Discussion proximation the following combinations: Slater-Dirac ( S ) exchange25 and Vosko, Wilk, and Nusair (VWN) correlation The geometries for these complexes were predicted using functionalz6(S-VWN); Slater-Dirac exchange (S-null) and the density functional methods. All the optimized geometries were compared to conventionalab initio methods. Hartree-Fock and Xuz7 exchange functional only (p4I3with the empirical coefficient of 0.7). Two other LSD approximationswere also tested. correlated methods have been extensively employed in many Gradient corrections were introduced using Becke28exchange chemical ~tudies.4~ They provide a very reliable way to calibrate density functional methods applied to chemical reactions and only (B-null) or using the Lee, Yang, and Parr29correlation (BLYP). The LYP correlation is the density functional energy especially to weakly bound complexes. Interaction energies formulation of the Colle-Salvetti correlation energy.30 We have with and without counterpoise corrections were compared with also tested the following combined functionals: B-VWN and Hartree-Fock and MPn ab initio results as well. S-LYP. Becke's three-parameter hybrid functional with gradient The calculated geometries for this series of hydrogen bond corrections provided by the LYP functional (B3-LYP)3' has also complexes and N H 3 are summarized in Tables 1-6. Systematic been tested in this study. studies of the basis sets and correlated effects of the HF-HF, All the calculations for the interaction energies have been HzO-HzO, and NH3-NH3 complexes have been reported in performed with the following: the correlation consistent polarthe literature.'& In this section we only summarize our results ized valence double-zeta ((9s4pld)/[3s2pld] for C, N, 0, and using DFT. Figures 1 and 2 illustrate all the parameters used F and (4slp)/[2slp] for H) denoted as C C - ~ V D and Z ~ ~b1 in throughout these calculations. this work; Pople's split-valence 6-3 l++G(2d,2p) (these basis Geometries. A. HF-HF Dimer. Table 1 summarizes all sets include diffuse functions on heavy and hydrogen atoms, the optimized geometries for the HF-HF dimer. All our results where the polarization function space consists of two d's on obtained using Hartree-Fock and MP2 methods are in good heavy atoms and two p's on hydrogen and is denoted bII in agreement with results of previous calculations at similar levels this and augmented correlation consistent polarized of theory. In the case of the HF-HF dimer (see Figure la), valence triple-zeta (( 1ls6p3d2f)/[%4p3d2fl for C, N, 0, and F DFT calculations tend to have a marked basis set effect. This and (6s3p2d)/[4s3p2d] denoted as aug-cc-pVTZ (denoted bIII is particularly true for the r(FF) distance. All the functionals in this work)32basis sets. tested with cc-pVDZ basis sets tend to underestimate the r(FF) Hartree-Fock and MPn (n = 2,3,4) orbital calculationswere distance by about 0.2 8, to almost 0.4 8,. On the other hand, performed with the Gaussian92AIFT program.22 Full geometry the r(H,F; x = a and c) bond lengths predicted with the different optimizations at the HF, MP2, and DFT levels were carried out functionalshave an error range 0.001-0.01 A. A more dramatic with the 6-3 l++G(2d,2p) basis sets only. Electron correlation basis set effect may be seen by looking at the predicted angles; single point energies were estimated by using Mgller-Plesset small basis sets (bI) clearly underestimate the f3fand e,, angles. perturbation theory34 up to fourth order, including all single, The errors are on the order of 60". This difference may be double, triple, and quadruple excitations (MP4SDTQ, frozen partially attributed to the need for using better basis sets to core). properly describe the density on fluorine because of the lack of diffuse functions in the cc-pVDZ basis sets. All the interaction energy calculations were carried out using Similar effects may be observed at the MP2 level as well. experimental geometries where available. For the HzO-H~O However, they are not as dramatic as those in the DFT case. In complex, instead of using the experimental geometry,35 we general, very good agreement has been found between the selected Feller's optimized MF'2/6-311++G(2d,2p) geometry?O 6-3 1++G(2d,2p) and the aug-cc-pVTZ results for all the forcing the 0-H-0 bond to be collinear. For the HF-HF complex, the dimer was located at its experimental e ~ m e t r y ~ ~ methods. with the H-F distance of each monomer at 0.917 . For the In addition to the basis set effects, it is important to note the effect of the different functionals on the structures of this dimer. CZHZ-HZOcomplex, the CzVconformation was selected (the effective conformation observed e~perimentally)~~ and posiIn general, it may be seen (S-null B-null) that adding gradient tioned the monomers at the experimentalH-0 distance of 2.229 corrections only to the exchange approximation tends to elongate The geometry of the fragments was optimized at the MP2/ the r(FF) distance by as much as 0.4 A. On the other hand,
x
-
Evaluation of the Density Functional Approximation
J. Phys. Chem., Vol. 99,No. 43, 1995 15839
TABLE 1: Optimized Geometries for the HF-HF CompleP methodb
r(FF)
bI bII bIII
2.756 2.800 2.819
bI bII bIII
2.682 2.725 2.733
bI bII bIII
2.392 2.618 2.619
bI bII bIII
2.661 3.023 3.015
bI bII bIII
2.367 2.575 2.574
bI bII bIII
2.351 2.550 2.553
bI bII bIII
2.617 2.877 2.871
bI bIII
2.304 2.504 2.514
bI bII bIII
2.517 2.764 2.755
bI bII bIII
2.523 2.729 2.726
Of
r(HaF)C r(RFY HF 0.905 0.906 0.900 0.902 0.901 0.903 MP2 0.924 0.925 0.923 0.926 0.925 0.928 S-null 0.967 0.967 0.950 0.960 0.951 0.961 B-null 0.954 0.954 0.944 0.946 0.945 0.948 Xa 0.955 0.955 0.939 0.949 0.939 0.950
Od"
113.6 115.9 119.7
9.0 7.4 6.7
97.9 110.1 110.3
12.8 6.2 6.4
42.5 103.4 103.0
42.1 7.6 7.4
46.3 108.7 110.6
47.4 9.4 6.9
TABLE 4: Optimized Geometries for the CJ&-HzO ComplelP
42.5 102.9 102.4
42.6 8.0 7.9
HF MP2 S-null B-null
42.8 101.9 104.5
41.8 8.2 6.5
46.4 110.3 110.7
47.5 7.1 6.8
41.5 97.5 103.6
41.1 10.1 6.6
44.8 110.5 110.7
45.4 5.7 5.5
48.3 112.4 112.0
45.1 5.6 6.1
S-VWN
bII
0.954 0.954 0.937 0.948 0.937 0.949 B-VWN 0.941 0.940 0.932 0.935 0.932 0.935 s-LYP 0.964 0.964 0.943 0.956 0.944 0.957 B-LYP 0.948 0.948 0.937 0.943 0.938 0.943 B3-LYP 0.935 0.935 0.926 0.931 0.927 0.932
Bond lengths in angstroms; angles in degrees. bI, bII, and bIII correspond to cc-pVDZ, 6-31++G(2d,2p), and aug-cc-pVTZ,respectively. See Figure la for parameter definition.
TABLE 2: Optimized Geometries for the HzO-HZO Compleg methodb HF MP2 S-null B-null
r(00) Or ed r(OHJ eh 3.020 48.9 3.6 0.945 0.940 0.941 106.2 2.902 57.5 5.0 0.968 0.959 0.961 104.5 2.795 65.7 6.2 1.003 0.985 0.987 105.1 3.229 59.3 5.3 0.988 0.982 0.984 104.1 5.5 0.991 0.972 0.975 105.4 2.746 63.6 x, S-VWN 2.715 65.6 5.9 0.989 0.969 0.972 105.4 B-VWN 3.063 59.4 5.1 0.974 0.967 0.968 104.5 S-LYP 2.668 67.5 7.3 0.998 0.975 0.978 106.0 B-LYP 2.943 59.8 5.1 0.982 0.972 0.974 104.7 B3-LYP 2.915 57.7 4.8 0.971 0.962 0.963 105.2 expt 2.946' 57 f 1W 2 f 10' a Bond lengths in angstroms; angles in degrees. All the optimizations were carried out with the 6-31++G(2d,2p) basis sets. See ref 14h.
-
TABLE 3: Optimized Geometries for the C ~ H Z - H ~ O ComDleP methodb r(C0) Of 0, r(CC) r(CH,) r(CHb) r(OH) HF 3.399 7.8 0.2 1.187 1.060 1.056 0.941 MP2 3.282 33.1 0.7 1.219 1.069 1.064 0.960 S-null 3.108 9.0 -1.1 1.104 1.227 1.092 0.986 B-null 3.759 7.8 0.0 1.083 1.227 1.081 0.983 1.080 0.973 0.0 1.092 1.214 x, 3.056 7.9 S-VWN 3.011 8.7 -0.8 1.089 1.211 1.075 0.971 B-VWN 3.527 7.9 0.1 1.068 1.212 1.065 0.968 S-LYP 2.939 8.0 -0.1 1.099 1.215 1.082 0.976 B-LYP 3.330 8.0 0.0 1.076 1.216 1.071 0.973 0.0 1.070 1.206 1.065 0.963 B3-LYP 3.285 7.9
the effect of introducing correlation (S-null S-VWN) has the opposite effect. The r(FF) distance is decreased by about 0.06 A. Gradient-corrected functionals (S-LYP)shorten the r(FF) distance even further. These approximations (S-null, B-null, S-VWN, and S-LYP) tend to either overestimate or underesti-
a Bond lengths in angstroms; angles in degrees. All the optimizations were carried out with the 6-31++G(2d,2p) basis sets.
4.113 3.765 3.506 5.974 3.414 Xa S-VWN 3.352 B-VWN 5.911 s-LYP 3.234 B-LYP 3.936 B3-LYP 3.891
28.7 27.7 32.6 29.6 41.7 50.7 29.6 50.6 37.5 37.3
3.1 4.5 1.1 3.9 1.5 1.1 3.0 1.0 6.3 6.6
1.083 1.087 1.120 1.111 1.106 1.102 1.092 1.110 1.098 1.091
1.084 1.087 1.118 1.112 1.104 1.099 1.092 1.105 1.098 1.092
1.084 1.OS8 1.119 1.112 1.105 1.099 1.093 1.106 1.098 1.092
0.941 0.960 0.986 0.983 0.974 0.971 0.968 0.977 0.973 0.963
a Bond lengths in angstroms; angles in degrees. All the optimizations were carried out with the 6-31++G(2d,2p) basis sets.
TABLE 5: Theoretical and Experimental Geometries for NH-P r e P CC-PVDZ MP2 1.025 104.2 1.872 QCISD 1.027 104.2 1.784 B3-LYP 1.027 104.6 1.733 6-31++G(2d,2p) MP2 10.15 106.7 1.777 QCISD 1.016 106.7 1.761 B3-LYP 1.018 106.9 1.745 cc-aug-pVTZ MP2 1.015 106.6 1.616 B3-LYP 1.016 107.0 1.599 exptb 1.012 106.7 1.47 Bond lengths in angstroms; angles in degrees; dipole moments in Debye. All the optimizations were canied out with the 6-31++G(2d,2p) basis sets. * See refs 50 and 51. mate the M E predicted r(FF) distance. This behavior may be observed in the family of molecules reported by Pople and coworkersasa Their results suggest that adding gradient corrections to only an exchange functional is to elongate bond lengths between two heavy atoms. The opposite effect may be observed when one of the atoms is hydrogen or lithium. Again, opposite effects are observed in both cases when only gradient corrections are added to the correlation functional. Our predicted r(FF) distance with the B3-LYP approximation using aug-cc-pVTZ basis set is 2.726 A. This value is in excellent agreement with MP2 results. It is also within the experimental range for the reported Re value of 2.72 A.36 Similarly, as in previous DFT studies,44d the bond lengths predicted with the local and non-local approximations tend to be shorter than ab initio results carried out at the CISD and ACCD leve1.44bb This difference may be attributed to use of
15840 J. Phys. Chem., Vol. 99, No, 43, 199s
Novoa and Sosa
TABLE 6: Optimized Geometries for the NHrNH? Complex" methodb r(") ed r(NHa) 4") HF MP2 S-null
3.485 3.240 3.093 3.757 3.034 2.985 3.494 2.919 3.325 3.295
B-null X, S-VWN
B-VWN s-LYP B-LYP B3-LYP (I
11.7 19.4 14.4 15.9 15.2 80.3 15.7 15.0 8.6 10.2
1.002 1.016 1.052 1.039 1.039 1.024 1.022 1.045 1.030 1.021
"c)
1.Ooo 1.014 1.041 1.037 1.027 1.037 1.019 1.029 1.025 1.016
1.Ooo
1.013 1.041 1.037 1.027 1.023 1.019 1.029 1.025 1.016
r(NHd)
eh
6,
Of
1.Ooo 1.011 1.041 1.037 1.028 1.023 1.019 I .030 1.025 1.016
107.9 107.0 107.5 106.1 107.8 107.3 106.6 108.2 106.9 107.4
121.9 125.7 123.1 124.1 123.7 115.4 124.0 123.8 120.0 121.0
85.3 73.8 85.3 83.9 80.6 15.1 82.8 79.8 95.0 90.7
Bond lengths in angstroms; angles in degrees. All the optimizations were carried out with the 6-31++G(2d,2p) basis sets.
geometries and the MP2 and experimental g e ~ m e t r i e s . 'Their ~~ basis sets are similar to the ones employed throughout this work. In order to systematically study the effect of different functionals on the geometries for these complexes (HzO-H~O, CZHZ-H~O, and CH4-H20), we have selected the 6-31++G(2d,2p) basis sets, as previously pointed out. Examination of the results in Tables 2-4 indicate considerable changes in the 4x0)(X = C or 0) distance at the DIT level as different functionals are compared. The trend observed for these water complexes is similar to the trends already pointed out for the HF-HF dimer. Adding gradient corrections to the exchange functional (S-null B-null) tends to elongate the 4 x 0 ) distance. In all the cases B-null overestimatesthe 4 x 0 ) distance when compared to MP2 values. In the case of CH4-H20, B-null overestimates the r(C0) distance by as much as 2.2 A. Fortuitously, S-null predicts r ( X 0 ) distances that are in close agreement with MP2 results. Again, the effect of adding correlation corrections within the electron-gas approximation (local approximation) S-VWN) is to predict 4 x 0 ) distances that are too (S-null short in comparison to MP2 and experimental results where available. The effect of including gradient corrections to the correlation functional is negligible and tends to shorten the 4x0)distance even further. On the other hand, adding gradient corrections on an equal footing (B-LYP and B3-LYP) tend to predict r(X0) distances that are in very good agreement with predicted MP2 distances. C. NH3-NH3. Tables 5 and 6 give parameters of the optimized geometries for NH3 and the NH3-NH3 dimer. The prediction of the structure of the ammonia dimer has proven to be difficult and particularly sensitive to the BSSE.47-49In order to calibrate the B3-LYF' approximationfor the NH3-NH3 dimer, Table 5 summarizes optimized and experimental geometries for the NH3 molecule. Pople and have reported optimized geometries for NH3 using most of the functionals employed in this study. Our best estimated Hartree-Fock-based theoretical parameters for the NH3 geometry are 1.021 4Z 0.007 8, and 105.5" 4Z 1" for the r(NH) and HNH angle, respectively. These results are within 0.01 8, and 1-2" from experimental values for the NH distance and HNH angle, r e s p e c t i ~ e l y . ~ ~ ~ ~ ' All the parameters predicted using the B3-LYP approximation with the 6-31++G(2d,2p) and aug-cc-pVTZ basis sets are within 0.004 8, and 0.4" when compared to experimental results.50 Similar results were found with our best estimated Hartree-Fock-based prediction. The largest discrepancy appears to be between the predicted dipole moment and the reported experimental ~ a l u e . ~Calculated ' dipole moments with cc-pVDZ and 6-3 1++G(2d,2p) overestimate the experimental dipole moment by about 0.2-0.3 D. Larger basis sets (augcc-pTVZ) predict a dipole moment of 0.1 D when compared to experimental results. Three structures have mainly been suggested the equilibrium geometry of this dimer. The first structure corresponds to the
-
/ H C
Figure 1. Hydrogen bond complexes: (a) HF-HF (b) NH3-NH3.
-
Ha.
Ha,
,H
Figure 2. Hydrogen bond complexes. (a) H20-H20 (b) HCCH-H20; (c) Cb-HzO.
different basis sets. We have optimized the HF-HF dimer at the QCISD/6-3l++G(2d,2p) level. The predicted r(FF) distance is 2.736 A. This value is in very good agreement with the predicted B3-LYP distance. B. H20, CzHz, and CH4 Water Complexes. Tables 2-6 summarize all the optimized geometries (see Figures 1 and 2) for these three complexes. All these structures were optimized using the 6-31++G(2d,2p) basis sets. This basis set was selected based on previous experience with the HF-HF dimer. None of the X-H-Y hydrogen bonds are collinear (X and Y are the two heavy atoms in each of the XH,-YH, complexes, where X and Y can each be 0 or C), although the deviation from colinearity is negligible for the C2H2-H20 complex (see angle &, Figure 2a). Recently, Feller20 and Kim and Jordanik have extensively studied the H20-H20 dimer using Hartree-Fock and postHartree-Fock methods as well as DlT. Kim and Jordan'& have reported good agreement between the predicted B3-LYP
J, Phys. Chem., Vol. 99, No. 43, 1995 15841
Evaluation of the Density Functional Approximation
TABLE 7: Total Energies Eb and Interaction Energies AEc and A E d Calculated for the Minimum Energy Structure of €IF-€IF
CC-pVDZ method" HF
MP2 MP3 MP4SDTQ S-null
B-null
x, S-VWN B-VWN s-LYP B-LYP B3-LYP expt
E
AE
-200.048 966 -5.84 -200.458 450 -6.80 -200.462 833 -6.45 -200.472 046 -6.68 -198.175 288 -10.87 -200.125 850 -5.85 -199.122 999 -10.89 -199.983 194 -11.12 -201.935 868 -6.24 -198.903 431 -13.37 -200.854 571 -8.37 -200.886 359 -8.10 -4.56 f 0.2gd
6-31+fG(2d32p)
A& -4.08 -3.93 -3.86 -3.78 -6.77 -1.81 -7.26 -7.81 -2.98 -9.46 -4.52 -4.92
E
-200.061 -200.514 -200.513 -200.531 -198.222 -200.171 -199.165 -200.021 -201.973 -198.949 -200.898 -200.921
221 111 931 846 714 557 133 581 080 398 887 386
aug-cc-pVTZ
AE
AEc,
E
AE
AEm
-4.07 -4.62 -4.69 -4.68 -6.59 - 1.45 -7.06 -7.67 -2.70 -9.33 -4.20 -4.68
-3.79 -3.88 -3.98 -3.87 -6.36 - 1.20 -6.84 -7.44 -2.45 -9.10 -3.94 -4.44
-200.128 447 -200.690 216 -200.688 131 -200.71 1 373 -198.291 517 -200.237 491 -199.234 419 -200.090 649 -202.039 419 -199.018 931 -200.965 447 -200.987 246
-3.81 -4.61 -4.76 -4.69 -6.12 -1.03 -6.60 -7.20 -2.29 -8.82 -3.72 -4.26
-3.67 -4.05 -4.19 -4.09 -6.02 -0.94 -6.49 -7.09 -2.18 -8.72 -3.63 -4.16
All calculations were canied out with the 6-31++G(2d,2p) basis sets. Total energy in au. Interaction energies in kcdmol. See ref 16.
TABLE 8: Total Energies Eb and Interaction Energies AEc and AE,$ Calculated for the Minimum Energy Structure of H20-H20 cc-pV DZ 6-31++G(2d,2p) aug-cc-pVTZ method" HF
MP2 MP3 MP4SDTQ S-null B-null
x,
s-VWN B-VWN
s-LYP B-LYP B3-LYP expt
E
AE
-152,062783 -5.66 -152.472 35i -7.30 -152.484 805 -6.83 -152.496 224 -7.09 -152.402 157 -11.08 -152.125 177 -5.47 -151.212 905 -11.13 -152.123 286 -11.43 -153.848 704 -5.98 -151.087 414 -13.67 -152.810 677 -8.06 -152.854 746 -7.89 -5.4 f 0.7d
-3.59 -3.85 -3.72 -3.65 -6.55 -1.05 -7.10 -7.83 -2.46 -9.43 -3.92 -4.42
E
AE
AEcn
E
AE
AEcD
-152.078969 -152.535 456 -152.543 296 -152.561 166 -150.446099 -152.169 124 -151.252 565 -152.159 220 -153.885 261 -151.129 564 -152.852 870 -152.889 369
-3.99 -5.27 -5.14 -5.30
-3.62 -4.44 -4.37 -4.38 -6.66 -1.12 -7.20 -7.93 -2.51 -9.53 -3.95 -4.48
-152.126973 -152.666 634 -152.671 984 -152.696 009 -150.500 275 -152.220 025 -151.307 134 -152.213 605 -153.936 499 -151.184454 -151.904 427 -152.939 345
-3.55 -5.16 -5.14 -5.20 -6.70 -1.24 -7.21 -7.91 -2.61 -9.51 -4.03 -4.52
-3.45 -4.57 -4.58 -4.61 -6.60 -1.12 -7.09 -7.78 -2.47 -9.40 -3.90 -4.26
-7.05
-1.48 -7.57 -8.27 -2.85 -9.91 -4.31 -4.82
All calculations were carried out with the 6-31++G(2d,2p) basis sets. Total energy in au. Interaction energies in kcdmol. See refs 17. classical structure as defined by Nelson et al.48b In this approach one of the N-H bonds is collinear with the C3 axis of the other NH3 molecule, one with Cssymmetry,Ik and the other structure corresponds to a cyclic structure with CZh ~ y m m e t r y . 4 Re~~ cently, Tao and Klempereflgchave investigated the equilibrium geometry of the NH3-NH3 dimer using 6-311+G(3d,2p) and [7s5p3d,4slp] extended with bond functions at the HartreeFock and MPn levels. They concluded that when the basis sets are enlarged with bond functions, the C Z structure ~ becomes more stable. On the other hand, basis sets without bond functions predict the C, structure to be more stable. Recently, Muguet and Robinson49have reported a new method to correct BSSE. In their study they carry out an extensive analysis of the BSSE for the NH3-NH3 dimer and presented a section discussing Tao and K l e m p e r e r ' ~work. ~ ~ ~ The objective of the present study is mainly to calibrate DFT against MPn results for this type of interaction. All the DFT geometries were obtained by reoptimizing Hartree-Fock geometries reported by Frisch et a l . I k Figure 1 illustrates all the optimized parameters for the NH3NH3 dimer. All the calculations carried out with DFT approximations predict a nonlinear hydrogen bond. As in the two previous cases, r(NN) distances predicted with the B-null approximation tend to overestimate the MP2 results by about 0.4-0.5 A. The S-VWN underestimates this distance by 0.282 A when compared to the MP2 distance. B-LYP and B3-LYP approximationspredict r(NN) values that are in close agreement with MP2. In particular, the B3-LYP predicted distance is within 0.03 A of MP2. All the r ( M ) (x = a, b, c, and d) parameters computed with DFT approximations are in good
agreement with MP2 r(NHx) values. The largest discrepancy is for r(NHa) predicted with the S-null approximation. The deviation for this parameter from the MP2 prediction is 0.035 A. Similar trends may be observed for all the angles. Energetics. The total energy (E), and the interaction energies with and without counterpoise corrections (AI3 and hEcp)for HF-HF, H20-Hz0, C ~ H Z - H ~ OC&-H20, , and NH3-NH3 complexes are collected in Tables 7-11. These series of calculations indicate a strong basis set effect for all the methods employed in this study. The basis set dependence of the interaction energy (AE), shown in Figure 3 for the HF dimer, follows the same pattem in all the complexes. The 6-31++G(2d,2p) basis set provides values of hE that are in excellent agreement with the aug-ccpVTZ values. On the other hand, this is not the case for results obtained with smaller basis sets (cc-pVDZ). Although Figure 3 shows similar pattems among all the curves, the interaction energy depends on the density functional approximation. The DFT energy difference between the curves generated using 6-31++G(2d,2p) and aug-cc-pVTZ and the cc-pVDZ basis sets is about twice the difference in the Hartree-Fock or MPn values. Part of this difference may be attributed to the large BSSE when using cc-pVDZ basis sets. DFT approximations tend to be more susceptible to the basis sets truncation than conventional methods for the systems presented in this study (see Figure 4). Although the BSSE of the DFT functionals is almost constant for the 6-3 l++G(2d,2p) or aug-cc-pVTZ basis sets, it oscillates with different functionals for the cc-pVDZ basis set and is larger in the DFT computations. The shape of the BSSE curve shown in Figure 4 follows the same pattem as in
15842 J. Phys. Chem., Vol. 99,No. 43, 1995
Novoa and Sosa
TABLE 9: Total Energies Eb and Interaction Energies Mcand MCpC Calculated for the Minimum Energy Structure of CIHI-H~O
CC-pVDZ method"
HF MP2 MP3 MP4SDTQ S-null B-null XLI S-VWN B-VWN s-LYP B-LYP B3-LYP @
E -152.858 -153.324 -153.341 -153.364 -150.916 -152.922 -151.812 -152.927 -154.936 -151.703 -153.709 -153.765
452 311 570 062 819 302 701 648 178 932 702 272
6-3 1++G(2d,2p)
aug-cc-pVTZ
AE
AEcp
E
AE
AEcp
E
AE
AEq
-4.16 -5.12 -4.81 -5.00 -7.77 -4.07 -7.64 -7.72 -4.16 -9.38 -5.66 -5.55
-2.72 -2.76 -2.70 -2.65 -4.41 -0.74 -4.65 -5.06 -1.52 -6.21 -2.52 -2.94
-152.869 171 -153.376968 -153.390496 -153.417 873 -150.941 894 -152.948 199 -151.834 863 -152.947 226 -154.957 041 -151.727 234 -153.733 854 -153.785 055
-2.49 -3.13 -3.08 -3.17 -4.08 -0.50 -4.28 -4.69 -1.26 -7.87 -2.32 -2.65
-2.30 -2.56 -2.55 -2.52 -3.76 -0.27 -3.98 -4.37 -0.98 -7.52 -2.03 -2.39
-152.912 774 -153.499 310
-2.31 -3.13
-2.26 -2.77
-150.995 479 -152.999 331 -151.888 745 -153.000952 -155.008 546 -151.782 208 -153.786 273 -153.834 697
-3.71 -0.29 -3.96 -4.36 -1.07 -5.44 -1.99 -2.46
-3.69 -0.24 -3.93 -4.34 - 1.02 -5.40 - 1.93 -2.41
All calculations were carried out with the 6-31++G(2d,2p) basis sets. Total energy in au. Interaction energies in kcal/mol.
TABLE 10: Total Energies Eb and Interaction Energies AEc and AEc{ Calculated for the Minimum Energy Structure of CH4-H20 CC-PVDZ methoda
HF MP2 MP3 MP4SDTQ S-null B-null
X, s-VWN B-VWN s-LYP B-LYP B3-LYP a
E -116.198 530 -116.558 572 -116.584 994 -116.595 396 -114.652 176 -116.194 886 -115.357 657 -116.303 121 -117.849 074 -115.295 306 -116.837 718 -1 16.904 642
AE -1.08 -1.99 -1.80 -1.97 -3.77 -1.80 -3.54 -3.50 -1.65 -4.53 -2.55 -2.28
6-3 1++G(2d,2p) AEcp 0.17 -0.07 -0.06 -0.07 -0.93 1.04 -1.01 -1.23 0.62 -1.88 0.09 0.08
E -116.216 541 -116.614 792 -116.638 527 -116.652 046 -114.686708 -116.232 547 -115.389 947 -116.333 115 - 117.882 704 -115.328 835 -1 16.874 423 -116.935 793
AE 0.1 1 -0.50 -0.48 -0.57 -1.07 1.03 -1.12 -1.36 0.58 -2.02 0.09 -0.12
aug-cc-pVTZ AEcu 0.20 -0.21 -0.22 -0.24 -0.81 1.27 -0.87 1.11 0.82 -1.76 0.33 0.08
E -116.250 819 -116.715 383 -116.737 170 - 116.754 740 -114.728 065 -116.273 015 -115.431 564 -116.374512 - 117.923 302 -115.370914 -116.915 611 -1 16.974 622
AE 0.17 -0.56 -0.57 -0.63 -0.84 1.16 -0.92 -1.15 0.72 -1.76 0.25 -0.02
AEcu 0.21 -0.36 -0.38 -0.43 -0.78 1.25 -0.85 - 1.07 0.82 - 1.70 0.33 -0.10
All calculal:ions were carried out with the 6-31++G(2d,2p) basis sets. Total energy in au. Interaction energies in kcal/mol.
TABLE 11: Total Energies Eb and Interaction Energies AEc and A&{ NH3-NH3 CC-pVDZ methoda
HF MP2 MP3 MP4SDTQ S-null B-null XCI S-VWN B-VWN s-LYP B-LYP B3-LYP
Calculated for the Minimum Energy Structure of
6-3 1++G(2d,2p)
aug-cc-pVTZ
E
AE
AEcp
E
AE
AEcp
E
AE
AEcp
-112.396 134 -112.772319 -112.798 576 -112.810 333 -110.889027 -112.411 981 -111.580302 -1 12.529 126 -114.055 041 -111.528491 -113.051 217 -113.114027
-3.17 -4.22 -3.98 -4.11 -6.41 -3.24 -6.30 -6.39 -3.37 -7.68 -4.50 -4.38
-1.88 -2.25 -2.17 -2.16 -3.76 -0.56 -3.94 -4.34 -1.29 -5.29 -2.08 -2.34
-112.413 166 -112.827 712 -112.850875 -112.865 304 -110.923906 -112.447906 -111.612075 - 112.557 836 -114.085 486 -111.561 755 -113.085 595 -113.143055
-2.11 -3.08 -2.97 -3.09 -4.18 -0.94 -4.32 -4.72 -1.66 -5.67 -2.43 -2.67
-1.96 -2.64 -2.57 -2.60 -3.86 -0.62 -4.00 -4.38 -1.32 -5.31 -2.08 -2.36
-112.444683 -112.926294 -112.947 153 -112.966 677 -110.965 166 -112.487 948 -111.653 831 - 112.599 402 -114.125 873 -111,603923 -113.126503 -113.181 327
-1.86 -2.99 -2.94 -3.01 -3.86 -0.74 -4.06 -4.42 -1.47 -5.30 -2.17 -2.42
-1.82 -2.76 -2.72 -2.79 -3.73 -0.60 -3.92 -4.26 -1.29 -5.16 -2.01 -2.28
All calculations were carried out with the 6-31++G(2d,2p) basis sets. Total energy in au. Interaction energies in kcal/mol.
the case of the interaction energy curves, that is, the largest gap between the 6-31++G(2d,2p) and aug-cc-pVTZ and the cc-pVDZ basis sets corresponds to all the points computed using DFT approximations. In general, the BSSE plays an important role in the DFT computations. In the case of the three basis sets employed here, as the basis set quality increases (cc-pVDZ 6-31++G(2d, 2p) aug-cc-pVTZ), the BSSE decreases, although, as in the MPn case, it does not become zero. When the BSSE is corrected, as in the AE,, curves plotted in Figure 5 for the HF dimer, the results become almost identical for the three basis sets and are very close to the uncorrected aug-cc-pVTZ curve (not shown in the plot). In general, it seems that for the selected dimers in this study, CP corrections may be applied to DFT approximations.
-
-
At this point, we focus our attention on the aug-cc-pVTZ results. Figure 6 shows AE plotted against all the methods employed here using the aug-cc-pTVZ basis sets. Within each curve the D l T approximations and MPn methods appear to be consistent, the largest discrepanciescorrespond to the interaction energies predicted with B-null, B-VWN, and S-LYP. Similar trends may be observed for the corrected interaction energy plot (see Figure 7). Although not represented here, the cc-pVDZ and 6-31++G(2d,2p) interaction energies follow a similar behavior to that plotted in Figures 6 and 7 for the aug-cc-pVTZ values. Therefore, it is possible to extract general conclusions on how D l T interaction energies compare to the MPn values. At this point, we note that our results are in good agreement with previous o b ~ e r v a t i o n s , ' ~ ~ - that ~ . ~ in - ~all * the complexes MP2 and MP4 results are remarkably similar and independent
J. Phys. Chem., Vol. 99,No. 43, 1995 15843
Evaluation of the Density Functional Approximation 01
Figure 3. AE interaction energy computed for the HF dimer with the (0)cc-pVDZ, (0)6-31++G(2d,2p), and (A) aug-cc-pVTZ basis sets at the Hartree-Fwk and MPn methods and indicated DFT approximations.
21
Figure 6. Variation of the AE values computed with aug-cc-pVTZ basis sets: (0)HF-HF; (0)HzO-HzO; (A) NH3-NH3; ( 0 )C2HzHzO; (B) C&-H20 L
'1 a
1
0
Figure 4. BSSE computed using the counterpoise method for the HF dimer with the (0)cc-pVDZ, (0)6-31++G(2d,2p), and (A) aug-ccpVTZ basis sets at the Hartree-Fock, and MPn methods and indicated DFT approximations.
Figure 5. AEcp interaction energy computed for the HF dimer with the (0)cc-pVDZ, (0)6-31++G(2d,2p), and (A) aug-cc-pVTZ basis sets at the Hartree-Fwk, MPn methods and indicated DFT approximations.
of the basis set employed. At the same time, the MPn results are very similar to the available experimental results. When all the DFT results are compared to the MP2 or MP4 values, we find some interesting tendencies. First of all, we see that independent of the correlation functional employed, all the local functionals (S and X,) provide AE and AEcp values larger than the MPn values, while the functionals that include gradient corrections give smaller values. The inclusion of the VWN local correlation functional increases the interaction energy in all the cases, and the interaction energy is further
Figure 7. Variation of the AEcpvalues computed with aug-cc-pVTZ basis sets: (0)HF-HF; (0)HzO-HzO; (A) NH3-NH3; (0)C2H2HzO;(B) CW-HzO.
increased when the LYP functional is used. The best agreement between the DFT functionals and the MPn interaction energies is found for the hybrid B3-LYP method and the B-LYP, in this order. Interestingly enough, however, the S-null results are also reasonable in some of the cases. The B3-LYP interaction energies are shifted in general by a value of + O S kcdmol to a smaller value, a fact that makes the CW-HzO practically unstable at this level of theory. This is probably one of the weakest hydrogen bonds available and by substitution of a hydrogen by another group is likely to greatly increase the stability of the C-H-0 contact, thus allowing the use of the B3-LYP functional for a correct description of the attractive nature of the contact. We mention here that we have studied the effect of increasing the precision of the numerical integration by increasing the number of points in the grid (by using the tight option), and it is found to be 0.02 kcdmol in the (2%H20 case, ruling out a precision problem as the reason for the negligible stabilization of this complex. We should also mention here that the size of the shift is in good agreement with that reported by Kim and Jordan for the water dimer in their computations.'" The BSSE does not change the presence or the average size of the shift, so it must be attributed to the functional itself.
Conclusion The computations presented here on the HF-HF,Hz0-Hz0, CzH2-H20, Ch-HzO, and NH3-NH3 complexes indicate that the B-LYP and B3-LYP functionals give results closer to experimental and M E values. For this series of systems, B3-
15844 J. Phys. Chem., Vol. 99, No. 43, 1995
LYP functionals predict geometries and interaction energies that are in better agreement with MP2 geometries and MP4 energetics. The results also show a strong basis sets effect. At least triple-zeta-type basis set are required to reduce the importance of the BSSE in the computed interaction energies. The counterpoise-correctedvalues in all the cases are very close to the best computed uncorrected results. They follow trends similar to the ones obtained at the MP2 level. The basis set is also important to obtain good optimized geometries. Finally, the interaction energies and the optimized geometries show the presence of a general trend for the exchange and correlation functionals with and without gradient corrections: the inclusion of only gradient-corrected exchange functionals tends to overestimate X-Y (X, Y = C, N, 0, and F) distances and underestimate interaction energies. On the other hand, the inclusion of only gradient-corrected correlation functionals shows the opposite effect. For these systems, it has been found that gradient corrections for exchange and correlation are important to compute geometries and interaction energies.
Acknowledgment. We thank the management of Cray corporate network for generous allocation of computer time on the Cray C-90 and Y-MP machines. We also thank Professor K. Jordan for providing us with a preprint of his work prior to publication. We especially thank Dr. Chengteh Lee for his valuable comments and for reading this manuscript. J.J.N. also thanks the “Fundacib Catalana per a la Recerca” (Catalunya, Spain) for making possible his visit to Cray Research and the Spanish DGICYT for continuous support under Project PB920655-C02-02. He also expresses his gratitude to the people and management of the chemistry applications group for support and a nice working environment while visiting Cray. References and Notes (1) (a) Desiraju, G. R. Crystal Engineering. The Design of Organic Solids; Elsevier: Amsterdam, 1989. (b) Structure Correlation; Dunitz, J. D. Burgi, H. B. Eds.; VCH: Heidelberg, 1993. (2) Jeffrey, G. A,; Saenger, W. Hydrogen Bonding in Biological Structures; Springer-Verlag: Berlin, 1991. (3) (a) Cohen, R. C.; Saykally, R. J. Annu. Rev. Phys. Chem. 1991, 42, 369. (b) D. J. Nesbitt, Chem. Rev. 1988, 88, 843. (c) A. C. Legon, Chem. Soc. Rev. 1990, 19, 197. (4) Buckingham, A. D.; Fowler, P. W.; Hutson, J. M. Chem Rev. 1988, 88,963, and references therein. (a) van Lenthe, J. H.; van Duijneveldt-van de Rijdt, J. G. C. M.; van Duijneveldt, R. F. In Ab Initio Methods in Quantum Chemistry;Lawley, K. P., Ed.; Wiley: New York, 1987; Vol. 11. (b) Hobza, P.; Zahradnik, R. Intermolecular Complexes; Elsevier: Amsterdam, 1988. (c) Hobza, p.; Zahradnik, R. Chem. Rev. 1988,88, 871. (d) Chalasinski, G.; Gutowski, M. Chem. Rev. 1988, 88, 943. (e) Liu, B.; McLean, A. D. J. Chem. Phys. 1989, 91, 2348, and references therein. (f) Scheiner, S. Reviews in Computational Chemistry; VCH: New York, 1991; Vol. 2, p 165. ( 5 ) (a) Parr, R. G.; Yang, W. Density-Functional Theory of Atoms and Molecules; Oxford University Press: New York, 1989. (b) M. Levy, Proc. Narl. Acad. Sei. U.S.A. 1979, 76, 6062. (c) Ziegler, T. Chem. Rev. 1991, 91, 651. (d) Jones, R. 0.;Gunnarsson, 0. Rev. Mod. Phys. 1990, 61, 689. (6) (a) Johnson, B. G.; Gill, P. M. W.; Pople, J. A. J. Chem. Phys. 1993, 98, 5612. (b) For a discussion of the scaling properties of several methods, see, for example, the following Scuseria, G. E.; Lee, T. J. J. Chem. Phys. 1990, 93, 5851. Strout, D. L.; Scuseria, G. E. J. Chem. Phys. 1995, 102, 8448. Foresman, J. B.; Frisch, E . In Exploring Chemistry with Electronic Structure Methods: A Guide to Using Gaussian; Gaussian Inc.: Pittsburgh, PA, 1993. (7) (a) Satoko, C. Chem. Phys. Lett. 1981.83, 111. (b) Satoko, C. Phys. Rev. 1984, B30, 1754. (c) Averill, F. W.; Painter, G. S. Phys. Rev. 1985, B32, 2141. (d) Averill, F. W.; Painter, G. S. Phys. Rev. 1986, 834, 2088. (e) Verluis, L.; Ziegler, T. J. Chem. Phys. 1988, 88, 322. (f) Foumier, R. Andzelm, J.; Salahub, R. J. Chem. Phys. 1989, 90, 6371. (8) Delley, B. In Density Functional Methods in Chemistry; Labanowski, J.; Andzelm, J. Eds.; Springer: New York, 1991; p 101. (h) Johnson, B. G.; Frisch, M. J. Chem. Phys. Lett. 1994,216, 133. (i) Johnson, B. G.; Frisch, M. J. J. Chem. Phys. 1994, 100, 7429. (8) (a) Stoll, H.; Pavlidou, C. M. E.; Preuss, H. Theor. Chim. Acta 1978, 49, 143. (b) Langreth, D. C.; Perdew, J. P. Phys. Rev. 1980, B21, 5469. (c) Langreth, D. C.; Mehl, M. J. Phys. Rev. 1983, 828, 1809. (d) Becke, A. D. Int. J. Quantum Chem. 1983, 27, 1915. (e) Perdew, J. P.
Novoa and Sosa Phys. Rev. Lett. 1985, 55, 1655. (f) Wilson, L. C.; Levy, M. Phys. Rev. 1990,841, 12930. (g) DePristo, A. E.; Kress, J. D. J. Chem. Phys. 1987, 86, 1425. (9) (a) Perdew, J. P. Phys. Rev. 1986, 833, 8822. (b) Lee, C.; Yang, W.; Parr, R. G. Phys. Rev. 1988, B37, 785. (c) Perdew, J. P.; Wang, Y. Phys. Rev. 1992, B45, 13244. (d) Perdew, J. P. Physica 1991, B172, 1. (e) Wang; Y.; Perdew, J. P. Phys. Rev. 1991, B44, 13298. (f) Perdew, J. P. In Electronic Structure of Solids ‘91,Ziesche, P. Eschrig, H., Eds.; Akademie Verlag: Berlin, 1991; p 11. (10) (a) Hohenberg, P.; Kohn, W. Phys. Rev. 1964, B136, 864. (b) Kohn, W.; Sham, L. J. Phys. Rev. 1965, A140, 1133. (11) (a) Baerends, E. J.; Ellis, D. E.; Ros, P. Chem. Phys. 1973, 2, 41. (b) Delley, B. J. Chem. Phys. 1990,92,508.(c) Salahub, D. R. In Ab initio Methods in Quantum Chemistry-II; Lawley, K. P., Wiley: New York, 1987; p 447. (d) Wimmer, E.; Andzelm, J. J. Chem. Phys. 1992, 96, 1280. (12) (a) Becke, A. D. J. Chem. Phys. 1986.84.4524. (b) Becke, A. D. In The Challenge of d and f Electrons: Theory and Computation; Salahub, D. R., &mer, M. C. Eds.; ACS Symposium Series No. 394; American Chemical Society: Washington D.C., 1989. (c) Becke, A. D. J. Chem. Phys. 1988.88, 2547. (13) (a) Andzelm, J.; Sosa, C.; Eades, R. A. J. Phys. Chem., 1993, 97, 4664. (b) Sosa, C.; Andzelm, J.; Lee, C.; Blake, J. F.; Chenard, B. L.; Butler, T. W. Int. J. Quantum Chem. 1994,49,511. (d) Sosa, C.; Lee, C. J. Chem. Phys. 1993, 98, 8004. (14) (a) Schwenke, D. W.; Truhlar, D. G. J. Chem. Phys. 1985, 82, 2418; 1986, 84, 4113; 1987, 86, 3760. (b) Collins, J. R.; Gallup, G. A. Chem. Phys. Lett. 1986,123, 56. (c) Frish, M. J.; Del Bene, J. E.; Binkley, J. S.; Schaeffer, H. F., 111. J. Chem. Phys. 1986, 84, 2279. (d) Olivares del Valle, F. J.; Tolosa, S.; Ojalvo, E. A.; Espinosa, J. J. Chem. Phys. 1987, 127, 343. (e) Kim, K.; Jordan, K. D. J. Phys. Chem. 1994, 98, 10089. (0 Olivares del Valle, F. J.; Tolosa, S.; Ojalvo, E. A,; Espinosa, J. Chem. Phys. 127, 343. (E) Tolosa, S.; Esuirilla, J. J.; Espinosa, J.; Olivares del Valle, F. J. Chem. phys. 1988, 127,-65. (h) Odutoia, J. A.; Dyke, T. R. J. Chem. Phys. 1980, 72, 5062. (15) Boys, S. F.; Bemardi, F. Mol. Phys. 1970, 19, 553. (16) Pine, A. S.; Howard, B. J. J. Chem. Phys., 1986, 84, 590. (17) (a) Curtiss, L. A.; Frurip, D. J.; Blander, M. J. Chem. Phys. 1979, 71, 2703. (b) Reimers, J.; Watts, R.; Klein, M. Chem. Phys. 1982, 64, 95. (18) Gutowski, M.; Chalasinski, G. J. Chem. Phys. 1993, 98, 5540. (19) (a) Gutowski, M.; van Duijneveldt-van de Rijdt, J. G. C. M.; van Lenthe, J. H.; van Duijneveldt, F. B. J. Chem. Phys. 1993, 98, 4728. (b) Tao, F.-M.; Pan, Y.-K. J. Phys. Chem. 1991, 95, 3582. (20) Feller, D. J. Chem. Phys. 1992, 96, 6104, and references therein. (21) (a) Novoa, J. J.; Planas, M.; Whangbo, M.-H. Chem. Phys. Lett. 1994, 225, 240. (b) Novoa, J. J.; Planas, M.; Rovira, C. Submitted for publication. (22) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Gill, P. M. W.; Johnson, B. G.; Wong, M. W.; ForesmanJ. B.; Robb, M. A,; Head-Gordon, M.; 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/DFT, revision G.3; Gaussian, Inc.: Pittsburgh PA, 1993. (23) Gill, P. M. W.; Johnson, B. G.; Pople, J. A. Chem. Phys. Lett. 1993, 209, 506. (24) Johnson, B. G.; Frisch, M. J. Chem. Phys. Lett. 1993, 216, 133. (25) (a) Dirac, P. A. M. Proc. Cambridge Philos. SOC.1930, 26, 376. (b) Slater, J. C. The Self-consistent Field for Molecules and Solids; Quantum Theory of Molecules and Solids, Vol 4; McGraw-Hill: New York, 1974. (26) Vosko, S. H.; Wilk, L.; Nusair, M.; Can J. Phys 1980, 58, 1200 . (27) Slater, J. C. Phys. Rev., 1951, 81, 385. (28) Becke, A. D. Phys. Rev. A 1988, 38, 3098. (29) Lee, C.; Yang, W.; Parr, R. G. Phys. Rev. B 1993, 37, 785. (30) Colle, R.; Salvetti, D. Theor. Chim. Acta 1975, 37, 329. (31) The B3-LYP is an exchange-conelation hybrid functional which includes a mixture of, The B3-LYP. Hartree-Fock exchange with VWN exchange and LYP correlation functional according to the following expression P 3 L y p = (1 - ao)ExLsD a a X H F axAExBss acEcLYP (1 - ac)ECwN(See the following. (a) GaussiadDFT supplement manual. (b) Becke, A. D. J . Chem. Phys. 1993,98, 1372; 98,5648. (c) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. J. Phys. Chem. 1994, 98, 1623. (32) (a) Dunning, T. H., Jr. J. Chem. Phys. 1989,90, 1007. (b) Kendall, R. A.; Dunning, T. H., Jr. J. Chem. Phys. 1992, 96, 6796. (33) (a) Hehre, W. J.; Ditchfield, R.; Pople, J. A. J. Chem. Phys. 1972, 56, 2257. (b) Hariharan, P. C.; Pople, J. A. Mol. Phys. 1974, 27, 209. (c) Frisch, M. J.; Pople, J. A.; Binkley, J. S. J. Chem. Phys. 1984, 80, 3265. (c) Clark, T.; Chandrasekhar, J.; Spitznagel, G. W.; Schleyer, P. v. R. J. Comp. Chem., 1983, 4, 294. (34) (a) Moller, C.; Plesset, M. S. Phys. Rev. 1934,46, 618. (b) Pople, J. A.; Seeger, R.; Krishnan, R. Inr. J. Quantum Chem. Symp. 1977, 11, 149. (35) (a) Dyke, T. R.; Mack, K. M.; Muenter, J. S. J. Chem. Phys. 1977, 66, 498. (b) Odutola, J. A.; Dyke, T. R. J. Chem. Phys. 1980, 72, 5062. (36) Howard, B. J.; Dyke, T. R.; Klemperer, W. J. Chem. Phys. 1984, 81, 5417.
+
+
+
+
Evaluation of the Density Functional Approximation (37) (a) Peterson, K. I.; Klemperer, W. J. Chem. Phvs. 1984.81, 3842. (b) Block, P. A.; Marshall, M. D:; Pedersen, L. G.; Miiler, R. E. J. Chem. Phys. 1992, 96, 7321. (38) (a) Seiler, P.; Weisman, G. R.; Glendening, E. D.; Weinhold, F.; Johnson, V. B.; Dunitz, J. D. Angew. Chem., lnt. Ed. Engl. 1987,26, 1175. (b) Seiler, P.; Dunitz, J. Helv. Chim Acta 1989, 72, 1125. (39) AE = E(AB)AB- E(A)A - E(B)B, where E(AB)AB,E(A)A, and E(B)B are the energies of A-B, A, and B obtained using the truncated basis sets located on A B, A, and B, respectively. AEcp= E(AB)ABE(A)AB- E(B)AB,where the truncated basis set of A B is used on all the computations. BSSE = AEcp- AE. (40) Racine, S. A.; Davidson, E. R. J. Phys. Chem. 1993, 97, 6367. (41) Knowles, P. J.; Somasundram, K.; Handy, N. C.; Hirao, K. Chem. Phys. Lett. 1985, 113, 8. (42) Frish, M. J.; Pople, J. A.; Del Bene, J. E. J. Phys. Chem. 1985,89, 3664. (43) Hehre, W. J.; Radom, L.; Schleyer, P. Pople, J. A. Ab Initio Molecular Orbital Theory; John Wiley & Sons: New York, 1986. (44) (a) Michael, D. W.: Dykstra, C. E.; Lisy, J. M. J. Chem. Phys. 1984, 81, 5998. (b) Gaw, J. F. Yamaguchi, Y. Vincent, M. A,; Schaefer, H. F. J. Am. Chem. SOC. 1984, 106, 3133. (c) Benzel, M. A.; Dykstra, C.
+
+
J. Phys. Chem., Vol. 99,No. 43, 1995 15845 E.J. Chem. Phys. 1993, 78,4052. 1984, 80,35lO(E).(d) Mile, F.; Mineva, T.; Russo,N.; Toscano, M. Theoret. Chim. Acta, accpeted for publication. (45) (a) Sosa, C.; Lee, C.; Fitzgerald, G.; Eades, R. A. Chem. Phys. Lett. 1993, 211, 265. (b) Sosa, C.; Novoa, J. J.; Lee, C. In preparation. (46) Del Bene, J. E. J. Chem. Phys. 1987, 86, 2110. (47) (a) Sagarik, K. D.; Ahlrichs, R.; Brode, S. Mol. Phys. 1986, 57, 1247. (b) Hassen, D. M.; Marsden, C. J.; Smith,B. J. Chem. Phys. Lerr. 1991,183,449.(c) Latajka, Z.; Scheiner, S. J. Chem. Phys. 1986,84,341. (48) (a) Nelson, D. D., Jr.; Fraser, G. T.; Klemperer, W. J. Chem Phys. 1985,83,6201. (b) Nelson, D. D., Jr.; Klemperer, W.; Fraser, G. T.; Lovas, F. J.; Suenram, R. D. J. Chem. Phys. 87, 6364. (c) Tao, F. M.; Klemperer, W. J. Chem. Phys. 99, 5976. (d) Olfiof, E. H. T.; van der Avoird, A.; Wormer, P. E. S. J. Chem. Phys. 1994, 101, 8430. (49) (a) Muguet, F. F.; Robinson, G. W. J. Chem. Phys. 1995, 102, 3648. (b) Muguet, F. F.; Robinson, G. W.; Basses-Muguet, M. P. J. Chem. Phys. 1995, 102, 3655. (50) Benedict, W. S.;Plyler, E. K.; Can J. Chem. 1957, 35, 1235. (51) Cohen, E. A.; Poynter, A. L.; J. Mol. Spectrosc. 1974, 53, 131. Jp95 10220
