6630
J . Phys. Chem. 1992, 96,6630-6636
A Local Density Functional Study of the Structure and Vibrational Frequencies of Molecular Transition-Metai Compoundst Carlos Sosa, Jan Andzelm,t Brad C. Elkin, Erich Wimmer,**t Cray Research, Inc., 655 E . Lone Oak Dr.. Eagan, Minnesota 55121
Kerwin D. Dobbs,* and David A. Dixon* Du Pont Central Research and Development, P.O. Box 80328, Experimental Station, Wilmington, Delaware 19880-0328 (Received: December 26, 1991)
Local density functional (LDF) theory has been used to calculate the geometry and vibrational frequencies of a set of transition-metal compounds in their molecular forms containing halogens, oxygens, alkyl groups, carbonyls, nitrosyls, and other substituents. The calculations were done with polarized doublezeta numerical and Gaussian basis sets, and the geometries were obtained by analytic gradient methods. The frequencies were evaluated by numerical differentiation of the analytic first derivatives. The results obtained with the numerical and the Gaussian basis sets were found to be in good agreement. The agreement with experiment for the geometries is quite good with an average mean deviation of 0.026 A. The largest errors involve dative bonds with the LDF method predicting the bonds to be too short. Nonlocal corrections were applied to some of the methyl-carbonyl and metal-nitrosyl bond lengths, and this correction was found to lengthen the bonds to give better agreement with experiment. The frequencies are also predicted quite accurately. The LDF results are in much better agreement with experiment as compared to HartreeFock results.
Introduction The structure of transition-metal systems has proven to be a challenge for traditional molecular orbital methods based on Hartree-Fock theory. Even molecular geometries which are usually handled well by Hartree-Fock theory are difficult to calculate accurately for transition-metal systems. Besides the problems inherent in applying the Hartree-Fock SCF method to transition-metal systems is the additional problem of the size of transition-metal systems. There are many possible hybridization schemes for transition metals, and the types of ligands, that are bonded to the metal may, in themselves, be of significant size. Hartree-Fock methods scale formally as iV'where N is the number of basis functions although in practice the scaling can be lower ( N 3 - W 5 ) . Methods that include electron correlation effects are often required to treat transition-metal systems, and these methods scale as N" with m 1 5. It would be of real benefit to have a method that can be used to treat transition-metal systems that gives good geometries and vibrational frequencies and scales with a lower exponent for the number of basis functions. One method that has been applied with great success by the solid-state physics community to treat the electronic structure of metallic systems (including transition metals) is the local density functional (LDF) method.'-s This method is based on rigorous theoretical grounds and has no empirical parameters. LDF theory has been applied to transition-metal systems to predict the energetics and structures along the reaction path for a number of organometallic processes with one set of appro~imations.6~~ We have recently shown that LDF theory with either numerical basis sets or Gaussian basis sets can be used to calculate the structure, vibrational spectra, and energetics of a number of first-row molecules including those such as FOOF which require a correlated treatmenL8-I0 Thus, it is appropriate to test the ability of these formulations of the LDF method to treat the structures and vibrational spectra of a set of transition-metal compounds involving first- and second-row transition metals. There is one consistent Hartree-Fock study"-'3 of the structures and vibrational frequencies of a range of such molecules far which experimental data are available. We have chosen this set for our study. Furthermore, this set of compounds contains most of the simple metal oxides 'Contribution No. 6033. 'Present address: BIOSYM Technologies, 10065 Barnes Canyon Rd., San Diego, CA 921 21-2717.
and halides which are important as models for active catalysts.
Metbods The calculations described below were done with the program systems DMo14 and D G a ~ s s ~ von ' ~Cray J ~ YMP computers in a single processor mode. LDF theory is based on a theorem of Hohenberg and KohnI7 which states that the total energy E, is a functional of the charge density, p Elbl = TIPI + UtPl + Exctpl
(1)
where T is the kinetic energy of the noninteracting electrons of density p, U is the classical Coulomb electrostatic energy, and E,, includes all of the many-body contributions to the energy. For a molecular system, the most important terms in E, are the exchange and correlation energies. The density is given as a sum over the squares of the occupied molecular orbitals as in standard molecular orbital treatments. The first term in eq 1 is the kinetic energy of the electrons. The second term can be written as
PI = ( P ~ N )+ ( p V e / 2 ) + VNN
(2)
where the first term in ( 2 ) corresponds to the nuclear attraction energy, the second term to the electron-electron repulsion, and the third term to the nuclear-nuclear repulsion. These quantities and the kinetic energy term can all be evaluated using straightforward techniques. Up to this point eq 1 is exact, and it is in the last term of eq 1 that the local density approximation is introduced. It has been found that a good approximation for the final term can be taken from the exchange-correlation energy of the uniform electron gas.*8.19The physical assumption for this approximation is that the charge density varies slowly on the scale of exchange and correlation effects and one can obtain E,, by integrating the uniform electron gas result Excbl
= Jdrlexc[P(r)l dr
(3)
with e,, being the exchangecorrelation energy per particle in the uniform electron gas. In the calculations described below, the form derived by von Barth and Hedido is used in DMol and that of Vosko, Wilk, and Nusair2' in DGauss. The value of E, is determined by optimizing the variations in E, with respect to variations on p subject to the usual orthonormality constraints for molecular orbitals. This leads to a set of coupled differential equations which'must be solved iteratively until the density and
0022-3654/92/2096-6630$03.00/00 1992 American Chemical Society
The Journal of Physical Chemistry, Vol. 96, No. 16, 1992 6631
Molecular Transition-Metal Compounds TABLE I: Gaussian Basis Sets for DGauss basis H C-F C1 Sc-Zn YCd
(4 1I/PI (621/41/ 1)/[3/2/11 (6321 /411/1)[4/3/1] (63321/531/41)/[5/3/2] (633321 15321 1/53 1)/ [6/5/3]
aux basis (4) (71313) (9/4/4) (10/5/5) (10/5/5)
the energy are converged just as in standard molecular orbital treatments. DMol employs numerical functions for the atomic basis sets. The atomic basis functions are given numerically on an atomcentered, spherical-polar mesh. The radial portion of the grid is obtained from the numerical solution of the atomic LDF equations. The use of exact spherical atom results offers some advantages. The molecule will dissociate exactly to its atoms within the LDF framework, although this does not guarantee correct dissociation energies. Furthermore, because of the quality of the atomic basis sets, basis set superposition effects should be minimized, and correct behavior at the nucleus is obtained. Since the basis sets are numerical, the various integrals arising from the expression for the energy equation are evaluated over TABLE II: Geometries for First-Row Transition-Metal Compounds" molecule parameter STO-3Gb r(ScF) 1.845 r(ScF) 1.852 1.701 r(TiF) r(TiC1) 2.167 r(TiC) r(TiC1) L(C1TiC) L(HCTi) (VO) r(VF) L(FVF) r(VO) r(VC1) L(C1VCI) r(VF*,) r(VF,) rWN) r(VC1) r(NC1) L(CINV) L(NVC1) L(C1VCI) r(Cr0) r(CrF) f(OCr0) f(FCrF) r(Cr0) r( CrCI) L(OCr0) L(CICrC1) r(Cr0)
r( MnO) r(NiC) r(NiN) r(CC) r(NO) r(CH) r(CuH) r(CuF) r(CuC1) r(CuC) r(CW f(HCH) r(ZnC) r(CH) L(HCH) r(ZnC) L(HCH)
1.494 1.619 108.7 1.465 2.108 110.1 1.648 1.616 1.526 2.111 1.715 179.9 109.9 109.0 1.442 1.584 109.6 106.8 1.425 2.091 106.6 110.2 1.520 1.455 2.08 1 1.424 1.420 1.268 1.078 1.595 2.07 1
1.682 1.083 105.6 1.774 1.1 11 100.2
a grid in DMol. The grid is generated in terms of radial functions and spherical harmonics. The number of radial points NR is given as
NR = 1.2
X 14(2
+ 2)1/3
(4) where 2 is the atomic number. The maximum distance for any function is 12 au. The angular integration points No are generated at the NRradial points to form shells around each nucleus. The value of N8 ranges from 14 to 302 depending on the behavior of the density." The Coulomb potential corresponding to the e l m o n repulsion term is found by solving Poisson's equation for the charge density numerically.
-V2Ve(r) = 4 7 d p ( r )
(5)
All of the DMol calculationswere done with a double numerical basis set augmented by polarization functions. For comparison to traditional molecular orbital methods, this can be considered in terms of size as a polarized double-zeta basis set. However, because exact numerical solutions for the atom are used, this basis set is of significantlyhigher quality than a normal molecular orbital double-zeta basis set. The multipolar fitting functions for the model density used to evaluate the effective potential were all
3-21GE 1.775 1.812 1.719 2.186 2.000 2.226 1.086 102.1 108.5 1.520 1.687 110.4 1.511 2.157 110.7 1.695 1.673 1.553 2.186 1.705 178.9 106.6 112.2 1.506 1.670 108.6 109.8 1SO2 2.160 107.3 112.7 1.604 1.552 2.132 1.530 1.419 1.168 1.066 1.382 1.605 2.093 1.901 1.098 105.5 1.944 1.088 106.9 2.154 1.104 105.6
DMol 1.795 1.833 1.747 2.167 2.012 2.174 1.105 106.3 108.1 1.574 1.721 110.7 1.573 2.131 110.8 1.747 1.712 1.639 2.138 1.618 175.6 107.0 112.6 1.567 1.702 108.3 110.4 1.568 2.099 109.2 110.2 1.661 1.611 2.077 1.607 1.411 1.179 1.082 1.465 1.730 2.035 1.927 1.115 106.2 1.912 1.107 107.4 2.124 1.117 106.4
DGauss 1.787 1.831 1.740 2.160 2.008 2.170 1.108 104.8 108.7 1.567 1.712 110.8 1.568 2.124 110.6 1.740 1.706 1.632 2.133 1.613 175.9 106.2 112.4 1.562 1.699 108.0 110.8 1.565 2.099 109.4 111.6 1.656 1.606 2.053 1.607 1.427 1.177 1.095 1.456 1.723 2.030 1.917 1.114 106.2 1.907 1.107 107.9 2.126 1.116 106.6
exptd 1.787 1.91 1.745 2.170 2.047 2.185 1.098 105.6 109.0 1.569 1.729 11 1.2 1.570 2.142 111.3 1.734 1.703 1.651 2.138 1.597 169.7 106.0 113.4 1.575 1.720 107.8 111.9 1.58 1 2.126 108.5 113.3 1.66 1.629 2.1 1 1.626 1.43 1.165 1.09 1.462 1.745 2.05 1 1.935 1.929 1.090 107.7 2.07
Bond distances in angstroms. Bond angles in degrees. Reference 13. CReference1 1. "Experimental results summarized in ref 1 1 unless noted.
6632 The Journal of Physical Chemistry, Vol. 96, No. 16, 199'2
generated by using an angular momentum number, C, one greater than that of the polarization function; thus, a value of C = 3 was used for the fitting function with d polarization functions, and C = 2 was used for H which had p polarization functions. In order to use an efficient analytical a p p r ~ a c h ? ~ 'Gaussian ~*'~~~ basis sets are required. This allows one to use the wealth of experience gained from Hartree-Fock molecular orbital calculations. Cartesian Gaussians are used as primitives in DGauss and are contracted in the same way as in Hartree-Fock methods24 although the actual form of the basis sets is somewhat different. Although the experience about basis sets gained from HartreeFock molecular calculations is invaluable in terms of defining the basis set size and the need for additional polarization and diffuse functions, Hartree-Fock basis sets used in an LDF calculation may exhibit large basis set superpositionerrors (BSSE). In order to minimize BSSE, LDF optimized basis sets have been develThe Gaussian basis sets employed in the DGauss calculations are summarized in Table I. This level of basis set is denoted as doublezeta valence plus polarization (DZVP) and does not include p functions on H.16 In the analytical approach used in DGauss, the total electron density is variationally fit to a set of auxiliary Gaussian-type basis functions leading to exact Coulomb forces.26a The remaining exchange-correlation energy term is a smooth function of the density and can be accurately fit to another auxiliary set of Gaussian-type functions26busing an adaptive set of grid points as described below. It is this numerical evaluation of the exchangecorrelation energy term that leads to the differences in the geometries obtained from the minimum in the energy and the one with zero gradient as discussed below. Auxiliary basis sets are required to fit the electron density as obtained from the molecular orbitals and the exchangecorrelation potential and energy. These additional basis sets are expanded in a set of atom-centered Gaussian-type basis functions with s, p, and d character. It is found that even-tempered expansions of triple-zeta and quadruplezeta quality are adequate to reproduce the electron density. The coefficients in this density expansion are determined from a variational condition that requires that the residual (second order) Coulomb energy term arising from the difference between the exact and the fitted density be minimized. In practice, this leads to analytic expressions for the fitting coefficients involving Coulomb-type three-index, two-electron integrals. This auxiliary basis set is also summarized in Table I. In DGauss, the exchange-correlation potential as well as the exchangmrrelation energy is expanded in Gaussian-type basis functions. However, the determination of the expansion coefficients is carried out numerically using a grid similar to the one described above for the DMol program. The grid selection in DGauss is derived from an adaptive procedure based on the values of the exchangecorrelation energy in angular shells around each atom.15 The radial distribution of points is accomplished by using the method of B ~ k e . ~Once ' the Gaussian expansions of the exchangmrrelation potential are determined, the matrix elements involving the three-index, two-electron integrals are calculated. To this end, the DGauss approach leads to an N3 method which in practice can be made close to M by using sparse, direct matrix algorithms for the integral calculations, numerical integration, and density synthesis. This can be compared to the four-index, tweelectron integrals required for the Hartree-Fock method which formally scale as N4. All three-index, two-electron integrals are calculated by using an algorithm based on the work of Obara and Saika.28 Geometries were optimized by using analytic gradient methods.5~15~25929~30 First derivatives in the LDF framework can be calculated efficiently and for the numerical basis sets only take on the order of 3-4 SCF iterations or 10-25% of the calculation of the energy. This is in contrast to the evaluation of the derivatives in traditional molecular orbital methods which usually takes at least a factor of 1-2 times the evaluation of the energy. There are two problems with evaluating gradients in the LDF framework which are due to the numerical methods that are used. The first
Sosa et al. TABLE I I I : Geometries for Seed-RowTnnsition-Metal Comporulds"
molecule YF YCI ZrF4 ZrCI4 NbCI5
uarameter STO-3Gb r(YF) 1.975 r(YC1) 2.486 r(ZrF) 1.870 r(ZrC1) 2.316 r(NbCI,,) r(NbC1,) MoF, ~(MoF) OMoF4 r(MoF) 1.776 r(MoO) 1.611 L(OMOF) 107.3 OM&& r(MoC1) r(MoO) ~(0MoCl) Ru04 r(Ru0) AgH r(AgH) AgF r(AgF) 1.633 AgCl r(AgC1) 2.083 Cd(CH3)Z r(CdC) 2.003 r(CH) 1.084 L(HCH) 106.2
3-21GC 1.909 2.495 1.876 2.387 2.362 2.338 1.814 1.830 1.645 105.9 2.359 1.638 103.5 1.691 1.782 2.019 2.485 2.220 1.087 108.2
DMol DGauss exdd 1.923 2.353 1.902 2.320 2.330 2.277 1.851 1.856 1.683 105.5 2.301 1.685 104.5 1.716 1.660 1.996 2.306 2.128 1.105 108.0
1.959 2.376 1.915 2.321 2.329 2.287 1.858 1.864 1.683 104.9 2.297 1.690 103.9 1.716 1.657 1.994 2.301 2.119 1.105 108.0
1.926 2.406 1.886 2.32 2.338 2.241 1.820 1.836 1.650 103.8 2.279 1.658 102.8 1.706 1.618 1.983 2.281 2.112 1.09 108.4
Bond distances in angstroms. Bond angles in degrees. bReference 13. cReference 12. dExperimental results summarized in ref 12.
is that the energy minimum does not necessarily correspond exactly to the point with a zero derivative. The second is that sum of the gradients may not always be zero as required for translational invariance. These tend to introduce errors on the order of 0.001 A in the calculation of the coordinates if both a reasonable grid and basis set are used. This gives bond lengths and angles with reasonable error limits. The frequencies were determined using numerical differentiation of the gradients. A two-point difference formula was used with a displacement of 0.01 au.
Ralllts We divide the results as follows. The geometries of compounds with first-row transition metals with ionic bonding are in Table I1 whereas those with second-row elements are in Table 111. The frequencies associated with these structures are in Tables IV and V, respectively. We present the calculated and experimental geometries for molecules with dative bonds in Table VI. The agreement between the LDF geometries and the experimental ones for the first- and second-row transition-metal compounds with strong metal-ligand bonds is excellent, especially when compared to the HF/3-21G geometries. At the HF/3-21G level, essentially all calculated metal-fluorine bond lengths are shorter than corresponding experimental values with the exception being the calculated bond length for AgF. There is no similar trend for the two LDF methods, but these methods better reproduce metal-fluorine bond distances than the HF/3-21G method. For example, HF/3-21G predicts the CuF bond length to be 0.14 A shorter than the experimental value whereas DMol and DGauss both underestimate the bond length by about 0.02 A. A closer look at the LDF calculated parameters reveals a large discrepancy for ScF3 where the calculated Sc-F bond is short by 0.08 A. Considering the excellent agreement found for the other structures, this difference is surprising. The HF/3-2lG method consistently calculated bond distances between chlorine and transition metals to be longer than experimental values. On the other hand, the trend for the LDF methods is to predict bond lengths between chlorine and first-row transition metals to be shorter than experimental ones while no analogous trend exists for the calculated metal-chlorine bond lengths of the second transition series. Again, the LDF methods outperform the HF/3-21G method in predicting metal-chlorine bond lengths. For example, the HF/3-21G bond length for AgCl is significantly longer (by about 0.20 A) than the experimental value whereas the corresponding LDF bond lengths are longer than the experimental one by about 0.02 A. Metal-oxygen bond distances are consistently underestimated at the HF/3-21G level for both transition series, errors ranging from a minimum of 0.005 A (in OMoF4) to a maximum of 0.079 A
Molecular Transition-Metal Compounds
The Journal of Physical Chemistry, Vol. 96, No. 16, 1992 6633
TABLE Iv: Vibratioaal Frequencies for Fitst-Row Trmeition-Metal
ComwPlds”
molecule ScF TiCl,
assignt
str a, str e bend t2 str t2 bend VOF3 a, VO str a, VF str a, FVF bend e VF str e VF30 rock e FVF bend VOCl, a, VO str a, VCI str a1 ClVCl bend e VCI str e VC130rock e ClVCl bend a,’ VF, str VF5 a,‘ VF2 str a; VF2 str a; VF3 bend e’ VF, str e’ VF, bend e’ FVF bend e’ FVF bend Cr02C12 a, CrO str a1 CrCl str a, Cr02 bend a, CrC12bend a2 CrO, twist b, CrO str b, Cr02 rock b2 CrCl str b2 Cr02 wag Ni(CO), a, CO str a, NiC str e NiCO bend e CNiC bend tl NiCO bend t2 CO str t2 NiCO bend t2 NiC str t2 CNiC bend CuH str str CuF str ’ CUCl
HFb
DMol
DGauss
exptC
860 364 130 500 143 1346 871 285 978 35 1 223 1350 412 175 496 279 135 886 769 991 370 994 302 126 387 1274 471 377 141 26 1 1320 228 523 285 2355 445 491 62 342 2270 527 421 122 2020 928 402
729 377 96 507 129 1108 727 249 820 304 201 1094 413 158 521 243 121 720 622 790 333 822 264 109 315 1096 480 360 138 222 1070 209 523 261 2135 436 485 laid 277 2068 515 479 33 1967 655 438
751 380 114 503 118 1121 732 249 827 303 196 1107 417 163 515 249 129 722 626 801 320 819 272 112 329 1072 479 366 144 236 1099 221 523 272 2155 512 486 49 290 2089 491 417 82 1996 652 438
736 389 114 498 136 1055 720 256 801 304 204 1042 408 163 502 246 125 718 608 784 33 1 810 282 110 336 991 470 356 139 224 1002 212 503 257 2154 371 380 62 300 2092 459 423 79 1941 623 415
Frequencies in cm-I. Assignments are approximate. Reference 11. ‘Experimental results summarized in ref 11. dImaginary frequency. Sct text.
(in Cr02C12). The absolute errors for the corresponding LDF bond lengths are smaller and range from a minimum of 0.001 A (in [CrO,]*-) to a maximum of 0.033 A (in OMoF4). The LDF methods calculate metal-carbon lengths which are in much better agreement with experiment (error of about 0.03 A) than those calculated at the HF/3-21G level (error of about 0.05 A). Overall, the mean absolute deviation of HG/3-21G from experimental bond distances (connecting heavy atoms) for the molecules listed in Tables I1 and I11 is 0.052 A. The corresponding deviation for the LDF methods is about 3 times smaller, at 0.019 A. The 3-21G basis set is double-zeta in the valence space for the first- and second-row atoms and hydrogen. For the metals, there are two sets of s, p, and d functions in the valence space which is of the same size as the LDF metal basis sets. Although the LDF basis sets have polarization functions on the substituents except for H in the Gaussian basis sets and the HF basis sets do not, it is unlikely that the significantly better agreement with experiment for the LDF results is due to this difference in basis sets. It is important to note here that somewhat larger differences between the LDF and experimental values are found in going from the first-row to the second-row transition-metal compounds. This is not surprising as we are not employing relativistic corrections which is especially evident in AgH where the bond is calculated
TABLE V Vibratioarl Frequencies for Second-Row Tnnsition-MeW Compoluldso
molecule
assignt str str a, str e bend t2 str bend
HFb DMol DGauss
exaC
739 349 348 99 417 105 393 302 431 177 418 187 42 144 1198 369 163 271 51 213 397 263 127 1005 339 958 328 3162 1354 457 42 3162 1353 520 3235 1656 757
631
a,’ NbC13 str a,’ NbC12 str a; NbC12 str a2“ NbC1, bend e‘ NbC13 str e’ NbC13 bend e’ NbC1 bend e” NbCl bend a, MOO str a, MoCl str a, MoCl, bend b, MoCl str bl OMoCl b2 MoC12 bend e MoCl str e OMoCl bend e ClMoCl bend a, str e bend t2 str t2 bend a,’ CH3 str a,’ CH3 def a,’ CCd str a,” torsion a; CH, str a? CH, def a; CCdC def e’ CH3 str e’ CH, def e’ CH, rock e’ CCdC def 113 e” CH, str 3233 e” CHI def 1634 e” CH3rock 705 str 1422 str 535 Str
248
671 395 371 85 428 87 388 314 398 162 432 171 54 136 1030 384 160 303 18 214 403 236 146 916 309 955 323 3035 1124 469 101 3036 1132 542 3036 1393 677 105 2944 1397 599 1688 511 336
644388 370 93 433 97 384 320 411 174 431 180 57 141 1013 390 166 308 35 223 414 261 151 916 313 953 333 2966 1137 472 121 3067 1134 539 3070 1401 694 129 3067 1396 616 1694 519 337
381 377 98 418 113 394 317 409 148 430 180 54 148 1015 450 143 400. 148* 220. 396 256 172 885 322 921 336 2903 1127 459 0 2923 1136 535 2980 1315 700 124 2859 1427 634 1760 513 343
Frequencies in cm-I. Assignments are approximate. Reference 12. ‘Experimental results summarized in ref 12. Starred values for OMoC1, from OWCIO,. a
to be about 0.04 A too long. In comparing the calculated metric parameters for the two LDF methods, we found very good agreement with most differences less than 0.01 A. The calculated results for the compounds with dative bonds (Table VI), notably the carbonyls, do not show as good agreement with experiment. This is especially apparent with the HF/3-21G geometries with a mean absolute deviation from experiment of 0.105 A. For DMol and DGauss, this deviation, 0.033 A, is much smaller but not as good as that for the molecules in Tables I1 and 111. The major problem is the description of the metal-ligand bond as the geometries of the ligands are treated properly. For the simple carbonyls, the LDF calculations consistently predict metal-carbon bonds that are too short with a mean absolute deviation of 0.042 A. Note that the Cr-C bond in Cr(C0)2(PH3)4 is only 0.006 A too short, but the Cr-P bonds are 0.089 and 0.060 A too short with the numerical basis sets but only 0.066 and 0.039 A with the Gaussian basis sets. Probably the worst agreement between the calculated and experimental structures is found for Mn(CO),NO, where the Mn-N bond distance is 0.131 A too short. The large difference in the calculated and experimental H-Fe-H bond angle for H2Fe(C0)4is probably exaggerated as the experiment does not determine the positions of the H’s well, and the angle has an error of loo associated with it.” A number of calculations of the structure of HCO(CO)~ have been made. A density functional study based on the Hartrec-
6634 The Journal of Physical Chemistry, Vol. 96, No. 16, 1992 TABLE VI: Geometries for First-Row T-itkMetaI molecule C r K O ," L
.
parameter r(CrC)
ricoj
Sosa et al.
Compomda with Lhtive Bonds' STO-3Gb 3-21GC 1.789 1.167
2.003 1.926 1.798 1.124 1.131 1.205 168.1 104.4
1.624 1.732 1.723 1.162 1.163 72.7
2.016 1.643 1.147 1.171
1.583 1.162
1.965 1.986 1.707 1.125 1.124 100.4 151.6 105.4 2.008 1.844 1.124 1.139 1.955 1.894 1.668 1.122 1.130 97.2
1.831 1.131
DMol 1.869 1.157 1.811 2.249 2.222 1.175 89.5 178.6 90.5 1.820 1.802 1.676 1.156 1.160 1.174 165.7 113.2 2.336 1.778 1.827 1.160 1.151 85.2 1.560 1.806 1.807 1.158 1.156 82.2 1.758 1.776 1.515 1.156 1.157 80.2 151.7 101.5 1.774 1.772 1.152 1.159 1.762 1.758 1.481 1.155 1.157 99.4 1.964 1.037 115.0 1.790 2.203 1.962 1.190 1.034 1.035 118.6 176.7 90.6 115.6 113.9 1.785 1.154 1.943 1.943 1.153 1.158
DGauss 1.867 1.159 1.812 2.272 2.243 1.180 91.5 176.8 87.8 1.812 1.802 1.666 1.158 1.161 1.170 158.0 121.8 2.353 1.774 1.827 1.162 1.153 85.2 1.570 1.804 1.806 1.160 1.158 82.2 1.756 1.771 1.520 1.158 1.158 81.9 149.9 100.3 1.769 1.767 1.157 1.160 1.755 1.754 1.48 1 1.156 1.157 99.8 1.962 1.041 114.7 1.790 2.210 1.962 1.189 1.036 1.037 118.7 176.7 90.6 113.5 113.8 1.780 1.155 1.943 1.944 1.155 1.159
exptd 1.909' 1.137 1.817f 2.338 2.282 1.15 90.4 179.6 89.3 1.886 1.851 1.797 1.135 1.167 1.152 179.6 118.9 2.3678 1.807 1.892 1.122 1.108 88.3 1.576h 1.854 1.856 1.142 1.142 85.5 1.832 1.802 1.556 1.145 1.145 100.0 148.5 96.0 1.807 1.827 1.152 1.152 1.764 1.818 1.556 1.141 1.141 99.7 1.96' 1,871'3' 2.220 1.981 1.154 0.90 0.90 119.0 169.7 88.5 109.5 109.5 1.838 1.141 1.950' 1.969 1.143 1.143
'Bond distances in angstroms. Bond anales in degrees. Reference 13. Reference 11 Excerimental results from ref 11 unless noted. 'Reference 48 for experiGent. /Reference 49 for experiment. #Reference 50 for experiment. hReferencc 31 for experiment. 'Reference 51 for = 1.950 A experiment. 'Reference 52 for experiment. 'Because the CoNO bond is not linear, average values are given. For example, r(Cc-N,) for the NH, groups adjacent to the N O group and 1.974 A for the nonadjacent groups. 'Reference 42 for experiment.
Fock-Slater method predicted r(Co-H) = 1.52 A, R(Co-C,,) = 1.78 A, r(Co-C,) = 1.75 A, r(C-0,) = 1.14 A, and r(C-OaX) = 1.16 A.32 These values are similar to our LDF results except that our calculated Co-H bond length is shorter by about 0.04
A. Calculations at the SCF level do not give a good prediction of the structure as shown in Table VI (see also ref 33). For H,Fe(CO),, the tabulated SCF geometry is not as good as the LDF geometry except for the H-FeH bond angle. Another
Molecular Transition-Metal Compounds geometry optimization of this molecule at the SCF level and with slightly larger basis sets yields a value of 92.7O for the H-FE-H bond angle.34 However, as noted above, the significant error associated with the experimental value invites further study from the experimental and theoretical communities on the structure of this iron compound. Calculations on cr(co)635and ferrocene36have shown that it is important to treat the nonbonded interactions between the ligands properly. The ligands are close enough that a reasonable description of the van der Waals potential between the ligands is needed. It has been established that the LDF method leads to overbinding and that nonbonded interactions are too attractive: for example, there is not enough repulsion between the molecules in a hydrogen bond.16 It has been found that inclusion of nonlocal corrections to the exchange and correlation energies led to more repulsion between the ligands, Le., to a weaker attractive interaction in the hydrogen-bond case. We have thus reoptimized the geometries pointwise with nonlocal corrections for the metalcarbon bonds in Ni(C0)4 and Fe(CO)5 and for the Mn-N bond in Mn(C0)4N0. We have used the Becke correction3's3*to the exchange energy and the Perdew correction39to the correlation energy. The nonlocal corrections are applied after a self-consistent LDF energy has been calculated. The nonlocal calculations were all done with Gaussian basis sets. For Ni(C0)4 the Ni-C bond lengthens to 1.827 A, only 0.01 1 A shorter than experiment. The Ni-C bond length is in excellent agreement with molecular orbital results at the CCSD(T) level which give a bond length of 1.83 A.40
A similar result is found for Fe(CO)5where the two F e C bonds lengthen to 1.8 1 1 and 1.8 10 A for the equatorialand axial bonds, respectively. The axial bond is now in excellent agreement, but the equatorial bond is still too short. Thus, the calculations do not reproduce the difference in the axial and equatorial bond lengths. However, it has been reported that the low-temperature crystal structure of Fe(CO)5 does not show a difference in bond lengths4] The structure of Fe(CO)5 is not well-determined by Hartree-Fock The Fe-C,, bond is calculated to be too long (2.05 A) by 0.24 A. The Fe-C, bond is calculated to be too long (1-88 A) by 0.05 A. Correlation significantly shortens the Fe-C,, bond to between 1.78 and 1.80 A depending on the basis set and the level of correlation. At the same level, the F e C bond is calculated to be 1.836 A and the FeC, bond is calculatd to be 1.798 It is interesting to note that the agreement between the LDF calculations and experiment43for R U ( C O ) is ~ much better than for Fe(CO)5. The same improvements in the M-C bond lengths for metal carbonyls due to the inclusion of nonlocal effects was recently reported by Fan and Ziegler.44These authors used the same nonlocal corrections as us but solved for the energy by using the Hartree-Fock-Slater method of Baerends et al.45 For the Mn-N bond in Mn(CO),NO, a similar increase of 0.052 A due to nonlocal corrections is found giving an Mn-N bond length of 1.718 A. This bond is still too short by 0.08 A as compared to experiment. A careful examination of the reported crystal structure46shows that the structure is disordered for the three ligands in the plane. A disordered structure with a 2/3 C and N occupancy is essentially the same as the ordered structure with an N O on the 2-fold axis. Both structures give a long Mn-N bond. We then searched the Cambridge Crystal Structure Data Base4' for structures containing the Mn-NO moiety. For a reasonable representative set, the NO bonds range from about 1.64 to 1.70 A. Thus, it seems likely that the calculated structure is a good representation of the structure of the isolated molecule. The large range of errors found for the HF/3-21G normal-mode vibrational frequencies of transition-metal compounds discouraged the authors from recommending this level of theory for frequency assignments."S1* However, the LDF range of errors for the same set of frequencies is clearly smaller than observed for the HF/ 3-21G method. The calculated frequencies for the first-row diatomics are in good agreement with the observed harmonic values, usually within 50 cm-I. There are some differences between
The Journal of Physical Chemistry. Vol. 96, No. 16, 1992 6635 the numerical and Gaussian results, but these are not large. The LDF method seems to underestimate the harmonic stretching frequencies for a number of these species as the calculated values are harmonic and the experimental ones include anharmonic effects which are usually on the order of 5%. This is especially noticeable for the M-Cl and M-F stretches. The numerical basis set results are very sensitive to the orientation of the tetrahedral molecules and to the grid that is used. This is most important for the bending modes where the calculations underestimate the lowest bending mode. This can be seen in the results for TiC14and ZrC14where the e bends are too low. Although this is only a small error for these compounds, the problem is more pronounced for Ni(C0)4. The t2 and e bends are too small with the numerical basis set, and the e bend actually becomes imaginary. In order to only obtain one set of imaginary frequencies, the grid had to be made tighter (GRID parameter in DMol set to XFINE). Thus, the method with numerical basis sets is very sensitive to these low-frequency bends in a tetrahedral environment whereas the method with Gaussian basis sets is not so sensitive. Since the largest errors in the geometry for the compounds whose frequencies were calculated are found for Ni(C0)4, it is not surprising that the largest errors in the frequencies were found for this compound. Futhermore, the differences between the numerical and Gaussian basis sets are largest for this compound. The CO stretches are treated adequately by the LDF methods, but the NiC stretches are predicted to be too high consistent with the too short predicted Ni-C bond length. The high-frequency e bend is predicted to be too high also. The calculated frequencies for the second-row diatomics are in reasonable agreement with the experimental values. The largest error is found for AgH, which has the largest error in the bond distance and as described above has the largest relativistic correction. The calculated results for the polyatomics are in quite good agreement with experiment except for the bl and b2 infrared-inactive species of OMoCI,. These frequencies were not observed and were assigned to the values from OWC14. For the b2 bend this is reasonable, but for the two bl bends this transfer is not appropriate, especially for the bend which is predicted to have a very low frequency.
Conclusions The local density functional method provides a useful method to predict the structures and vibrational frequencies of a range of compounds containing transition metals. This is in contrast to the results found at the Hartree-Fock level where much larger errors in the predicted structures are found. The local density functional approach apparently incorporates the dominant correlation effects required for predicting these structures. The results suggest that the LDF method should be of use in predicting the structures of compounds that are important in catalytic systems. The use of the method is limited at present in the search for catalytic pathways by the lack of analytic second derivatives that have been incorporated in a computer program. This means that searching for transition states will be difficult as the second derivatives must be calculated numerically. Furthermore, the demonstration of local minima for large systems is also very expensive if this can only be done with numerical second derivatives. Finally, it is clear that nonlocal corrections play a small but important role in predicting the structure of compounds with dative bonds. Again, the ability to obtain analytical derivatives with nonlocal corrections will be an important advantage. Acknowledgment. We thank Nathalie Godbout for helpful discussions. We thank D. Thorn and A. J. Arduengo for help with the Cambridge Crystal package. Registry No. [Cu(CH,),]-, 47942-29-0; [Zn(CH&] *-,142423- 16-3.
References and Notes (1) Parr, R. G.; Yang, W. Density Functional Theory of
Atoms and Molecules; Oxford Press: New York, 1989 and references therein. ( 2 ) Salahub, D.R In Ab Initio Methods in Quantum Chemistry-II; LawIcy, K. P., Ed.; J. Wiley & Sons: New York, 1987; p 447.
6636
J. Phys. Chem. 1992,96, 6636-6640
(3) Wimmer, E.; Freeman, A. J.; Fu, C.-L.; Cao, P.-L.; Chou, S.-H.; Delley, B. In Supercomputer Research in Chemistry and Chemical Engineerins Jensen, K. F., Truhlar, D. G., Eds.; ACS Symposium Series; American Chemical Society: Washington, DC,1987; p 49 and references therein. (4) Jones, R. 0.;G u n n a m n , 0. Rev. Mod. Phys. 1989,61,689. (5) Andzelm, J.; Wimmer, E.; Salahub, D. R. In The Challenge of d and f Electrons: Theory and Computation; Salahub, D. R., Zerner, M. C., Eds.; ACS Symposium Series No. 394; American Chemical Society: Washington, DC, 1989; p 228 and references therein. (6) Ziegler, T.; Tschinke, V. In Density Functional Methods in Chemistry; Labanowski, J., A n h l m , J., Ed%; Springer-Verlag: New York, 1991; p 139. (7) Ziegler, T. Chem. Reo. 1991, 91, 651. (8) Dixon, D. A,; Andzelm, A.; Fitzgerald, Wimmer, E.; Delley, B. In Science and Engineering on Supercomputers; Pitcher, E. J., Ed.; Computational Science Publications: Southhampton, England, 1990; p 285. 19) Dixon. D. A.: Andzelm. J.: Fitznerald, G.; Wimmer. E.: Jasien. P. In &&ty FunciionaIMerhods inchemis&, L&nowslti, J., Andzelm, J.; Eds.; Springer-Verlag: New York, 1991; p 33. (LO) Dixon, D. A.; Andzelm, J.; Fitzgerald, G.; Wimmer, E. J . Phys. Chem. 1991,95,9197. (11) Dobbs, K. D.; Hehre, W. J. J. Comput. Chem. 1987,8, 861. (12) Dobbs, K. D.; Hehre, W. J. J. Comput. Chem. 1987,8, 880. (13) Pietro, W. J.; Hehre, W. J. J . Comput. Chem. 1983, 4, 241. (14) Delley, B. J . Chem. Phys. 1990, 92, 508. Dmol is available commercially from BIOSYM Technologies, San Diego, CA. (IS) Andzelm, J. In Density Functional Methods in Chemistry; Labanowski, J., Andzelm, J., Eds.;Springer-Verlag: New York, 1991; p 155. DGauss is available as part of the UniChem software environment from Cray Research, Eagan, MN. (16) Andzelm, J. W.; Wimmer, E. J . Chem. Phys. 1992, 96, 1280. (17) Hohenberg, P.; Kohn, W. Phys. Rev. B 1964, 136, 184. (18) Hedin, L.; Lundquist, B. I. J. Phys. C 1971, 4, 2064. (19) Ceperley, D. M.; Alder, B. J. Phys. Reu. Lett. 1980, 45, 566. (20) von Barth, U.; Hedin, L. J. Phys. C 1972, 5, 1629. (21) Vosko, S. H.; Wilk, L.; Nusair, M. Can. J . Phys. 1980, 58, 1200. (22) This grid can be obtained by using the FINE parameter in DMol. (23) Dunlap, B.; Connolly, J.; Sabin, J. J. Chem. Phys. 1979, 71, 3396. (24) Huzinaga, S.;Andzelm, J.; Klobukowski, M.; Radzio, E.; Sakai, Y.; Tatasaki, H. Gaussian Basis Sets of Molecular Calculations;Elsevier: Amsterdam, 1984. (25) (a) Andzelm, J.; Radzio, E.; Salahub, D. R. J. Comput. Chem. 1985, 6. 520. (b) Godbout.. N.:. Salahub. D.: Andzelm. J. W.: Wimmer, E. Can. J . Chem., h‘press. (26) (a) Foumier, R.; Andzelm, J.; Salahub, D. R. J . Chem. Phys. 1989, 90,6371. (b) Dunlap, B. I.; Andzelm, J.; Mintmire, J. W. Phys. Rev. A 1990,
42, 6354. (27) Becke, A. J. Chem. Phys. 1988, 88, 2547. (28) Obara, S.;Saika, A. J. Chem. Phys. 1986,84, 3963. (29) Delley, B. In &miry Functional Methods in Chemistry; Labanowski, J., Andzelm, J., Eds.; Springer-Verlag: New York, 1991; p 101. (30) Versluis, I.; Ziegler, T. J . Chem. Phys. 1988, 88, 3322. (31) McNeill, E. A.; Scholer, F. R. J. Am. Chem. Soc. 1977,99,6243. (32) Veduis, L.; Ziegler, T.; Baerends, E. J.; Ravenek, W. J. Am. Chem. Soc. 1989, I l l , 2018. (33) Antolovic, D.; Davidson, E. R. J . Am. Chem. Soc. 1987, 109, 977; J. Chem. Phys. 1988,88,4967. See also: Veillard, A,; Daniel, C.; Rohmer, M.-M. J. Phys. Chem. 1990, 94, 5556. (34) Dedieu, A.; Nakamura, S.;Sheldon, J. C. Chem. Phys. Lett. 1987, 141, 323. (35) (a) Davidson, E. J. Am. Chem. SOC.,submitted. (b) Davidson, E. Presented at North Carolina ACS Symposium, Sept 1990. (36) Park, C.; Almlbf, J. J. Chem. Phys. 1991, 95, 1829. (37) Becke, A. D. In The Challenge of d and f Electrons: Theory and Computation; Salahub, D. R., Zerner, M. C., Eds.;ACS Symposium Series No. 394; American Chemical Society: Washington, DC, 1989; p 166. (38) Becke, A. D. Int. J. Quantum Chem., Quuntum Chem. Symp. 1989, 23, 599. (39) Perdew, J. P. Phys. Rev. B 1986, 33, 8822. (40) Blomberg, M. R. A,; Siegbahn, P. E. M.; Lee, T. J.; Rendell, A. P.; Rice, J. E. J. Chem. Phys. 1991, 95, 5898. (41) Unpublished results quoted by: Scriver, D. F.; Whitmire. K. H.In Comprehensiue Organometallic Chemistry; Willtinson, G., Stone, F. G. A., Abel, E. W., Eds.; Pergamon: New York, 1982; Vol. 4, p 245. (42) (a) Lothi, H.P.; Siegbahn, P. E. M.; Almbf, J. J. Phys. Chem. 1985, 89, 2156. (b) Demuynck, J.; Strick, A.; Veillard, A. Nouu. J. Chim. 1977, I, 217. (43) Huang, J.; Hedberg, K.; Davis, H.B.;Pomeroy, R. K. Inorg. Chem. 1990,29,3923. (44)Fan, L.; Ziegler, T. J. Chem. Phys. 1991, 95, 7401. This work appeared after completion of our work. (45) Baerends, E. J.; Ellis, D. E.; Ros, P. Chem. Phys. 1973, 2, 41. (46)Frenz, B. A.; Enemark, J. H.; Ibers, J. A. Inorg. Chem. 1969,8,1288. (47) For Cambridge Crystal Structure Data Base, see: Allen, F. H.; Kennard, 0.;Taylor, R. Acc. Chem. Res. 1983, 16, 146-153. (48) Whitaker, A.; Jeffery, J. W. Acta Crysrallogr. 1967, 23, 977. (49) Huttner, G.; Schelle, S.J. Cryst. Mol. Struct. 1971, 1, 69. (50) Greene, P. T.; Bryan, R. F. J. Chem. SOC.A 1971, 1559. (51) Barnet, M. T.; Craven, B. M.; Freeman, H. C.; Kime, N. E.; Ibers, J. A. J. Chem. Soc., Chem. Commun. 1966, 307. (52) Pratt, C. S.;Coyle, B. A,; Ibers, J. A. J . Chem. SOC.A 1971, 2146.
Ab Initio Post-Hartree-Fock Study of Molecular Structures and Vibrational Spectra of Phosphine Oxide, Phosphinous AcM, and Their Thio Analogst Jiizef S. Kwiatkowski**% and Jerzy Leszczyiiski* Department of Chemistry, Jackson State University, Jackson, Mississippi 3921 7 (Received: February 28, 1992) The vibrational IR spectra (harmonic frequencies, integrated intensities) of phosphine oxide (H,PO), two conformers (cis, trans) of phosphinous acid (H,POH), and their thio analogs are predicted by ab initio secondsrder Mallel-Plesset perturbation theory (MP2) with the 6-31G** and 6-31 1G** basis sets. Also the molecular parameters of these species (bond lengths and bond angles, dipole moments, rotational constants) computed at the MP2/6-311G** level are presented. The calculated vibrational spectra compared well with the available experimental data, and the observed shifts of the IR wavenumbers upon isotopic substitution of the molecules are correctly predicted by the calculations. Reported data might be applied to identifkation of the sulfur species for which experimental data are still incomplete.
Introduction For a long time, information concerning the molecular structures and properties of simple phosphines H3PX (phosphine oxide (H3PO), phosphine sulfide (H3FS)) and their acidic forms (phopshinous acid (H2POH), phosphinothious acid (H2PSH)) have been available only from quantum mechanical calculations. ‘Part 4 in the series “The Quantitative Prediction and Interpretation of Vibrational Spectra of Organophosphorus Compounds”. For part 3 see ref 8.
f Permanent address: Insiytut Fizyki, Uniwersytet M. Kopemika, ul. Grudziadzka 5 , 87-100 Torufi, Poland.
0022-365419212096-6636$03.00/0
The main attention in these studies has been devoted to understand the character of the PO and PS bonds in these model phosphines with less interest directed toward prediction of their other parameters including vibrational IR spectra. Fortunately, the recent and infrared spectroscopic studies on phosphorus phosphorus sulfide^^.^ characterized vibrational IR spectra of several small phosphorus and sulfur species, providing the basic (reference) data for testing the reliability of different levels of a b initio quantum mechanical methods. However, there is still a lack of experimental geometrical parameters of H3PX and H2PXH (X = 0, S), though the vibrational IR spectra of both H3P0and H2POH species are well characterized.’ Unfortunately, 0 1992 American Chemical Society
