J . Phys. Chem. 1991,95,9811-9816 Applications (NCSA), Urbana, IL, under a research and development grant from Cray Research, Inc., Mendota Heights, Minnesota. Molecular mechanics calculations and data analysis were performed on an Ardent Titan at the Theoretical Chemistry Group, Argonne National Laboratory, Argonne, IL. We are grateful for the support of DowElanco, Cray Research, Inc.,
9811
Argonne National Laboratory, NCSA, and The Laboratory of Computer-Aided Molecular Modeling and Design a t UIC for support of this work. Registry No. H S 0 2 N H 2 , 25310-87-6; CH~SOZNHCHI,1184-85-6; CHIS02NH2, 3144-09-0; H S 0 2 N H C H , , 32797-25-4.
Maller-Plesset Perturbative ab Initio Calculations of the Ne, Potential Fu-Ming Taot and Yuh-Kang Pan* Department of Chemistry, Boston College, Chestnut Hill, Massachusetts 021 67 (Received: May 20, 1991)
The Msller-Plesset perturbation approach, up to the complete fourth order, is employed to calculate the energy potential curve of the neon dimer interaction (R = 2.0-6.0 A). A basis set of the energy-optimized 6 - 3 1 1 G set extended with the sp-, d-, and f-type diffuse and polarization functions is used, and the the full function counterpoise correction for the basis set superposition error is applied. HartreeFock interaction energies and correlation interaction energies from the fourth-order Mdler-Plesset perturbative calculations in the region of the van der Waals minimum are studied and compared to the available ab initio and other nonempirical results in the literature. The role of triple exctitations for the correlation interaction energy is stressed and confirmed by the calculations. Converged and nearly converged interaction energies are achieved at the distances (R = 3.2-6.0 A) beyond the van der Waals minimum, where about 90-100% of the total interaction energies are recovered by the complete fourth-order calculations. At shorter distances (R = 2.0-3.1 A), the interaction energies are still considerably underestimated even at the complete fourth-order level, due mainly to the lack of higher than f polarization functions in the present basis set, instead of the truncation of the Msller-Plesset series. This work demonstrates that the complete fourth-order Msller-Plesset perturbation may be used as a reliable method for reasonably accurate calculations of van der Waals potentials.
1. Introduction The study of van der Waals interactions has experienced a golden age of discovery in the past two decades.lJ Fast advances in experimental technologies, mostly by spectroscopic means,h5 and growing experimental data have greatly stimulated development of a theoretical understanding of the character of van der Waals interactions. Ab initio computation, as a major theoretical tool for molecular study, has been applied during the past 10 years to the study of van der Waals interactions among very small molecular specie^.^*^*' In contrast to its widespread applications in conventional molecular studies, however, ab initio studies of van der Waals molecules are still in a very limited stage. Special difficulties exist in the case of nonpolar-nonpolar pairs such as He2, Ne2, or (H2)2,whose interactions are predominanted by the dispersion energy. The effects of dispersion energy can only be fully accounted for by use of effective calculational methods with sufficiently large basis sets. Therefore, the search for effective methods and efficient basis sets is especially important at the present stage for the full establishment of ab initio method routinely used for study of van der Waals molecules. In this work, the Msller-Plesset perturbation approach through the fourth-order (MP2, MP3, MP4)8-'0 will be applied for the determination of the a b initio neon dimer potential. The study will focus primarily on the performance and convergence of the Molller-Plesset theory in calculations of van der Waals interaction energies, which has long been of particular interest to theoreticians as well as to experimental researchers." Several ab initio and other nonempirical calculations of the Ne-Ne interaction energies near the potential minimum have appeared in the literature. Chalasinski, van Lenthe, and GroenI2 analyzed the dispersion energy contributions for Ne2 by using partial wave expansions through the h-h term. Sauer, Hobza, Carsky, and Zahradnik') studied the performance of the MP2 on Ne2 with various extended basis sets. Krauss, Regan, and K o n o ~ a l o wdetermined '~ the Ne2 interaction energy curve by ab *To whom correspondence should be addressed. 'Current address: Department of Chemistry, Brown University, Providence, RI 0291 2.
0022-3654/91/209S-9811%02.50/0
initio Hartree-Fock and damped dispersion energy calculations. To the best of our knowledge, however, no complete a b initio computations at higher theoretical levels have been performed on the Ne2 interaction potential. For similar but smaller systems, such as H e H e or H2-He, the situation is changing. The MR-CI potential by van Lenthe et aI.I5 and the LM-2 potential by Liu and McLean,16 both for the He dimer, have reached near-complete agreement with those of the best experimental potential of Aziz et a1.I' In addition to the H e dimer potential, an accurate a b ( I ) Novick, S. E.; Leopold, K. R.; Klemperer, W. The Structures of Weakly Bound Complexes as Elucidated by Microwave and Infrared Spectroscopy. In Atomic and Molecular Clusters; Bernstein, E. R., Eds.; Elsevier: Amsterdam, 1990; p 359. (2) van der Waals interactions: Chem. Reu. 1988, 88. (3) Dyke, T. R. Top. Curr. Chem. 1984, 120,85. Legon, A. C.; Millen, D. J. Chem. Rev. 1986,86,635. (4) Fraser, A. T.; Nelson, D. D., Jr.; Charo, A.; Klemperer, W. J . Chem. Phys. 1985,82, 2535. (5) Dyke, T. R.; Howard, B. J.; Klemperer, W. J . Chem. Phys 1972,56, 2442; 1984,81, 5417. (6) Le Roy, R. J.; Carley, J. S. Adu. Chem. Phys. 1980, 42, 353. (7) van Lenthe, J. H.; van Duijneveldt-van der Rijdt. J. G. C. M.; van Duijneveldt, F. B. Adu. Chem. Phys. 1987, 69, 521. (8) Maller, C.; Plesset, M. S.Phys. Reu. 1934, 46, 618. (9) Pople, J. A. Faraday Discuss. Chem. Soc. 1982, 73, 7.
(IO) Pople, J. A.; Krishnan, R.; Schlegel, H. 8.; Binkley, J. S. In?. J .
Quantum Chem. 1978, 14, 545. Krishnan, R.; Frisch, M. J.; Pople, J. A. J . Chem. Phys. 1980, 72, 4244. Frisch, M. J.; Krishnan, R.; Pople, J. A. Chem. Phys. Lett 1980, 75, 60.
(1 1) Wilson, S. Diagrammatic Many-Body Perturbation Theory ofAtomic and Molecular Electronic Structure; North-Holland: Armsterdam, 1985. (12) Chalasinski, G.; van Lenthe, J. H.; Groen, Th. P. Chem. Phys. Lert. 1984, 110, 369. (13) Sauer, J.; Hobza, P.; Carsky, P.; Zahradnik, R. Chem. Phys. Lert. 1987, 134, 553. (14) Krauss, M.; Regan. R. M.; Konowalow, D. D. J . Phys. Chem. 1988, 92, 4329. ( 1 5) Vos, R. J.; van Lenthe, J. H.; van Duijneveldt, F. B. J . Chem. Phys. 1990, 93, 643. (16) Liu, B.; McLean, A. D. J . Chem. Phys. 1989, 91, 2348. (17) Aziz, R. A,; McCourt, F. R. W.; Wong. C. C. K. Mol. Phys. 1987, 61. 1487.
0 1 9 9 1 American Chemical Society
!MI2 The Journal of Physical Chemistry, Vol. 95, No. 24, 1991 initio potential for He-H2 has been determined by Meyer et a1.18 All these accurate ab initio potentials reported so far, however, were performed under extremely demanding conditions and with the use of approximate C1 methods instead of a perturbation approach for electron correlation calculations. These approximate methods, not being size consistent, have limited their wide popularity among most researchers in the area of van der Waals molecular studies. A number of experiments have determined the empirical Ne-Ne potential in various forms and parameters. One of the most important forms is the HFD-C, a modification of the original HFD form, which considers a potential to be a sum of the HartreeFock interaction energy and the correlation interaction (or dispersion) energy. Two Ne2 potentials of the HFD-C form, HFD-C1 and HFD-C2,'9'z' have been reported with the use of two different sets of parameters. The HFD-C1 potential has a well depth of t = 136.17 phartrees at the equilibrium distance Re = 3.060 A and a zero crossing point of u = 2.746 A. The HFD-C2 potential has the corresponding parameters of t = 133.80 phartrees Re = 3.087 A, and u = 2.758 A, respectively. The potential well depth t was reported to be smaller in other forms of the potential, i.e., 128.34 or 123.39 phartrees in the n - 6 form,22132.54, 130.67, or 133.54 phartrees in the XC mode1,2"z5 and 130.47 phartrees in the original HFD form.26 The experimental equilibrium distance was'also reported to be larger, Re = 3.102 A, by Farrar et al.27 The discrepancies of the well depth and of the equilibrium distance reflect the levels of accuracy of current experimental Ne, potentials. 11. Calculations A. Basis Set. In general, there are three problems that limit the ab initio approach in the determination of interaction potentials: ( I ) the size of the requisite basis set, (2) the demanding level of the theoretical method, and (3) the effects of basis set superposition error (or BSSE).2-8 Dealing with these three problems has become a major task in the ab initio determination of van der Waals interaction potentials. We consider that the size of requisite basis set is the first and most prominent problem among all the three. This is because first the supersized basis set can make even a very simple van der Waals molecule a huge problem to deal with, taking unusually long computer time and an extraordinarily large core memory for integrals. For example, for the He, system, Liu and M c h n l 6 used a primitive Slater-type basis of [9s, 7p, 4d,3f, 2g, 1 h] in the determination of the LM-2 potential. If one wishes to achieve the same accuracy for the Ne2 system as of the LM-2 potential, one would have to face a basis set that is much larger than Liu and McLean's set and should also include symmetry polarization functions higher than f-type. Not only are the integral codes that can handle these basis functions not easily accessible but also the dimension of the basis set soon becomes prohibitively large. The second reason is that solving the other two problems depends largely on the resolution of the first problem, i.e., the size of requisite basis set. Wells and WilsonzE3 indicated that as the size of the basis set systematically increases, the function counterpoise correction, which is a measure of the BSSE, increases and passes through a maximum before decreasing to its limiting value of zero. This means that if a complete or nearly complete basis set is used, there will be no need (18)Meyer, W.; Hariharan, P. C.; Kutzclnigg, W. J . Chem. Phys. 1980, 73, 1880. (19)Aziz, R. A. High Temp. High Press. 1980, 12,565. (20)Aziz, R. A.; Chen, H. H. J . Chem. Phys. 1977,67,5719. (21)Aziz, R. A. Mol. Phys. 1979, 38, 177. (22)Aziz, R. A. Chem. Phys. Leu. 1976, 40. 57. (23) Aziz, R.A.; Meath, W.J.; Allnatt, A. R. Chem. Phys. 1983.78.295. (24) Ng, K. C.; Meath, W. J.; Allnatt, A. R. Chem. Phys. 1978.32, 175. (25)Ng, K. C.; Meath, W. J.; Allnatt, A. R. Mol. Phys. 1979, 37,237. (26)Hepburn, J.; Scoles, G.; Penco, R. Chem. Phys. Lerr. 1975,36,451. Ahlrichs, R.; Penco, R.; Scoles, G. Chem. Phys. 1977, 19, 119. (27) Farrar, J. M.; et al. Chem. Phys. Lerr. 1973, 19, 359. (28)Wilson, S.Ab Inirio in Quantum Chemistry; Lawley, K. P., Ed.; Wiley: New York, 1987;p 439. (29)Wells, B. H.;Wilson, S.Mol. Phys. 1983, 50, 1295;1985, 54, 787.
Tao and Pan
extending function +(sp)
4d 3f
exwnent 0.108 163' 9.216. 2.304, 0.576, 0.144 2.5, 0.625, 0.156 25
aI.81" 0.108 163, 1.915962 0.144, 2.0 0.15625, 2.0
"See eq 1. 'The exponents of the other two sp functions (part of the 6-311G) are 0.397 06 and 1.457 56.
for any BSSE corrections. But this scheme to eliminate BSSE is far from being practical at the present time, and the BSSE corrections are still needed to get any meaningful calculated interaction energies. Among all the BSSE correction methods, the function counterpoise schemeMis the most widely used, but very controversial since it is believed by m a n ~ ~tol underestimate - ~ the real BSSE. However, from recent calculations we found that35 the overcorrection of the real BSSE by the counterpoise method vanishes rapidly with the increase of the basis set. Consequently, the BSSE can be virtually fully removed by the full counterpoise method if sufficiently large basis sets, which are practically accessible at present for most of the small van der Waals systems, are used. In summary, the basis set problem is the dominant factor in determining the accuracy of ab initio van der Waals potentials. A practical choice of basis set adopted by us is to take a given energy-optimized basis set such as 6-31 IG,% which is capable of giving a fairly good description of the core and valence orbitals near the nuclei, and to extend it with certain polarization functions and diffuse functions in order to improve the long-range polarizabilities and the tail regions of the orbitals. The 6-31 1G set is an especially good starting energy-optimized set for electron correlation calculations because, first, it already has the triple-{ quality valence orbitals for increasing flexibility of the representation for the outer valence region and, second, the orbital exponents and expansion coefficients are chosen by energy op timization at the second-order Mlaller-Pleset (MP2) perturbation level rather than at the Hartree-Fock level as usually done for other energy-optimized sets. The extended functions chosen by us are one set of sptype Gaussian functions with the exponent of 0.108 163 (called the sp diffuse function), four sets of d-type functions with the exponents of 9.216,2.304,0.576, and 0.144, and three sets of f-type functions with the exponents of 2.5,0.625, and 0.156 25. The larger of these Gaussian exponents are taken directly from the standard basis set 6-31 lG(2df3, and the rest are determined by the even-scaling rule; i.e., the exponents of the extended functions of the same symmetry are given as37
where NI is the number of basis functions of symmetry type 1. The al and ,!I1in eq 1 are constants determined by the existing exponents of the 6-31 lG(2df). Note that there is only one f function in the 6-31 1G(2df), and so the same value of & for the d functions is also used as Br for the f functions. Such an extended basis set is termed 6-31 lG+(4d3f) according to the conventions adopted by Pople. All the exponents of the extended functions in the M i set and the Corresponding al, values are summarized in Table I. It should be pointed out that only the lower exponent (30)Boys, S. F.;Bernardi, F. Mol. Phys. 1970. 19. 553. (31) Johansson, A.; Kollman, P.; Rothenberg, S. Theor. Chim. Acra (Berlin) 1973,29, 167. (32)Olivarcs del Valle, F. J.; Tolosa, S.;Ojalvo, E. A.; Espinosa, J. Chem. Phys. 1988, 127,343. (33)Schwenkc, D.W.; Truhlar, D. G. J . Chem. Phys. 1985.82, 2418; 1986,84,4113. (34) Collins, J. R.; Gallup, G. A. Chem. Phys. Leu. 1986,123,56; 1986, 129,329. (35)Tao, F. M.; Pan, Y. K. J . Phys. Chem. 1991, 95,3582. (36)Krishnan, R.; Frisch, M. J.: Pople, J. A. J . Chem. Phys. 1980, 72, 4244. (37)Ruedenberg, K.; Raffenetti. R. C.; Bardo, R. Energy Srrucrure and Reacrioiry; Proceedings of the Boulder Research Conference on Theoretical Chemistry, 1972;Wiley: New York, 1973.
The Journal of Physical Chemistry, Vol. 95, No. 24, 1991 9813
Calculation of the Ne2 Potential values are considered for the extended functions after the corresponding functions in the 6-31 lG(2df). This is primarily because orbitals with the upper exponent values, e.g., 36.864 for d and 10.0 for f functions, are so confined to the nuclei that they have virtually no contributions to the calculated interaction energies at the van der Waals distances. Upon finishing our computation with the 6-31 1G+(4d3f) basis set, we also tried to test for the saturation of the basis set (within the spdf quality) by adding one more of the even-tempered sp, d, or f functions with even lower exponents to the basis set in each of the additional calculations at several critical distances. The results, however, generally did not show any noticeable improvement on the potential curve except that an additional more diffuse sp function resulted in some improvement of the Hartree-Fock interaction energies at the long range. This means that the 6-31 1+G(4d3f) is close to the basis set saturation on the condition that higher than f-symmetry polarization functions are excluded from consideration. It should be pointed out that when the complete fourth-order M~ller-Plesset (MP4) calculations are carried out, the regular f function (with the exponent of 2.5) is dropped because of the limitation of the computer memory. The basis set for the MP4 calculations is termed 6-31 1+G(4d2f). B. Correlation Energy by the Mailer-Plesset Perturbation Approach. The total interaction energy AEINTis comprised of the Hartree-Fock interaction energy AEHFand the correlation interaction energy A S o R : Accurate calculations of AEHFis relatively easy and straightforward while accurate calculation of &OR necessitates the use of computational methods which adequately account for the effect of intra- and intermolecular correlations. In this work the hEcoR will be calculated by use of the second-, third-, and fourth-order Mdler-Plesset perturbation theory (MP2, MP3, MP4) as an approximate approach for electron correlation calculations. Therefore, the aE“R at the MP2, MP3, and MP4 approximation can be expressed respectively as ASoR(MP2) = A,??)
(3)
+
(4)
G o R ( M P 3 ) = AEcZ) AE(3)
+ @3) + U(4)
G o R ( M P 4 ) = AE(*)
(5) where AP2), AE(’), and AE(4)are contributions of the second-, third-, and fourth-order Mdler-Plesset perturbation energies, respectively. The M0ller-Plesset approach for electron correlation calculations is a preferred choice by us for its unique and very important feature in calculations of van der Waals interaction energies: Jt is size consistent.’* For other approximate CI methods, which are usually not size consistent, such as the correlated electron pair approximation (CEPA) and multireference CI (MR-CI), additional complicated procedures and arbitrary approximations must be employed solely to eliminate the sizeinconsistency error. Moreover, any correction procedure so far does not guarantee the complete removal of such an error.39 In this regard, the Moiler-Plesset perturbation approach is advantageous over other CI methods. The Morller-Plesset approach is particularly attractive also because it has recently become easily accessible (GAUSSIAN 82 and -86)@s4’ and its applications are very straightforward for anyone even those without specific theoretical backgrounds. (38) Pople, J. A.; Binkley, J. S.;Seeger, R. Inr. J. Quunrum Chem. 1976,
SIO. 1.
(39) Shavitt, 1.; Brown, F. B.; Burton, P. G. Int. 3. Quantum Chem. 1987, 31, 507. (40)Binkley, J. S.;Frisch, M. H.; DeFrccs, D. J.; Rahgavachari, K.; Whiteside. R. A.; Schlegal, H. B.; Fluder, E. M.; Pople, J. A. GAUSSIAN 82, Department of Chemistry, Carnegie-Mellon University, Pittsburgh, PA. (41) Frisch, M. J.; Binkley, J. S.;Schlegal, H. E.; Rahgavachari, K.; Melius. C. F.; Martin, R. L.; Steward, J. J. P.; Bobrowicz, F. W .Rohlfing C. M.; Kahn, L. R.; DeFras, D. J.; Seeger, R.; Whiteside, R. A.;’kox, D. J.: Fluder, E. M.;Pople, J. A. GAUSSIAN 86, Carnegie-Mellon Quantum Chemistry Publishing Unit, Pittsburgh, PA, 1984.
The major disadvantage of the M~rller-Plessetapproach that limits its applications to calculations of van der Waals potentials, as claimed by many a~thors,4*-~’ is the slow convergence of the theory, which causes an inadequate account for the electron correlation energies and therefore underestimations of interaction energies at low levels of approximation. In contrast to its wide applications in regular molecular calculations, the Moller-Plesset approach is still not so popularly used in accurate calculations of van der Waals interaction potentials. However, we consider that the slow convergence of the Moller-Plesset theory results from the use of a severely truncated basis set instead of the theory itself. This is because most of the convergence studies of the Merller-Plesset theory have employed basis sets too small to be used properly for calculations of van der Waals interactions. Not only do these small basis sets provide inadequate descriptions of van der Waals interactions and lead to serious underestimations of the interaction energies but they may also severely affect the convergence of the theory itself. Another factor that also possibly limits the performance of the Merller-Plesset approach is that in the MP2, MP3, and MP4SDQ (the fourth-order Merller-Plesset theory with single, double, and quadruple excitations) calculations the triply excited configurations are not included in the CI expansion, which results in the incomplete description of the influence of intramolecular correlation upon the intermolecular (dispersion energy) correlation, or negligence of the coupling between the intra- and intermolecular correlations. The configurations considered in the MP2 and MP3 are at most doubly excited, while in the MP4SDQ the quadruply excited configurations are also included but no triple excitations (on GAUSSIAN 82, MP4 defaults to MP4SDQM) are considered. The lowest level to include the triple excitations is the complete fourth-order Mdler-Plesset perturbation (MP4). Therefore, there is a very good chance to improve the performance of the Maller-Plesset approach at the fourth-order level in calculation of van der Waals potentials by including the triply excited configurations in the CI expansion. In this study, the complete MP4 calculations are performed with the use of a slightly smaller basis set, 6-31 1+G(4d2f), due to the limitation of computer memory as already pointed out in an earlier section. At the levels of MP2, MP3, and MP4SDQ, the role of the f function absent in the 6-311+G(4d2f) basis set has turned out to be insignificant in contributions (
