Open Submitted: 01 August 2009 Accepted: 10 November 2009 Published Online: 08 December 2009
J. Chem. Phys. 131, 224104 (2009); https://doi.org/10.1063/1.3269802
more...View Affiliations
View Contributors
  • Takeshi Sato
  • Hiromi Nakai

A new method to calculate the atom-atom dispersion coefficients in a molecule is proposed for the use in density functional theory with dispersion (DFT-D) correction. The method is based on the local response approximation due to Dobson and Dinte [Phys. Rev. Lett. 76, 1780 (1996)], with modified dielectric model recently proposed by Vydrov and van Voorhis [J. Chem. Phys. 130, 104105 (2009)]. The local response model is used to calculate the distributed multipole polarizabilities of atoms in a molecule, from which the dispersion coefficients are obtained by an explicit frequency integral of the Casimir–Polder type. Thus obtained atomic polarizabilities are also used in the damping function for the short-range singularity. Unlike empirical DFT-D methods, the local response dispersion (LRD) method is able to calculate the dispersion energy from the ground-state electron density only. It is applicable to any geometry, free from physical constants such as van der Waals radii or atomic polarizabilities, and computationally very efficient. The LRD method combined with the long-range corrected DFT functional (LC-BOP) is applied to calculations of S22 weakly bound complex set [Phys. Chem. Chem. Phys. 8, 1985 (2006)]. Binding energies obtained by the LC-BOP+LRD agree remarkably well with ab initio references.
The importance of weak interactions, such as hydrogen-bonding and van der Waals interactions, cannot be overemphasized in view of roles they play in diverse area of chemical physics, biochemistry, and material science.1,21. A. J. Stone, Theory of Intermolecular Forces (Clarendon, Oxford, 1996).2. J. Šponer, K. E. Riley, and P. Hobza, Phys. Chem. Chem. Phys. 10, 2595 (2008). https://doi.org/10.1039/b719370j Such interactions have different physical origins: the electrostatic, exchange-repulsion, charge-transfer, induction, and dispersion interactions, although all of them can ultimately be attributed to the Coulomb interaction between electrons and nuclei constituting a system.11. A. J. Stone, Theory of Intermolecular Forces (Clarendon, Oxford, 1996). Among them, the dispersion interaction is distinguished from others by ubiquitous and growing importance in larger systems.
Theoretical treatment of the dispersion interaction has faced difficulty because it is entirely an electron correlation effect. In ab initio wave function theories, correlated methods are necessary to include the dispersion energy. The second-order Møller–Plesset perturbation (MP2) method performs well for interactions between saturated molecules, although is problematic for aromatic interactions.3–53. V. Špirko, O. Engkvist, P. Soldán, H. L. Selzle, E. W. Schlag, and P. Hobza, J. Chem. Phys. 111, 572 (1999). https://doi.org/10.1063/1.4793384. S. Tsuzuki, K. Honda, T. Uchimaru, M. Mikami, and K. Tanabe, J. Am. Chem. Soc. 124, 104 (2002). https://doi.org/10.1021/ja01052125. M. O. Sinnokrot, E. F. Valeev, and C. D. Sherrill, J. Am. Chem. Soc. 124, 10887 (2002). https://doi.org/10.1021/ja025896h The coupled-cluster with single, double, and perturbative triple excitations [CCSD(T)] is established as a powerful tool to predict accurate interaction energies of weakly bound systems. The major drawbacks of these methods are the slow convergence with respect to the number of basis functions and unfavorable scaling of the computational cost against the system size. A method that can capture the essential part of the dispersion interaction with computational cost significantly lower than these accurate but demanding methods is required.
Density functional theory66. P. Hohenberg and W. Kohn, Phys. Rev. 136, B864 (1964). https://doi.org/10.1103/PhysRev.136.B864 (DFT) has gained attention as an ideally cost-effective method for calculations of large systems. In particular, the Kohn–Sham (KS) method77. W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965). https://doi.org/10.1103/PhysRev.140.A1133 has a distinct advantage that it can include the electron correlation essentially with the mean-field equation. The exchange-correlation (XC) functional is the key to such formal simplicity. The generalized gradient approximation (GGA),8–118. J. P. Perdew and W. Yue, Phys. Rev. B 33, 8800 (1986). https://doi.org/10.1103/PhysRevB.33.88009. A. D. Becke, Phys. Rev. A 38, 3098 (1988). https://doi.org/10.1103/PhysRevA.38.309810. C. Lee, W. Yang, and R. G. Parr, Phys. Rev. B 37, 785 (1988). https://doi.org/10.1103/PhysRevB.37.78511. J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). https://doi.org/10.1103/PhysRevLett.77.3865 meta-GGA,12–1912. E. I. Proynov, A. Vela, and D. R. Salahub, Chem. Phys. Lett. 230, 419 (1994). https://doi.org/10.1016/0009-2614(94)01189-313. A. D. Becke, J. Chem. Phys. 104, 1040 (1996). https://doi.org/10.1063/1.47082914. T. Van Voorhis and G. E. Scuseria, J. Chem. Phys. 109, 400 (1998). https://doi.org/10.1063/1.47657715. A. D. Becke, J. Chem. Phys. 109, 2092 (1998). https://doi.org/10.1063/1.47672216. M. Filatov and W. Thiel, Phys. Rev. A 57, 189 (1998). https://doi.org/10.1103/PhysRevA.57.18917. J. P. Perdew, S. Kurth, A. Zupan, and P. Blaha, Phys. Rev. Lett. 82, 2544 (1999). https://doi.org/10.1103/PhysRevLett.82.254418. A. D. Boese and N. C. Handy, J. Chem. Phys. 116, 9559 (2002). https://doi.org/10.1063/1.147630919. J. Tao, J. P. Perdew, V. N. Staroverov, and G. E. Scuseria, Phys. Rev. Lett. 91, 146401 (2003). https://doi.org/10.1103/PhysRevLett.91.146401 and the Hartree–Fock (HF)/DFT hybrid functionals20–2220. A. D. Becke, J. Chem. Phys. 98, 1372 (1993). https://doi.org/10.1063/1.46430421. A. D. Becke, J. Chem. Phys. 98, 5648 (1993). https://doi.org/10.1063/1.46491322. P. J. Stephens, F. J. Devlin, C. F. Chabalowski, and M. J. Frisch, J. Phys. Chem. 98, 11623 (1994). https://doi.org/10.1021/j100096a001 have been proposed as useful tools for broad chemical applications. However, these XC functionals have problems in describing weak van der Waals interactions.23,2423. S. Kristyán and P. Pulay, Chem. Phys. Lett. 229, 175 (1994). https://doi.org/10.1016/0009-2614(94)01027-724. J. M. Pérez-Jordá and A. D. Becke, Chem. Phys. Lett. 233, 134 (1995). https://doi.org/10.1016/0009-2614(94)01402-H The problem stems from the nearsightedness of both exchange and correlation functionals.
The problem of exchange functionals was first noted by Lacks and Gordon.2525. D. J. Lacks and R. G. Gordon, Phys. Rev. A 47, 4681 (1993). https://doi.org/10.1103/PhysRevA.47.4681 They applied local density approximation (LDA) and several GGA exchange functionals to calculations of helium and neon dimers, and revealed that the behavior of the exchange energy density at the low-density-high-gradient region had the crucial impact on calculated interaction energies. Wesolowski et al.2626. T. A. Wesolowski, O. Parisel, Y. Ellinger, and J. Weber, J. Phys. Chem. A 101, 7818 (1997). https://doi.org/10.1021/jp970586k and Zhang et al.2727. Y. Zhang, W. Pan, and W. Yang, J. Chem. Phys. 107, 7921 (1997). https://doi.org/10.1063/1.475105 derived the similar conclusion from calculations of other systems. These observations motivated several studies which aimed to adjust parameters for curing exchange functionals for such region.28–3028. C. Adamo and V. Barone, J. Chem. Phys. 108, 664 (1998). https://doi.org/10.1063/1.47542829. N. Kurita and H. Sekino, Chem. Phys. Lett. 348, 139 (2001). https://doi.org/10.1016/S0009-2614(01)01089-230. X. Xu and W. A. Goddard, Proc. Natl. Acad. Sci. U.S.A. 101, 2673 (2004). https://doi.org/10.1073/pnas.0308730100 Kamiya et al.3131. M. Kamiya, T. Tsuneda, and K. Hirao, J. Chem. Phys. 117, 6010 (2002). https://doi.org/10.1063/1.1501132 suggested that the long-range correction (LC) scheme,3232. H. Iikura, T. Tsuneda, T. Yanai, and K. Hirao, J. Chem. Phys. 115, 3540 (2001). https://doi.org/10.1063/1.1383587 which is the range separated HF/DFT hybrid method, can solve this problem by showing that the LC-applied functionals give closely resembling repulsive potential curves for rare-gas interactions irrespective of an employed exchange functional. Several authors also used LC (Refs. 33–4133. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/0026897041233133347434. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.201139635. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.274724336. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.307968437. O. A. Vydrov and T. van Voorhis, Phys. Rev. Lett. 103, 063004 (2009). https://doi.org/10.1103/PhysRevLett.103.06300438. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b39. J. Gräfenstein and D. Cremer, J. Chem. Phys. 130, 124105 (2009). https://doi.org/10.1063/1.307982240. B. G. Janesko, T. M. Henderson, and G. E. Scuseria, J. Chem. Phys. 130, 081105 (2009). https://doi.org/10.1063/1.309081441. B. G. Janesko, T. M. Henderson, and G. E. Scuseria, J. Chem. Phys. 131, 034110 (2009). https://doi.org/10.1063/1.3176514) or so-called range-separation42,4342. J. G. Ángyán, I. C. Gerber, A. Savin, and J. Toulouse, Phys. Rev. A 72, 012510 (2005). https://doi.org/10.1103/PhysRevA.72.01251043. J. Toulouse, I. C. Gerber, G. Jansen, A. Savin, and J. G. Ángyán, Phys. Rev. Lett. 102, 096404 (2009). https://doi.org/10.1103/PhysRevLett.102.096404 techniques with different treatments of the short-range exchange and the dispersion interactions.
Various approaches have also been suggested to treat the dispersion interaction within the DFT formalism. Here, we review some of these approaches restricting our focus on such methods that are essentially within the KS formalism and computationally not much more demanding than the KS method with conventional XC functionals. For comprehensive review of below-mentioned and other approaches, we guide readers to the recent paper by Gräfenstein and Cremer,3939. J. Gräfenstein and D. Cremer, J. Chem. Phys. 130, 124105 (2009). https://doi.org/10.1063/1.3079822 in which the authors categorized various DFT based approaches and analyzed the strengths and limitations of each method.
According to the philosophy of KS theory, the desired route may be improving accuracy of an XC functional to describe weak interactions. In this respect, recent progress in the (hybrid) meta-GGA functionals,44–4644. Y. Zhao, N. E. Schltz, and D. G. Truhlar, J. Chem. Theory Comput. 2, 364 (2006). https://doi.org/10.1021/ct050276345. Y. Zhao and D. G. Truhlar, Theor. Chem. Acc. 120, 215 (2008). https://doi.org/10.1007/s00214-007-0310-x46. Y. Zhang, A. Vela, and D. R. Salahub, Theor. Chem. Acc. 118, 693 (2007). https://doi.org/10.1007/s00214-007-0347-x especially those due to Zhao and Truhlar,44,4544. Y. Zhao, N. E. Schltz, and D. G. Truhlar, J. Chem. Theory Comput. 2, 364 (2006). https://doi.org/10.1021/ct050276345. Y. Zhao and D. G. Truhlar, Theor. Chem. Acc. 120, 215 (2008). https://doi.org/10.1007/s00214-007-0310-x have to be mentioned. Another recent progress has been made on the construction of dispersion-corrected atom-centered potentials by the group of Rothlisberger47–5047. O. A. von Lilienfeld, I. Tavernelli, U. Rothlisberger, and D. Sebastiani, Phys. Rev. Lett. 93, 153004 (2004). https://doi.org/10.1103/PhysRevLett.93.15300448. O. A. von Lilienfeld, I. Tavernelli, U. Rothlisberger, and D. Sebastiani, Phys. Rev. B 71, 195119 (2005). https://doi.org/10.1103/PhysRevB.71.19511949. I. -C. Lin, M. D. Coutinho-Neto, C. Felsenheimer, O. A. von Lilienfeld, I. Tavernelli, and U. Rothlisberger, Phys. Rev. B 75, 205131 (2007). https://doi.org/10.1103/PhysRevB.75.20513150. P. C. Aeberhard, J. S. Arey, I. -C. Lin, and U. Rothlisberger, J. Chem. Theory Comput. 5, 23 (2009). https://doi.org/10.1021/ct800299y and dispersion correcting potentials due to DiLabio and co-workers.51,5251. G. A. DiLabio, Chem. Phys. Lett. 455, 348 (2008). https://doi.org/10.1016/j.cplett.2008.02.11052. I. D. Mackie and G. A. DiLabio, J. Phys. Chem. A 112, 10968 (2008). https://doi.org/10.1021/jp806162t A practical advantage of these methods is that they can benefit from the sophisticated computational code for the KS calculation without large modifications; gradients and Hessians as well as various properties are readily available. The limitation of these methods is the inability to describe the long-range part of the dispersion interaction. (Here, the term “dispersion” is loosely used.) The physics of dispersion interaction, the synchronized fluctuation of electron distributions, cannot be described by the (semi-) local XC functionals or the atom-centered potentials.5353. J. F. Dobson, K. McLennan, A. Rubio, J. Wang, T. Gould, H. M. Le, and B. P. Dinte, Aust. J. Chem. 54, 513 (2001). https://doi.org/10.1071/CH01052
Probably the most widely used approach is explicitly adding damped atom-atom corrections of the form C6R6fdamp(R) to the DFT total energy (DFT-D),38,54–5838. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b54. X. Wu, M. C. Vargas, S. Nayak, V. Lotrich, and G. Scoles, J. Chem. Phys. 115, 8748 (2001). https://doi.org/10.1063/1.141200455. Q. Wu and W. Yang, J. Chem. Phys. 116, 515 (2002). https://doi.org/10.1063/1.142492856. S. Grimme, J. Comput. Chem. 25, 1463 (2004). https://doi.org/10.1002/jcc.2007857. S. Grimme, J. Comput. Chem. 27, 1787 (2006). https://doi.org/10.1002/jcc.2049558. P. Jurečka, J. Černý, P. Hobza, and D. R. Salahub, J. Comput. Chem. 28, 555 (2007). https://doi.org/10.1002/jcc.20570 with dispersion coefficients determined empirically. The DFT-D approaches have been quite successful in describing weak interactions among organic systems or carbonic materials.38,56–5838. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b56. S. Grimme, J. Comput. Chem. 25, 1463 (2004). https://doi.org/10.1002/jcc.2007857. S. Grimme, J. Comput. Chem. 27, 1787 (2006). https://doi.org/10.1002/jcc.2049558. P. Jurečka, J. Černý, P. Hobza, and D. R. Salahub, J. Comput. Chem. 28, 555 (2007). https://doi.org/10.1002/jcc.20570 However, the highly empirical nature of the method is unsatisfactory from a theoretical point of view; it might eventually suffer from the limited applicability and transferability. A promising alternative to these empirical approaches was proposed by Becke and Johnson,59,6059. A. D. Becke and E. R. Johnson, J. Chem. Phys. 122, 154104 (2005). https://doi.org/10.1063/1.188460160. A. D. Becke and E. R. Johnson, J. Chem. Phys. 123, 154101 (2005). https://doi.org/10.1063/1.2065267 in which the exchange-hole dipole moment (XDM) is used to model the dispersion coefficients. Subsequent studies by Geerlings and co-workers,61,6261. A. Olasz, K. Vanommeslaeghe, A. Krishtal, T. Veszprémi, C. V. Alsenoy, and P. Geerlings, J. Chem. Phys. 127, 224105 (2007). https://doi.org/10.1063/1.280539162. A. Krishtal, K. Vanommeslaeghe, A. Olasz, T. Veszprémi, C. V. Alsenoy, and P. Geerlings, J. Chem. Phys. 130, 174101 (2009). https://doi.org/10.1063/1.3126248 Kannemann and Becke,6363. F. O. Kannemann and A. D. Becke, J. Chem. Theory Comput. 5, 719 (2009). https://doi.org/10.1021/ct800522r and Kong et al.6464. J. Kong, Z. Gan, E. Proynov, M. Freindorf, and T. R. Furlani, Phys. Rev. A 79, 042510 (2009). https://doi.org/10.1103/PhysRevA.79.042510 showed that the XDM model can be a practically useful tool.
Yet another theoretically more sound approaches are based on the nonlocal dispersion functionals derived from the first principles.36,37,65–6836. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.307968437. O. A. Vydrov and T. van Voorhis, Phys. Rev. Lett. 103, 063004 (2009). https://doi.org/10.1103/PhysRevLett.103.06300465. K. Rapcewicz and N. W. Ashcroft, Phys. Rev. B 44, 4032 (1991). https://doi.org/10.1103/PhysRevB.44.403266. Y. Andersson, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 76, 102 (1996). https://doi.org/10.1103/PhysRevLett.76.10267. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.178068. M. Dion, H. Rydberg, E. Schröder, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 92, 246401 (2004). https://doi.org/10.1103/PhysRevLett.92.246401 The older Andersson–Langreth–Lundqvist (ALL) functional,6666. Y. Andersson, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 76, 102 (1996). https://doi.org/10.1103/PhysRevLett.76.102 which was also derived by Dobson and Dinte6767. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780 independently, is for the asymptotic interaction between nonoverlapping fragments, while the newer seamless van der Waals density functional (vdW-DF) (Ref. 6868. M. Dion, H. Rydberg, E. Schröder, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 92, 246401 (2004). https://doi.org/10.1103/PhysRevLett.92.246401) is applicable to any geometries. Kamiya et al.3131. M. Kamiya, T. Tsuneda, and K. Hirao, J. Chem. Phys. 117, 6010 (2002). https://doi.org/10.1063/1.1501132 and Sato et al.33–3533. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/0026897041233133347434. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.201139635. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243 successfully described various weak interactions by applying a damping function to the ALL functional combined with the long-range corrected XC functional. A similar approach with more efficient and rigorous implementation was recently reported by Gräfenstein and Cremer.3939. J. Gräfenstein and D. Cremer, J. Chem. Phys. 130, 124105 (2009). https://doi.org/10.1063/1.3079822 Silvestrelli69,7069. P. L. Silvestrelli, Phys. Rev. Lett. 100, 053002 (2008). https://doi.org/10.1103/PhysRevLett.100.05300270. P. L. Silvestrelli, J. Phys. Chem. A 113, 5224 (2009). https://doi.org/10.1021/jp811138n and co-workers7171. P. L. Silvestrelli, K. Benyahia, S. Grubisiĉ, F. Ancilotto, and F. Toigo, J. Chem. Phys. 130, 074702 (2009). https://doi.org/10.1063/1.3077288 cleverly simplified the ALL-based dispersion correction using the Wannier functions. On the other hand, the seamless vdW-DF (Ref. 6868. M. Dion, H. Rydberg, E. Schröder, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 92, 246401 (2004). https://doi.org/10.1103/PhysRevLett.92.246401) has been applied to both extended7272. S. D. Chakarova-Käck, E. Schröder, B. I. Lundqvist, and D. C. Langreth, Phys. Rev. Lett. 96, 146107 (2006). https://doi.org/10.1103/PhysRevLett.96.146107 and finite73,7473. A. Puzder, M. Dion, and D. C. Langreth, J. Chem. Phys. 124, 164105 (2006). https://doi.org/10.1063/1.218922974. V. R. Cooper, T. Thonhauser, and D. C. Langreth, J. Chem. Phys. 128, 204102 (2008). https://doi.org/10.1063/1.2924133 systems mainly with revPBE (Ref. 7575. Y. Zhang and W. Yang, Phys. Rev. Lett. 80, 890 (1998). https://doi.org/10.1103/PhysRevLett.80.890) XC functional and plane-wave basis. An implementation with Gaussian orbitals was first made by Vydrov et al.7676. O. A. Vydrov, Q. Wu, and T. van Voorhis, J. Chem. Phys. 129, 014106 (2008). https://doi.org/10.1063/1.2948400 Vydrov and van Voorhis3636. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.3079684 introduced a modified construction of the seamless functional (vdW-DF-09), achieving accuracy improvement with the LC-ωPBE functional.7777. O. A. Vydrov and G. E. Scuseria, J. Chem. Phys. 125, 234109 (2006). https://doi.org/10.1063/1.2409292 Recently, Vydrov and van Voorhis also derived a new analytical form of the seamless functional (VV09) (Ref. 3737. O. A. Vydrov and T. van Voorhis, Phys. Rev. Lett. 103, 063004 (2009). https://doi.org/10.1103/PhysRevLett.103.063004) based on the physically motivated dielectric model and simpler wavevector screening. A problem shared by nonlocal functional (NLF) approaches is the high computational cost for the explicit evaluation of the double spatial integral. See Refs. 78 and 7978. A. Gulans, M. J. Puska, and R. M. Nieminen, Phys. Rev. B 79, 201105 (2009). https://doi.org/10.1103/PhysRevB.79.20110579. G. Román-Pérez and J. M. Soler, Phys. Rev. Lett. 103, 096102 (2009). https://doi.org/10.1103/PhysRevLett.103.096102 for efficient implementations of vdW-DF (Ref. 6868. M. Dion, H. Rydberg, E. Schröder, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 92, 246401 (2004). https://doi.org/10.1103/PhysRevLett.92.246401) with pseudoatomic orbitals. A more serious limitation relating to the ALL-type asymptotic functionals31,33–35,3931. M. Kamiya, T. Tsuneda, and K. Hirao, J. Chem. Phys. 117, 6010 (2002). https://doi.org/10.1063/1.150113233. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/0026897041233133347434. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.201139635. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.274724339. J. Gräfenstein and D. Cremer, J. Chem. Phys. 130, 124105 (2009). https://doi.org/10.1063/1.3079822 is the lack of general applicability: It requires a priori fragmentation of a system to non- (or only weakly) overlapping subsystems.
In the present article, we propose an approach similar to the DFT-D methods, but with dispersion coefficients calculated nonempirically and directly in KS calculations. They are subsequently used for correcting DFT energy with an (exchange-corrected/dispersion-lacking) XC functional via the following form:
Edisp[ρ]=a>bn6Cnab[ρ]Rabnfdamp(n)(Rab),(1)
where a and b label atoms. To this end, the atom-atom dispersion coefficients in a molecule are derived based on the same theoretical model employed to derive the ALL-type NLF.66,6766. Y. Andersson, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 76, 102 (1996). https://doi.org/10.1103/PhysRevLett.76.10267. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780 One can expect that this approach complements the deficiencies of both empirical DFT-D and NLF approaches mentioned above, allowing efficient and accurate calculations with minimal empiricism. In subsequent sections, we will show that this is in fact the case. Our theory and implementation are detailed in Sec. II, its numerical applications are described in Sec. III, and the concluding remarks with prospects to the future work are given in Sec. IV.
Let us begin with the second-order dispersion energy expression between isolated (distinguishable) molecules A and B,
EdispAB=mAnB|mAnB|V̂AB|0A0B|2ωmA+ωnB,(2)
where V̂AB is the electrostatic interaction operator between molecules A and B, |0A and |mA are the ground and the mth excited states of the molecule A, and ωmA is the corresponding excitation energy, with similarly defined quantities for B. The atomic units are used throughout in this section. Applying the following integral transformation:
1a+b=2π0aa2+u2bb2+u2du,(3)
Eq. (2) is equivalently expressed as80–8280. H. C. Longuet-Higgins, Discuss. Faraday Soc. 40, 7 (1965). https://doi.org/10.1039/df965400000781. E. Zaremba and W. Kohn, Phys. Rev. B 13, 2270 (1976). https://doi.org/10.1103/PhysRevB.13.227082. R. McWeeny, Methods of Molecular Quantum Mechanics (Academic, London, 1992).
EdispAB=12πdr1dr1dr2dr2×0duχA(r1,r1,iu)χB(r2,r2,iu)|r1r2||r1r2|(4)
in terms of the dynamic density response functions χA and χB defined by
χA(r,r,ω)=2mAωmA(ωmA)2ω2×0A|ρ̂A(r)|mAmA|ρ̂A(r)|0A,(5)
with B analogs, where ρ̂A(r)=iAδ3(rir) is the density operator. The response function is a complicated quantity, of which calculations require solution of the time-dependent response equation. In order to avoid such a demanding computation, a simple approximation to the response function is required.
A. Review of local response model and its application to nonlocal functional
In 1996, Dobson and Dinte6767. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780 proposed a local approximation to the response function. In this approximation, the real-space density response is locally expressed in terms of the total electron density ρ as
χ(r,r,ω)=[ρ(r)ω02(r)ω2δ3(rr)],(6)
where ω0 is the long wavelength limit of the dispersion relation of the dielectric medium. Equation (6) was obtained by forcing the charge conservation and the reciprocity constraints on the response for the homogeneous electron density, which gives the long wavelength limit of the plasmon-pole approximated response in the Fourier space, χ(q0,ω)=q2ρ/(ω02ω2). By utilizing Eq. (6) in the response functions in Eq. (4) and reversely applying the Casimir–Polder transformation of Eq. (3), Dobson and Dinte6767. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780 derived the following NLF for the asymptotic dispersion energy:
EdispAB=332π2dr1dr2ω0A(r1)ω0B(r2)ω0A(r1)+ω0B(r2)1|r1r2|6.(7)
The resulting functional was identical to the ALL functional mentioned in Sec. I.6666. Y. Andersson, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 76, 102 (1996). https://doi.org/10.1103/PhysRevLett.76.102 Derivation by ALL (Ref. 6666. Y. Andersson, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 76, 102 (1996). https://doi.org/10.1103/PhysRevLett.76.102) was based on the approximation to the screened long-range interelectron interactions, which is equivalent to the approximation to the response function by Dobson and Dinte.6767. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780
In ALL and the original work of Dobson and Dinte,66,6766. Y. Andersson, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 76, 102 (1996). https://doi.org/10.1103/PhysRevLett.76.10267. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780 ω0 was taken to be the local-plasma frequency, ωP=4πρ. This formulation requires the real-space cutoff introduced by Rapcewicz and Ashcroft6565. K. Rapcewicz and N. W. Ashcroft, Phys. Rev. B 44, 4032 (1991). https://doi.org/10.1103/PhysRevB.44.4032 in the integration of the response function,
|ρ|6ρωPkF,(8)
where kF=(3π2ρ)1/3 is the local Fermi wave vector. The cutoff has been incorporated in the direct31,33–3531. M. Kamiya, T. Tsuneda, and K. Hirao, J. Chem. Phys. 117, 6010 (2002). https://doi.org/10.1063/1.150113233. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/0026897041233133347434. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.201139635. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243 or the smoothed39,69,7039. J. Gräfenstein and D. Cremer, J. Chem. Phys. 130, 124105 (2009). https://doi.org/10.1063/1.307982269. P. L. Silvestrelli, Phys. Rev. Lett. 100, 053002 (2008). https://doi.org/10.1103/PhysRevLett.100.05300270. P. L. Silvestrelli, J. Phys. Chem. A 113, 5224 (2009). https://doi.org/10.1021/jp811138n manner in practical calculations. However, from our experience,33–3533. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/0026897041233133347434. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.201139635. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243 results obtained with this cutoff are sensitive to the implementation at least with the direct approach.31,33–3531. M. Kamiya, T. Tsuneda, and K. Hirao, J. Chem. Phys. 117, 6010 (2002). https://doi.org/10.1063/1.150113233. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/0026897041233133347434. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.201139635. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243 This is because condition (8) has equality at the spatial points where the response is not small.
In the vdW-DF approach,6868. M. Dion, H. Rydberg, E. Schröder, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 92, 246401 (2004). https://doi.org/10.1103/PhysRevLett.92.246401 the above difficulty has been avoided by the use of a modified dispersion relation and the introduction of the gradient correction. The latter has been the key to the seamless theory.8383. D. C. Langreth, M. Dion, H. Rydberg, E. Schröder, P. Hyldgaard, and B. I. Lundqvist, Int. J. Quantum Chem. 101, 599 (2005). https://doi.org/10.1002/qua.20315 In the present study, we focus on the simplified dispersion relation recently proposed by Vydrov and van Voorhis,3636. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.3079684 of which the long wavelength limit reads
ω0(r)=q02(r)3,q0(r)=kF(1+λs2),(9)
where s=|ρ|/(2kFρ) is the reduced density gradient and λ is an empirical parameter introduced for adjusting decay of the response at the density tail.
B. Derivation of atom-atom expression of dispersion energy
A key to the atom-atom expression of the dispersion energy is a multicenter multipole expansion of the intermolecular Coulomb operator.84–8784. A. J. Stone and C. -S. Tong, Chem. Phys. 137, 121 (1989). https://doi.org/10.1016/0301-0104(89)87098-385. G. J. Williams and A. J. Stone, J. Chem. Phys. 119, 4620 (2003). https://doi.org/10.1063/1.159472286. C. Hättig, G. Jansen, B. A. Hess, and J. G. Ángyán, Mol. Phys. 91, 145 (1997). https://doi.org/10.1080/00268979717184187. J. G. Ángyán, J. Chem. Phys. 127, 024108 (2007). https://doi.org/10.1063/1.2749512 A molecular volume has to be divided or distributed into constituent atoms somehow. For this purpose, we adopt the atomic partition function used in XC numerical quadrature.8888. A. D. Becke, J. Chem. Phys. 88, 2547 (1988). https://doi.org/10.1063/1.454033 Thus, Coulomb interaction between electrons at r1 in molecule A and r2 in molecule B reads
1|r1r2|=aAbBwa(r1)wb(r2)|r1r2|=aAbBwa(r1)wb(r2)|Rabra1+rb2|,(10)
where we insert the normalization conditions aAwa(r1)=1 (also for b/B) of Becke-type function,8888. A. D. Becke, J. Chem. Phys. 88, 2547 (1988). https://doi.org/10.1063/1.454033 and rewrite vectors using Rab=RbRa, ra1=r1Ra, and rb2=r2Rb, with Ra and Rb being position vectors of atoms a and b. Here and in the following, labels {a, a} and {b, b} are used for atoms belonging to molecules A and B, respectively. All vectors, multipole moments, and polarizabilities are in the space-fixed global coordinate system.
Being scaled by the weight functions inside summations of Eq. (10), the Coulomb operator is reasonably expanded on its centers a and b. Such an expansion in terms of spherical harmonics is a well-established matter. We follow Stone’s book11. A. J. Stone, Theory of Intermolecular Forces (Clarendon, Oxford, 1996). within the required generality for our purpose. We have
1|Rabra1+rb2|=l1=0l2=0Rabl1l21m1m2Sm1m2ab(l1l2)×Rm1(l1)(ra1)Rm2(l2)(rb2),(11)
where Rab is the internuclear distance, Rm(l) is the regular solid harmonics, and the angular factor Sm1m2ab(l1l2) is defined by
Sm1m2ab(l1l2)=(1)l1[(2l1+2l2+1)!(2l1)!(2l2)!]1/2δL,l1+l2×(l1l2Lm1m2M)CM(L)(R̂ab),(12)
where the 2×3 matrix is Wigner’s 3-j symbol and CM(L) is the spherical harmonics with Racah’s normalization factor evaluated at the polar angles of the unit vector R̂ab=Rab/Rab.
Using Eqs. (10)–(12) for intermolecular Coulomb operators in Eq. (4), the dispersion energy between molecules A and B is expanded in terms of multipole moments centered on atoms {a} in A and {b} in B,
EdispAB=12πaaAbbBl1l1l2l2Rabl1l21Rabl1l21×m1m1m2m2Sm1m2ab(l1l2)Sm1m2ab(l1l2)×0αm1m1aa(l1l1)(iu)αm2m2bb(l2l2)(iu)du.(13)
Here, αmmaa(ll) is given by
αmmaa(ll)(iu)=drdrwa(r)wa(r)χA(r,r,iu)×Rm(l)(rRa)Rm(l)(rRa),(14)
which is interpreted as distributed polarizability. Equation (13) has the nonlocal dependence82,8482. R. McWeeny, Methods of Molecular Quantum Mechanics (Academic, London, 1992).84. A. J. Stone and C. -S. Tong, Chem. Phys. 137, 121 (1989). https://doi.org/10.1016/0301-0104(89)87098-3 on the atom-atom distances involving Rab and Rab. This is because the second-order energy involves two intermolecular operators. In order to obtain the convenient additive atom-atom expression of Eq. (1), the polarizabilities have to be localized. Here, this is achieved by relying on the local response model discussed in Sec. II A. Substituting Eq. (6) into Eq. (14), we have
αmmaa(ll)(iu)=drdrwa(r)wa(r)×[ρ(r)ω02(r)+u2δ3(rr)]×Rm(l)(rRa)Rm(l)(rRa)=drwa(r)wa(r)ρ(r)ω02(r)+u2×Rm(l)(rRa)Rm(l)(rRa).(15)
In deriving the second equation, we integrate by parts on each space variable and then use the delta function, essentially following the derivation of Eq. (7) by Dobson and Dinte.6767. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780 The exponential decay of the polarization density ρ/(ω02+u2) at the density tail3636. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.3079684 is also utilized. Note that the weight function derivative terms are ignored, discarding some contributions from the overlapping region of the weight functions on different centers. If we further assume wa(r)wa(r)δaawa2(r), different weight functions do not overlap significantly and/or contributions from such region are small, the off-diagonal blocks of polarizabilities can be omitted as αaa=δaaαa with
αmma(ll)(iu)=drwa2(r)ρ(r)ω02(r)+u2×Rm(l)(rRa)Rm(l)(rRa).(16)
Equation (13) is then simplified to the desired atom-atom dispersion expression,
EdispAB=aAbBn6CnabRabn,(17)
where the nth order interaction consists of contributions from various multipole-multipole interactions,
Cnab=l1l1l2l2Cab(l1l1,l2l2),(18)
with l1+l1+l2+l2+2=n. Finally, Cab(l1l1,l2l2) is given in terms of imaginary frequency integrals of the product of atomic polarizabilities,
Cab(l1l1,l2l2)=12πm1m1m2m2Sm1m2ab(l1l2)Sm1m2ab(l1l2)×0αm1m1a(l1l1)(iu)αm2m2b(l2l2)(iu)du.(19)
At this moment, we note the crucial difference between Eqs. (1) and (17) (apart from the damping function). The latter is for the dispersion energy between isolated two molecules. It does not apply to intramolecular interactions from the beginning of derivation: Eq. (2). Nevertheless, we dare to extend it to intramolecular interactions by eliminating the restriction for atomic labels in Eq. (17) and using the total (rather than monomer) weight functions and the electron density in Eq. (16). This is inspired by the success of empirical DFT-D approaches,38,55–5838. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b55. Q. Wu and W. Yang, J. Chem. Phys. 116, 515 (2002). https://doi.org/10.1063/1.142492856. S. Grimme, J. Comput. Chem. 25, 1463 (2004). https://doi.org/10.1002/jcc.2007857. S. Grimme, J. Comput. Chem. 27, 1787 (2006). https://doi.org/10.1002/jcc.2049558. P. Jurečka, J. Černý, P. Hobza, and D. R. Salahub, J. Comput. Chem. 28, 555 (2007). https://doi.org/10.1002/jcc.20570 where such an intramolecular expression is intrinsically assumed. It seems possible to justify this extension for a well-separated pair of atoms for which the antisymmetricity of the wave function gets irrelevant. For a pair of nearby atoms, such as those directly connected by a chemical bond, the treatment is clearly inappropriate. We have to rely on the damping function discussed in Sec. II C.
C. Damping function
The most widely used damping function in DFT-D approaches is the Fermi-type function first proposed by Wu and Yang,5454. X. Wu, M. C. Vargas, S. Nayak, V. Lotrich, and G. Scoles, J. Chem. Phys. 115, 8748 (2001). https://doi.org/10.1063/1.1412004
fdamp(Rab)=11+exp[d(Rab/R¯1)],(20)
where parameters d and R¯ determine the shape and position of the switching region, respectively. The function converges to zero quickly at the short distance. In spite of success in many DFT-D studies,55–5855. Q. Wu and W. Yang, J. Chem. Phys. 116, 515 (2002). https://doi.org/10.1063/1.142492856. S. Grimme, J. Comput. Chem. 25, 1463 (2004). https://doi.org/10.1002/jcc.2007857. S. Grimme, J. Comput. Chem. 27, 1787 (2006). https://doi.org/10.1002/jcc.2049558. P. Jurečka, J. Černý, P. Hobza, and D. R. Salahub, J. Comput. Chem. 28, 555 (2007). https://doi.org/10.1002/jcc.20570 we found in our preliminary test that this function was not suited for damping higher-order interactions beyond C6/R6. This is because the function is too steep with a reasonable d value (d=2025; necessary for a good behavior of C6 terms). This is shown in Fig. 1, in which the Fermi-type function as well as other two types of functions introduced shortly (optimized for neon dimer interaction) are compared. On the other hand, the older HF-D approaches89–9289. J. Hepburn and G. Scoles, Chem. Phys. Lett. 36, 451 (1975). https://doi.org/10.1016/0009-2614(75)80278-890. R. Ahlrichs, R. Penco, and G. Scoles, Chem. Phys. 19, 119 (1977). https://doi.org/10.1016/0301-0104(77)85124-091. K. T. Tang and J. P. Toennies, J. Chem. Phys. 66, 1496 (1977). https://doi.org/10.1063/1.43411392. K. T. Tang and J. P. Toennies, J. Chem. Phys. 68, 5501 (1978). https://doi.org/10.1063/1.435678 often incorporated higher-order terms for intermolecular interactions and various damping functions have been proposed. Among them is the one proposed by Tang and Toennies,9393. K. T. Tang and J. P. Toennies, J. Chem. Phys. 80, 3726 (1984). https://doi.org/10.1063/1.447150
fdamp(n)(Rab)=1[k=0n1k!(RabR¯)k]exp(RabR¯).(21)
In our test calculations, the Tang and Toennies (TT) function was found to work well for rare-gas interactions. However, in contrast to the Fermi-type function, the TT function decays too slowly at the short distance (see Fig. 1), hence is inappropriate for damping intramolecular interactions. In order to supplement the above two types of damping functions, we make an arbitrary choice of the following exponential form:
fdamp(n)(Rab)=exp[m(RabR¯)6],n=2m+4,(22)
which is, with n=6, identical to the function proposed by Kamiya et al.3131. M. Kamiya, T. Tsuneda, and K. Hirao, J. Chem. Phys. 117, 6010 (2002). https://doi.org/10.1063/1.1501132 As shown in Fig. 1, the new damping function varies slowly at long distance and vanishes quickly inside the critical point. We notice that the function proposed by Chai and Head-Gordon3838. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b shows a similar behavior as the present damping function.
For the effective radii, R¯, predetermined free atomic radii are often used. However, the damping function should take into account the effects of molecular environment. Tkatchenko and Schffler9494. A. Tkatchenko and M. Scheffler, Phys. Rev. Lett. 102, 073005 (2009). https://doi.org/10.1103/PhysRevLett.102.073005 realized this by noting the relation between polarizability and volume9595. T. Brinck, J. S. Murray, and P. Politzer, J. Chem. Phys. 98, 4305 (1993). https://doi.org/10.1063/1.465038 and employing the Hirshfeld partitioning of the electron density for the latter. A similar idea was also proposed by Johnson and Becke.9696. E. R. Johnson and A. D. Becke, J. Chem. Phys. 124, 174104 (2006). https://doi.org/10.1063/1.2190220 For the same purpose, we use the information from the local response of Sec. II B,
R¯=κ[(αaeff)1/3+(αbeff)1/3]+R0,(23)
where αaeff=1/3Tr[αa(11)(0)] is from Eq. (16) evaluated at zero frequency, which should correspond to the average dipole polarizability of an atom in a molecule. κ and R0 are global constants which should be determined empirically. These parameters are discussed in Sec. III.
We call the method described thus far the local response dispersion (LRD) method. It is defined by the use of Eq. (1) and followings:
(1)
local response approximation of Eq. (6) with the dielectric model of Eq. (9);
(2)
dispersion coefficients determined by Eqs. (16), (18), and (19) using Becke-type atomic partition function;
(3)
damping by Eqs. (22) and (23).
D. Implementation
Now we briefly describe our implementation of the LRD method. In the present study, we restrict atomic polarizabilities of Eq. (16) to those having two identical multipole rank, αmma(ll). Thus Eqs. (18) and (19) are brought together as
Cnab=l1l2m1m2m1m2Sm1m2ab(l1l2)Im1m2,m1m2ab(l1l2)Sm1m2ab(l1l2),(24)
where l1+l2+1=n/2 and the Casimir–Polder-type integral Im1m2,m1m2ab(l1l2) is given by
Im1m2,m1m2ab(l1l2)=12π0αm1m1a(l1l1)(iu)αm2m2b(l2l2)(iu)du.(25)
The restriction is reasonable since the terms in Eq. (16) make by far the dominant contribution. Due to this simplification, the number of polarizability elements for each atom is significantly reduced.
To facilitate the evaluation of Eq. (25) we change the variable to t=u/1+u2 and note that the integrand is an even function, then
Im1m2,m1m2ab(l1l2)=14π11dt1t2α¯m1m1a(l1)(t)α¯m2m2b(l2)(t),(26)
and
α¯mma(l)(t)(1t2)1/2αmma(ll)(iu)=drwa2(r)(1t2)1/2ρ(r)(1t2)ω02(r)+t2×Rm(l)(rRa)Rm(l)(rRa),(27)
which is calculated by means of a usual single-center integration technique used in the evaluation of XC functional. Explicitly applying the Gauss–Chebyshev rule to Eq. (26) gives
Im1m2,m1m2ab(l1l2)=14Nk=1Nα¯m1m1a(l1)(tk)α¯m2m2b(l2)(tk),(28)
with N quadrature points, tk=cos[(2k1)/4Nπ].
In the present study, the LRD method is implemented in a post-KS manner. The computational cost for the method is expected to be similar to that for a single XC energy calculation since the atomic polarizabilities of Eq. (27) at N required frequencies are calculated simultaneously by (partition function-aided) single-center integrations as in (grid-based) XC energy calculations. The evaluation of the integrand of Eq. (27), given the density and density gradient, is an O(n0) procedure at each spatial grid point, where n measures the system size, and inexpensive relative to the formally O(n2) density evaluation. Importantly, the density formation can be shared with XC energy and potential matrix evaluation in realizing a self-consistent field (SCF) implementation of the LRD method. Costs for Eqs. (1), (24), and (28) are negligibly small. In consequence, the LRD method is capable of calculating the dispersion energy with small additional cost to a usual KS calculation.
E. Comparison with similar approaches
At the end of this section, we compare the LRD method with other similar approaches for the dispersion problem of DFT. First, the LRD method shares much philosophy with the XDM model of Becke and Johnson,59,60,96–9859. A. D. Becke and E. R. Johnson, J. Chem. Phys. 122, 154104 (2005). https://doi.org/10.1063/1.188460160. A. D. Becke and E. R. Johnson, J. Chem. Phys. 123, 154101 (2005). https://doi.org/10.1063/1.206526796. E. R. Johnson and A. D. Becke, J. Chem. Phys. 124, 174104 (2006). https://doi.org/10.1063/1.219022097. A. D. Becke and E. R. Johnson, J. Chem. Phys. 124, 014104 (2006). https://doi.org/10.1063/1.213966898. A. D. Becke and E. R. Johnson, J. Chem. Phys. 127, 154108 (2007). https://doi.org/10.1063/1.2795701 which also obtains the system-dependent dispersion coefficients from the electron density6060. A. D. Becke and E. R. Johnson, J. Chem. Phys. 123, 154101 (2005). https://doi.org/10.1063/1.2065267 or from the occupied orbitals.5959. A. D. Becke and E. R. Johnson, J. Chem. Phys. 122, 154104 (2005). https://doi.org/10.1063/1.1884601 Theoretical investigations of the XDM model by Ángyán8787. J. G. Ángyán, J. Chem. Phys. 127, 024108 (2007). https://doi.org/10.1063/1.2749512 and Heßelmann9999. A. Heßelmann, J. Chem. Phys. 130, 084104 (2009). https://doi.org/10.1063/1.3077939 clarified the closely relating physics underlying the local response6767. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780 and XDM (Ref. 5959. A. D. Becke and E. R. Johnson, J. Chem. Phys. 122, 154104 (2005). https://doi.org/10.1063/1.1884601) models.
Next, we have to be back to the original work by ALL,6666. Y. Andersson, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 76, 102 (1996). https://doi.org/10.1103/PhysRevLett.76.102 in which the authors already used the local response model to compute the average C6 coefficients between numbers of molecules. Recently Vydrov and van Voorhis3636. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.3079684 also reported C6 values based on exactly the same model employed in the present work, Eq. (9). In this sense, the LRD method is only an extension of the previous work to (i) general multipole interactions (ii) between atoms in a molecule. However, these extensions are crucial for the method to be practically useful, as shown in Sec. III.
It is worth mentioning that the imaginary frequency integrals in Eq. (19) are easily avoidable by a reverse use of the Casimir–Polder integral transformation, Eq. (3), just as in deriving ALL functional.66,6766. Y. Andersson, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 76, 102 (1996). https://doi.org/10.1103/PhysRevLett.76.10267. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780 The average dipole-dipole interaction, for example, is then given by
C6ab=32dr1dr2wa2(r1)wb2(r2)×ρ(r1)ρ(r2)ω0(r1)ω0(r2)[ω0(r1)+ω0(r2)],(29)
which is essentially what Vydrov and van Voorhis3636. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.3079684 noted for intermolecular interactions. The proposition of the current work is the calculation of multipole polarizabilities for each atom followed by the explicit Casimir–Polder integration, which is more efficient than the double spatial integrals such as Eq. (29) while not sacrificing the accuracy. Silvestrelli and co-workers69,7069. P. L. Silvestrelli, Phys. Rev. Lett. 100, 053002 (2008). https://doi.org/10.1103/PhysRevLett.100.05300270. P. L. Silvestrelli, J. Phys. Chem. A 113, 5224 (2009). https://doi.org/10.1021/jp811138n took a different route. Their use of the maximally localized Wannier function formalism100100. N. Marzari and D. Vanderbilt, Phys. Rev. B 56, 12847 (1997). https://doi.org/10.1103/PhysRevB.56.12847 allows analytical representation of the localized polarization density, reducing the computational cost for the analogs of Eq. (29) without reintroducing the frequency integral.
In this section, the LRD method described thus far is numerically applied to calculations of several weakly bound systems. Although the LRD method can be combined with any XC functionals, the LC-BOP (Ref. 3232. H. Iikura, T. Tsuneda, T. Yanai, and K. Hirao, J. Chem. Phys. 115, 3540 (2001). https://doi.org/10.1063/1.1383587) functional [short-range Becke-88 (Ref. 99. A. D. Becke, Phys. Rev. A 38, 3098 (1988). https://doi.org/10.1103/PhysRevA.38.3098) plus long-range HF exchanges with one-parameter progressive (OP) correlation101101. T. Tsuneda, T. Suzumura, and K. Hirao, J. Chem. Phys. 110, 10664 (1999). https://doi.org/10.1063/1.479012 functional] is employed throughout this study (LC-BOP+LRD). The range-separation parameter of μ=0.47a.u. in the LC scheme is taken from Ref. 102102. J. -W. Song, T. Hirosawa, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 154105 (2007). https://doi.org/10.1063/1.2721532. The NLF of Eq. (7) is also examined (LC-BOP+NLF) for comparison. The NLF is numerically calculated as described in Ref. 3535. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243 using the gridwise damping function with cutoff parameters optimized as C1=0.4294 and C2=3.226a.u. following the procedure detailed therein.3535. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243
All calculations are performed with the modified version of GAUSSIAN 03 program package,103103. M. J. Frisch, G. W. Trucks, H. B. Schlegel et al., GAUSSIAN03, Revision D.02, Gaussian, Inc., Pittsburgh, PA, 2003. in which the LC scheme, the OP correlation functional, the LRD method, and NLF of Eq. (7) are locally implemented. The KS-SCF calculations are converged with the “tight” option.103103. M. J. Frisch, G. W. Trucks, H. B. Schlegel et al., GAUSSIAN03, Revision D.02, Gaussian, Inc., Pittsburgh, PA, 2003. The reduced 99×590 grids104,105104. P. M. W. Gill, B. G. Johnson, and J. A. Pople, Chem. Phys. Lett. 209, 506 (1993). https://doi.org/10.1016/0009-2614(93)80125-9105. V. I. Lebedev, Dokl. Math. 45, 587 (1992). are used unless otherwise noted for spatial integrals in XC, LRD, and NLF computations with Becke’s atomic partition function.8888. A. D. Becke, J. Chem. Phys. 88, 2547 (1988). https://doi.org/10.1063/1.454033 All the reported interaction energies are corrected for the basis set superposition error by counterpoise method.106106. S. F. Boys and F. Bernardi, Mol. Phys. 19, 553 (1970). https://doi.org/10.1080/00268977000101561 Further details of computations are presented separately.
A. Rare-gas diatoms: Parameter fitting
Before proceeding to calculations of real interests, several parameters left unfixed are determined in this subsection based on calculations of rare-gas diatoms: He–He, He–Ne, He–Ar, Ne–Ne, Ne–Ar, and Ar–Ar. The rare-gas diatoms are ideal systems for setting up the LRD method since the only source of attractions is the dispersion interaction. The aug-cc-pVQZ basis set is used in this subsection.
Following Vydrov and van Voorhis,3636. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.3079684 the parameter λ was fixed to give accurate C6 coefficients for the six rare-gas diatoms. To avoid errors due to finite N and weight functions in Eq. (27), the average C6 coefficients were calculated using Eq. (29) with a large (100 Å) internuclear distance. The smallest mean absolute percent error of 3.7% relative to the experimental values was obtained with λ=0.232a.u. This is near λ=0.22a.u. obtained using the LC-ωPBE functional and a larger benchmark set.3636. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.3079684 The C6 coefficients calculated with this value of λ are listed in Table I.
Table icon
Table I. C6 dispersion coefficients (a.u.) of rare-gas diatoms calculated using Eq. (29) (using the LC-BOP/aug-cc-pVQZ SCF density). Experimental results (from Ref. 107107. A. Kumar and W. J. Meath, Mol. Phys. 54, 823 (1985). https://doi.org/10.1080/00268978500103191) are also shown for comparison.
 He–HeHe–NeHe–ArNe–NeNe–ArAr–Ar
This work1.5433.0479.5396.11918.6059.77
Experimenta1.4583.0299.5386.38319.5064.30
aDOSD values from Ref. 107107. A. Kumar and W. J. Meath, Mol. Phys. 54, 823 (1985). https://doi.org/10.1080/00268978500103191.
Next, the C6, C8, and C10 coefficients between rare-gas diatoms (100 Å parted) were calculated using Eqs. (24) and (28), with various values of N. The results for neon dimer are shown in Table II. As seen in the table, calculated dispersion coefficients converge rapidly to the values obtained with Eq. (29) and its analogs for C8 and C10, which correspond to the exact integration of Eq. (26). In the following calculations, N=12 is used unless otherwise noted.
Table icon
Table II. Dispersion coefficients (a.u.) of neon dimer calculated using Eqs. (24), (27), and (28) (using the LC-BOP/aug-cc-pVQZ SCF density), varying N in Eq. (28).
NC6C8C10
66.1312163.203264.3
86.1228163.133263.7
106.1204163.113263.6
126.1196163.103263.6
146.1193163.103263.6
166.1191163.103263.6
186.1191163.103263.6
a6.1189163.103263.5
Referenceb6.552790.3441535.6
aFrom Eq. (29) and its analogs to C8 and C10.
bMany-body perturbation theory values from Ref. 108108. A. J. Thakkar, H. Hettema, and P. E. S. Wormer, J. Chem. Phys. 97, 3252 (1992). https://doi.org/10.1063/1.463012. DOSD value reads C6=6.383 in Ref. 107107. A. Kumar and W. J. Meath, Mol. Phys. 54, 823 (1985). https://doi.org/10.1080/00268978500103191.
Table II reveals disappointing failure in reproducing ab initio C8 and C10. This was also found for other rare-gas diatoms. In contrast, Becke and Johnson97,9897. A. D. Becke and E. R. Johnson, J. Chem. Phys. 124, 014104 (2006). https://doi.org/10.1063/1.213966898. A. D. Becke and E. R. Johnson, J. Chem. Phys. 127, 154108 (2007). https://doi.org/10.1063/1.2795701 obtained reasonable values for the higher-order coefficients. We note that the problem is not due to the multipole expansion, rather it is reflecting the intrinsic accuracy of the local response approximation of Eqs. (6) and (9). This point is clearly shown in Fig. 2, in which the potential energy curves of neon dimer obtained by the LC-BOP, the LC-BOP+LRD without damping, and the LC-BOP+NLF are compared to the CCSD(T) potential.109109. T. J. Giese, V. M. Audette, and D. M. York, J. Chem. Phys. 119, 2618 (2003). https://doi.org/10.1063/1.1587684 As shown in the inset of Fig. 2, at the long distance, the LC-BOP+LRD potential curves converge to the LC-BOP+NLF curve in increasing the interaction order. This reflects the fact that the LRD corresponds to the multipole expansion of the NLF of Eq. (7) for interactions between free atoms. The overestimation of higher-order coefficients in Table II is consistent with the overstabilized LC-BOP+NLF potential compared to the CCSD(T) curve. However, the overall agreement of LC-BOP+NLF and CCSD(T) potential curves is reasonably good. This fact, with previous success of similar LC-BOP+ALL approach,31,33–3531. M. Kamiya, T. Tsuneda, and K. Hirao, J. Chem. Phys. 117, 6010 (2002). https://doi.org/10.1063/1.150113233. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/0026897041233133347434. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.201139635. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243 leaves us optimistic about the above disagreement.
For practical applications of the LRD method, the multipole expansion has to be truncated at a certain order. As seen in Fig. 2, inclusion of higher-order interactions results in better agreement with the NLF at the long range but encounters severe divergence at the short range. Considering these aspects, in the present article, we truncate the expansion after C10 terms and use the damping function of the exponential form discussed in Sec. II. In order to fix the parameters in Eq. (23), we first optimized R¯ values in Eq. (22) for each rare-gas interaction to reproduce the experimental equilibrium distances.110,111110. J. F. Ogilvie and F. Y. H. Wang, J. Mol. Struct. 273, 277 (1992). https://doi.org/10.1016/0022-2860(92)87094-C111. J. F. Ogilvie and F. Y. H. Wang, J. Mol. Struct. 291, 313 (1993). https://doi.org/10.1016/0022-2860(93)85053-W Then the optimal values were least square fitted to the linear form of Eq. (23) with the polarizabilities of constituent atoms calculated at the equilibrium distance. The fitting was nicely done with an R2 value greater than 0.989, giving κ=0.64192a.u. and R0=3.2925a.u.
Potential energy curves of neon dimer obtained by the LC-BOP+LRD including the damping function are shown in Fig. 3. The LC-BOP+LRD curve with up to C10 contributions well reproduces the CCSD(T) curve. Other rare-gas diatoms were also described accurately. In Table III, the binding energies and equilibrium bond distances of the rare-gas diatoms calculated by the LC-BOP+LRD are listed with the CCSD(T) results.109109. T. J. Giese, V. M. Audette, and D. M. York, J. Chem. Phys. 119, 2618 (2003). https://doi.org/10.1063/1.1587684 The mean absolute deviations (MADs) from the CCSD(T) results for binding energies and bond distances are 0.013 kcal/mol and 0.04 Å, respectively.
Table icon
Table III. Binding energies (kcal/mol) and equilibrium bond lengths (Å) of rare-gas diatoms calculated by the LC-BOP+LRD/aug-cc-pVQZ. The CCSD(T) results (from Ref. 109109. T. J. Giese, V. M. Audette, and D. M. York, J. Chem. Phys. 119, 2618 (2003). https://doi.org/10.1063/1.1587684) are also shown for comparison.
 He–HeHe–NeHe–ArNe–NeNe–ArAr–Ar
Binding energy
LC-BOP+LRD0.0250.0470.0740.0860.1470.306
CCSD(T)a0.0210.0420.0590.0820.1290.277
 
Equilibrium bond length
LC-BOP+LRD3.033.073.503.143.523.85
CCSD(T)a2.983.033.493.103.503.78
aReference 109109. T. J. Giese, V. M. Audette, and D. M. York, J. Chem. Phys. 119, 2618 (2003). https://doi.org/10.1063/1.1587684.
B. S22 database: Initial application
In this subsection, the LC-BOP+LRD method is applied to interaction energy (ΔE) calculations of S22 database complexes introduced by Jurečka et al.112112. P. Jurečka, J. Šponer, J. Černý, and P. Hobza, Phys. Chem. Chem. Phys. 8, 1985 (2006). https://doi.org/10.1039/b600027d The database has been constructed to be a reliable accuracy checker of a newly developed method. It includes 22 weakly bound complexes conveniently grouped into hydrogen-bonded, dispersion-dominated, and mixed complexes. Figure 4 shows the geometries of the complexes with their point group symmetries. The estimated complete basis set limit CCSD(T) ΔE’s reported in Ref. 112112. P. Jurečka, J. Šponer, J. Černý, and P. Hobza, Phys. Chem. Chem. Phys. 8, 1985 (2006). https://doi.org/10.1039/b600027d are accurate enough (with errors within few tenths of kcal/mol) for judging the overall performance of a tested method. The 6-311++G(3df,3pd) basis set is used in this subsection. Single-point energy calculations are performed at the reference geometries. The LC-BOP+LRD ΔE calculations are performed within a supermolecular approach: The total energies including dispersion contributions from Eq. (1) calculated for a dimer and monomers are subtracted. This is different from the LC-BOP+NLF calculations, in which the NLF is explicitly applied between preindicated monomers and added to the supermolecular LC-BOP ΔE.33–3533. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/0026897041233133347434. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.201139635. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243
In Table IV, ΔE’s of S22 molecules calculated by the LC-BOP, the LC-BOP+LRD with up to C6, C8, and C10 contributions, and the LC-BOP+NLF are listed with the reference CCSD(T) results. The LC-BOP functional totally underestimates ΔE’s, showing the general importance of the dispersion interaction. Inclusion of the C6 terms roughly halves the MAD and mean absolute percentage deviation (MAPD) from the CCSD(T) ΔE’s, but still leaves large underestimations especially for the dispersion-dominated complexes. This is unlike the good performance with only C6 terms in empirical DFT-D approaches,38,56,5738. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b56. S. Grimme, J. Comput. Chem. 25, 1463 (2004). https://doi.org/10.1002/jcc.2007857. S. Grimme, J. Comput. Chem. 27, 1787 (2006). https://doi.org/10.1002/jcc.20495 in which the strength56,5756. S. Grimme, J. Comput. Chem. 25, 1463 (2004). https://doi.org/10.1002/jcc.2007857. S. Grimme, J. Comput. Chem. 27, 1787 (2006). https://doi.org/10.1002/jcc.20495 and/or the damping38,55–5838. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b55. Q. Wu and W. Yang, J. Chem. Phys. 116, 515 (2002). https://doi.org/10.1063/1.142492856. S. Grimme, J. Comput. Chem. 25, 1463 (2004). https://doi.org/10.1002/jcc.2007857. S. Grimme, J. Comput. Chem. 27, 1787 (2006). https://doi.org/10.1002/jcc.2049558. P. Jurečka, J. Černý, P. Hobza, and D. R. Salahub, J. Comput. Chem. 28, 555 (2007). https://doi.org/10.1002/jcc.20570 of the C6 interactions are optimized without higher-order terms. It is encouraging that inclusions of C8 and C10 terms steadily reduce the deviations from the CCSD(T) ΔE’s for all types of interactions. The MAD and MAPD for the full dataset obtained by the LC-BOP+LRD (implying inclusion of up to C10 terms from here on) are only 0.27 kcal/mol and 5.7%, respectively. The accuracy is comparable to those reported for other recent approaches. For instance, among the work on empirical dispersion correction, MAD (MAPD) values obtained for the S22 set by wB97XD (Ref. 3838. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b) and B97D (Ref. 5757. S. Grimme, J. Comput. Chem. 27, 1787 (2006). https://doi.org/10.1002/jcc.20495) functionals have been reported3838. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b to be 0.22 kcal/mol (5.4%) and 0.50 kcal/mol (6.4%), respectively [with 6-311++G(3df,3pd) basis, counterpoise corrected]. The series of Minnesota hybrid meta-GGA functionals, e.g., M05-2X (Ref. 113113. Y. Zhao and D. G. Truhlar, J. Chem. Theory Comput. 3, 289 (2007). https://doi.org/10.1021/ct6002719) and M06-2X,4545. Y. Zhao and D. G. Truhlar, Theor. Chem. Acc. 120, 215 (2008). https://doi.org/10.1007/s00214-007-0310-x have shown similar performance.
Table icon
Table IV. Interaction energies (kcal/mol) of S22 complexes calculated by the LC-BOP, the LC-BOP+LRD, and the LC-BOP+NLF methods [with 6-311++G(3df,3pd) basis set]. Estimated CCSD(T)/complete basis set values (from Ref. 112112. P. Jurečka, J. Šponer, J. Černý, and P. Hobza, Phys. Chem. Chem. Phys. 8, 1985 (2006). https://doi.org/10.1039/b600027d) are also listed for comparison.
  LC-BOPLC-BOP+LRDLCBOP+NLFCCSD(T)d
6a6+8 b6+8+10 c
Hydrogen-bonded complexes 
1(NH3)2−2.50−2.82−2.95−3.01−3.17−3.17
2(H2O)2−4.88−5.09−5.19−5.23−5.33−5.02
3(HCOOH)2−18.94−19.77−20.11−20.24−20.58−18.61
4(HCONH2)2−14.89−15.75−16.10−16.23−16.65−15.96
5Uracil dimer−18.70−19.94−20.42−20.59−20.96−20.65
62-pyridoxine-2-aminopyridine−14.30−15.80−16.41−16.64−17.07−16.71
7A-T WC−13.64−15.32−16.01−16.29−16.71−16.37
 
Dispersion-dominant complexes
8(CH4)20.10−0.31−0.54−0.63−0.61−0.53
9(C2H4)2−0.15−0.95−1.33−1.47−1.65−1.51
10Benzene-CH40.11−0.81−1.27−1.45−1.78−1.50
11PD-benzene dimmer2.20−0.18−1.86−2.86−2.89−2.73
12Pyrazine dimmer0.73−1.72−3.39−4.33−4.18−4.42
13Stacked uracil dimmer−3.29−6.64−8.78−9.93−9.98−10.12
14Stacked indole-benzene2.55−0.88−3.28−4.51−4.57−5.22
15Stacked A-T−2.14−7.15−10.33−11.93−11.90−12.23
 
Mixed complexes
16Ethene-ethyne−0.86−1.30−1.54−1.63−1.64−1.53
17Benzene-H2O−2.08−2.95−3.43−3.66−3.74−3.28
18Benzene-NH3−0.88−1.78−2.27−2.49−2.68−2.35
19Benzene-HCN−3.17−4.17−4.65−4.87−5.24−4.46
20T-shaped benzene dimmer−0.05−1.49−2.20−2.48−3.02−2.74
21T-shaped indole-benzene−1.95−3.84−4.74−5.13−5.79−5.73
22Phenol dimer−4.29−5.93−6.75−7.11−7.39−7.05
 
Overall
MDe(kcal/mol)2.771.240.38−0.04−0.26 
MADf(kcal/mol)2.801.350.580.270.38 
Rangeg(kcal/mol)10.426.243.442.342.62 
MPDh(%)61.727.78.2−0.9−5.7 
MAPDi(%)61.828.410.45.77.6 
aWith C6 terms.
bWith C6 and C8 terms.
cWith C6, C8, and C10 terms.
dFrom Ref. 112112. P. Jurečka, J. Šponer, J. Černý, and P. Hobza, Phys. Chem. Chem. Phys. 8, 1985 (2006). https://doi.org/10.1039/b600027d.
eMean deviation from the CCSD(T) values.
fMean absolute deviation.
gThe largest positive minus the largest negative deviations.
hMean percentage deviation.
iMean absolute percentage deviation.
We emphasize the notable agreement between LC-BOP+LRD and LC-BOP+NLF ΔE’s, especially for the stacking interactions (Nos. 11–15). This suggests the adequacy of the seemingly drastic assumption of the local atomic polarizability; see discussions above Eq. (16). The larger difference [relative to the CCSD(T) ΔE] between LRD and NLF dispersion corrections are found for so-called XH-π interactions (Nos. 19–21). The difference is mainly due to the different short-range damping: The damping has a stronger impact on calculated ΔE’s for these spatially closer interactions. On the other hand, the hydrogen-bonded complexes, which involve even shorter interactions, are described similarly with both methods, reflecting the minor importance of the dispersion interaction in these complexes.
Efficiency of the LRD method is highlighted in Table V, which summarizes the calculation of stacked adenine-thymine complex (No. 15). In the third column of the table, ΔE’s calculated by the LC-BOP+LRD (with various number N of Gauss–Chebyshev quadrature points) and the LC-BOP+NLF are listed. In increasing N, the LC-BOP+LRD ΔE’s come close to the LC-BOP+NLF one. The default value of N=12 is justified here. The fourth and fifth columns show central processor unit (CPU) times for two computational bottlenecks in the LRD and NLF computations, respectively. These two steps for the LRD procedure are the calculation of atomic polarizabilities of Eq. (27) [Pol] and evaluation of Eqs. (1), (24), and (28) [Sum], while those for the NLF are the preparation of the integrand [Pol] and the doubly numerical evaluation [Sum] of Eq. (7). Comparison of the “Sum” parts of CPU times for the LRD and NLF makes clear the benefit of replacing the computationally demanding double spatial integral with atomic pair summations. It is also noted that the “Pol” timings increase only marginally against N, which indicates that the cost for Eq. (27) is dominated by the N-independent procedure, i.e., the density evaluation at the spatial grid points. This numerically confirms the theoretical arguments made at the end of Sec. II D.
Table icon
Table V. Interaction energies (kcal/mol) of stacked adenine-thymine dimer calculated by the LC-BOP+LRD and the LC-BOP+NLF methods [with 6-311++G(3df,3pd) basis set], and CPU time (on power 6/4.7GHz machine) (second) for the LRD and NLF dispersion corrections, varying N in Eq. (28).
MethodNΔECPU time
PolaSumb
LC-BOP+LRD6−11.63193.20.3
 8−11.82193.20.3
 10−11.90193.40.4
 12−11.93193.60.4
 14−11.95193.80.5
 16−11.95194.00.5
 18−11.96194.20.6
LC-BOP+NLF −11.90192.41014.5
aCPU time for calculating atomic polarizabilities of Eq. (27) for the LRD or preparing the integrand of Eq. (7) for the NLF.
bCPU time for numerical evaluation of Eqs. (1), (24), and (28) for the LRD or that of Eq. (7) for the NLF.
In this article, we propose a new method for calculating system-dependent dispersion coefficients between atoms in a molecule, called the LRD method. Based on the local response approximation due to Dobson and Dinte,6767. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780 the modified dispersion relation of Vydrov and van Voorhis,3636. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.3079684 and the Becke-type atomic partition function,8888. A. D. Becke, J. Chem. Phys. 88, 2547 (1988). https://doi.org/10.1063/1.454033 the multicenter multipole expanded expression of the second-order dispersion energy is derived. The distributed multipole polarizabilities for atoms in a molecule are first calculated from the local response model, from which the dispersion coefficients are obtained by an explicit frequency integral of the Casimir–Polder type. This procedure replaces the costly double spatial integral with the simpler double atomic summation without sacrificing accuracy. A new damping function is also introduced which possesses desirable features of both Fermi-type5555. Q. Wu and W. Yang, J. Chem. Phys. 116, 515 (2002). https://doi.org/10.1063/1.1424928 and Tang–Toennis9393. K. T. Tang and J. P. Toennies, J. Chem. Phys. 80, 3726 (1984). https://doi.org/10.1063/1.447150 functions, and takes into account the effective volume of atoms in molecular environment.
Thus formulated LRD method is combined with the LC-BOP functional via the damped atom-atom form and applied to calculations of weakly bound systems. Three empirical parameters included in the LRD method are first determined based on calculations for rare-gas diatoms. The method with the fixed parameters is then applied to interaction energy calculations of S22 database.112112. P. Jurečka, J. Šponer, J. Černý, and P. Hobza, Phys. Chem. Chem. Phys. 8, 1985 (2006). https://doi.org/10.1039/b600027d The mean absolute percent deviation of the calculated interaction energies from the CCSD(T) results is obtained as 5.7%, which is comparable to those obtained by other latest approaches. More importantly, it is confirmed that the LRD method performs similarly to the NLF approach with much less computational cost and no system fragmentation.
The present LRD method has a number of advantages in theoretical and practical points of view. It computes dispersion energy directly from the ground-state electron density. It is applicable to any geometries, free from outer sources of physical constants, and computationally no more demanding than the usual KS calculations.
The known problem is the overestimation of the higher-order dispersion coefficients observed in rare-gas interactions. This may indicate the need of further sophistication of the local response model. However, we are currently optimistic about this problem since the LRD method roughly corresponds to the multipole expansion of the parent NLF, which has been confirmed fairly accurate.31,33–3531. M. Kamiya, T. Tsuneda, and K. Hirao, J. Chem. Phys. 117, 6010 (2002). https://doi.org/10.1063/1.150113233. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/0026897041233133347434. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.201139635. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243 A more fundamental defect of the LRD and any DFT-D type approaches is the need of damping functions. We manage to increase its transferability by introducing a system dependence to the damping. It is clear that such a prescription is far from satisfactory from a purely theoretical stand point; the seamless functional approaches36,37,6836. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.307968437. O. A. Vydrov and T. van Voorhis, Phys. Rev. Lett. 103, 063004 (2009). https://doi.org/10.1103/PhysRevLett.103.06300468. M. Dion, H. Rydberg, E. Schröder, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 92, 246401 (2004). https://doi.org/10.1103/PhysRevLett.92.246401 or more rigorous treatments40–4340. B. G. Janesko, T. M. Henderson, and G. E. Scuseria, J. Chem. Phys. 130, 081105 (2009). https://doi.org/10.1063/1.309081441. B. G. Janesko, T. M. Henderson, and G. E. Scuseria, J. Chem. Phys. 131, 034110 (2009). https://doi.org/10.1063/1.317651442. J. G. Ángyán, I. C. Gerber, A. Savin, and J. Toulouse, Phys. Rev. A 72, 012510 (2005). https://doi.org/10.1103/PhysRevA.72.01251043. J. Toulouse, I. C. Gerber, G. Jansen, A. Savin, and J. G. Ángyán, Phys. Rev. Lett. 102, 096404 (2009). https://doi.org/10.1103/PhysRevLett.102.096404 of the long-ranged electron correlation are desired. However, the LRD method has an opportunity to provide an ideal cost/performance balance in large-scale calculations, which we believe is one of the most important features of KS DFT.
Another strength of KS method is the ease of calculating analytical energy derivatives with respect to arbitrary perturbations. In particular, the geometry optimization with analytical gradients is probably the most common application of DFT. In this respect, the development of analytical derivatives for the LRD method based on the SCF treatment is indispensable, which should be our next work.
This research was supported in part by a Grant-in-Aid for Scientific Research on priority areas “Molecular Theory for Real Systems,” under Grant No. KAKENHI 18066016 from the Ministry of Education, Culture, Sports, Science and Technology (MEXT), Japan, by the Global Center Of Excellence (COE) Practical Chemical Wisdom from the MEXT, and by a Grant-in-Aid for Young Scientist (Grant No. WAKATE-B-21750025) from the MEXT. We were also partly supported by the project research grant “Development of high-performance computational environment for quantum chemical calculations and its assessment” from the Research Institute for Science and Engineering (RISE), Waseda University. Calculations were performed at the Research Center for Computational Science (RCCS), Okazaki Research Facilities, National Institutes of Natural Sciences (NINS).
  1. 1. A. J. Stone, Theory of Intermolecular Forces (Clarendon, Oxford, 1996). Google Scholar
  2. 2. J. Šponer, K. E. Riley, and P. Hobza, Phys. Chem. Chem. Phys. 10, 2595 (2008). https://doi.org/10.1039/b719370j, Google ScholarCrossref, ISI
  3. 3. V. Špirko, O. Engkvist, P. Soldán, H. L. Selzle, E. W. Schlag, and P. Hobza, J. Chem. Phys. 111, 572 (1999). https://doi.org/10.1063/1.479338, Google ScholarScitation
  4. 4. S. Tsuzuki, K. Honda, T. Uchimaru, M. Mikami, and K. Tanabe, J. Am. Chem. Soc. 124, 104 (2002). https://doi.org/10.1021/ja0105212, Google ScholarCrossref, ISI
  5. 5. M. O. Sinnokrot, E. F. Valeev, and C. D. Sherrill, J. Am. Chem. Soc. 124, 10887 (2002). https://doi.org/10.1021/ja025896h, Google ScholarCrossref, ISI
  6. 6. P. Hohenberg and W. Kohn, Phys. Rev. 136, B864 (1964). https://doi.org/10.1103/PhysRev.136.B864, Google ScholarCrossref, ISI
  7. 7. W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965). https://doi.org/10.1103/PhysRev.140.A1133, Google ScholarCrossref, ISI
  8. 8. J. P. Perdew and W. Yue, Phys. Rev. B 33, 8800 (1986). https://doi.org/10.1103/PhysRevB.33.8800, Google ScholarCrossref, ISI
  9. 9. A. D. Becke, Phys. Rev. A 38, 3098 (1988). https://doi.org/10.1103/PhysRevA.38.3098, Google ScholarCrossref, ISI
  10. 10. C. Lee, W. Yang, and R. G. Parr, Phys. Rev. B 37, 785 (1988). https://doi.org/10.1103/PhysRevB.37.785, Google ScholarCrossref, ISI
  11. 11. J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). https://doi.org/10.1103/PhysRevLett.77.3865, Google ScholarCrossref, ISI
  12. 12. E. I. Proynov, A. Vela, and D. R. Salahub, Chem. Phys. Lett. 230, 419 (1994). https://doi.org/10.1016/0009-2614(94)01189-3, Google ScholarCrossref
  13. 13. A. D. Becke, J. Chem. Phys. 104, 1040 (1996). https://doi.org/10.1063/1.470829, Google ScholarScitation, ISI
  14. 14. T. Van Voorhis and G. E. Scuseria, J. Chem. Phys. 109, 400 (1998). https://doi.org/10.1063/1.476577, Google ScholarScitation, ISI
  15. 15. A. D. Becke, J. Chem. Phys. 109, 2092 (1998). https://doi.org/10.1063/1.476722, Google ScholarScitation, ISI
  16. 16. M. Filatov and W. Thiel, Phys. Rev. A 57, 189 (1998). https://doi.org/10.1103/PhysRevA.57.189, Google ScholarCrossref
  17. 17. J. P. Perdew, S. Kurth, A. Zupan, and P. Blaha, Phys. Rev. Lett. 82, 2544 (1999). https://doi.org/10.1103/PhysRevLett.82.2544, Google ScholarCrossref, ISI
  18. 18. A. D. Boese and N. C. Handy, J. Chem. Phys. 116, 9559 (2002). https://doi.org/10.1063/1.1476309, Google ScholarScitation, ISI
  19. 19. J. Tao, J. P. Perdew, V. N. Staroverov, and G. E. Scuseria, Phys. Rev. Lett. 91, 146401 (2003). https://doi.org/10.1103/PhysRevLett.91.146401, Google ScholarCrossref, ISI
  20. 20. A. D. Becke, J. Chem. Phys. 98, 1372 (1993). https://doi.org/10.1063/1.464304, Google ScholarScitation, ISI
  21. 21. A. D. Becke, J. Chem. Phys. 98, 5648 (1993). https://doi.org/10.1063/1.464913, Google ScholarScitation, ISI
  22. 22. P. J. Stephens, F. J. Devlin, C. F. Chabalowski, and M. J. Frisch, J. Phys. Chem. 98, 11623 (1994). https://doi.org/10.1021/j100096a001, Google ScholarCrossref, ISI
  23. 23. S. Kristyán and P. Pulay, Chem. Phys. Lett. 229, 175 (1994). https://doi.org/10.1016/0009-2614(94)01027-7, Google ScholarCrossref, ISI
  24. 24. J. M. Pérez-Jordá and A. D. Becke, Chem. Phys. Lett. 233, 134 (1995). https://doi.org/10.1016/0009-2614(94)01402-H, Google ScholarCrossref, ISI
  25. 25. D. J. Lacks and R. G. Gordon, Phys. Rev. A 47, 4681 (1993). https://doi.org/10.1103/PhysRevA.47.4681, Google ScholarCrossref
  26. 26. T. A. Wesolowski, O. Parisel, Y. Ellinger, and J. Weber, J. Phys. Chem. A 101, 7818 (1997). https://doi.org/10.1021/jp970586k, Google ScholarCrossref, ISI
  27. 27. Y. Zhang, W. Pan, and W. Yang, J. Chem. Phys. 107, 7921 (1997). https://doi.org/10.1063/1.475105, Google ScholarScitation, ISI
  28. 28. C. Adamo and V. Barone, J. Chem. Phys. 108, 664 (1998). https://doi.org/10.1063/1.475428, Google ScholarScitation, ISI
  29. 29. N. Kurita and H. Sekino, Chem. Phys. Lett. 348, 139 (2001). https://doi.org/10.1016/S0009-2614(01)01089-2, Google ScholarCrossref, ISI
  30. 30. X. Xu and W. A. Goddard, Proc. Natl. Acad. Sci. U.S.A. 101, 2673 (2004). https://doi.org/10.1073/pnas.0308730100, Google ScholarCrossref, ISI
  31. 31. M. Kamiya, T. Tsuneda, and K. Hirao, J. Chem. Phys. 117, 6010 (2002). https://doi.org/10.1063/1.1501132, Google ScholarScitation, ISI
  32. 32. H. Iikura, T. Tsuneda, T. Yanai, and K. Hirao, J. Chem. Phys. 115, 3540 (2001). https://doi.org/10.1063/1.1383587, Google ScholarScitation, ISI
  33. 33. T. Sato, T. Tsuneda, and K. Hirao, Mol. Phys. 103, 1151 (2005). https://doi.org/10.1080/00268970412331333474, Google ScholarCrossref, ISI
  34. 34. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 123, 104307 (2005). https://doi.org/10.1063/1.2011396, Google ScholarScitation, ISI
  35. 35. T. Sato, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 234114 (2007). https://doi.org/10.1063/1.2747243, Google ScholarScitation, ISI
  36. 36. O. A. Vydrov and T. van Voorhis, J. Chem. Phys. 130, 104105 (2009). https://doi.org/10.1063/1.3079684, Google ScholarScitation, ISI
  37. 37. O. A. Vydrov and T. van Voorhis, Phys. Rev. Lett. 103, 063004 (2009). https://doi.org/10.1103/PhysRevLett.103.063004, Google ScholarCrossref, ISI
  38. 38. J. -D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys. 10, 6615 (2008). https://doi.org/10.1039/b810189b, Google ScholarCrossref, ISI
  39. 39. J. Gräfenstein and D. Cremer, J. Chem. Phys. 130, 124105 (2009). https://doi.org/10.1063/1.3079822, Google ScholarScitation, ISI
  40. 40. B. G. Janesko, T. M. Henderson, and G. E. Scuseria, J. Chem. Phys. 130, 081105 (2009). https://doi.org/10.1063/1.3090814, Google ScholarScitation, ISI
  41. 41. B. G. Janesko, T. M. Henderson, and G. E. Scuseria, J. Chem. Phys. 131, 034110 (2009). https://doi.org/10.1063/1.3176514, Google ScholarScitation, ISI
  42. 42. J. G. Ángyán, I. C. Gerber, A. Savin, and J. Toulouse, Phys. Rev. A 72, 012510 (2005). https://doi.org/10.1103/PhysRevA.72.012510, Google ScholarCrossref, ISI
  43. 43. J. Toulouse, I. C. Gerber, G. Jansen, A. Savin, and J. G. Ángyán, Phys. Rev. Lett. 102, 096404 (2009). https://doi.org/10.1103/PhysRevLett.102.096404, Google ScholarCrossref, ISI
  44. 44. Y. Zhao, N. E. Schltz, and D. G. Truhlar, J. Chem. Theory Comput. 2, 364 (2006). https://doi.org/10.1021/ct0502763, Google ScholarCrossref, ISI
  45. 45. Y. Zhao and D. G. Truhlar, Theor. Chem. Acc. 120, 215 (2008). https://doi.org/10.1007/s00214-007-0310-x, Google ScholarCrossref, ISI
  46. 46. Y. Zhang, A. Vela, and D. R. Salahub, Theor. Chem. Acc. 118, 693 (2007). https://doi.org/10.1007/s00214-007-0347-x, Google ScholarCrossref
  47. 47. O. A. von Lilienfeld, I. Tavernelli, U. Rothlisberger, and D. Sebastiani, Phys. Rev. Lett. 93, 153004 (2004). https://doi.org/10.1103/PhysRevLett.93.153004, Google ScholarCrossref, ISI
  48. 48. O. A. von Lilienfeld, I. Tavernelli, U. Rothlisberger, and D. Sebastiani, Phys. Rev. B 71, 195119 (2005). https://doi.org/10.1103/PhysRevB.71.195119, Google ScholarCrossref, ISI
  49. 49. I. -C. Lin, M. D. Coutinho-Neto, C. Felsenheimer, O. A. von Lilienfeld, I. Tavernelli, and U. Rothlisberger, Phys. Rev. B 75, 205131 (2007). https://doi.org/10.1103/PhysRevB.75.205131, Google ScholarCrossref, ISI
  50. 50. P. C. Aeberhard, J. S. Arey, I. -C. Lin, and U. Rothlisberger, J. Chem. Theory Comput. 5, 23 (2009). https://doi.org/10.1021/ct800299y, Google ScholarCrossref
  51. 51. G. A. DiLabio, Chem. Phys. Lett. 455, 348 (2008). https://doi.org/10.1016/j.cplett.2008.02.110, Google ScholarCrossref
  52. 52. I. D. Mackie and G. A. DiLabio, J. Phys. Chem. A 112, 10968 (2008). https://doi.org/10.1021/jp806162t, Google ScholarCrossref, ISI
  53. 53. J. F. Dobson, K. McLennan, A. Rubio, J. Wang, T. Gould, H. M. Le, and B. P. Dinte, Aust. J. Chem. 54, 513 (2001). https://doi.org/10.1071/CH01052, Google ScholarCrossref, ISI
  54. 54. X. Wu, M. C. Vargas, S. Nayak, V. Lotrich, and G. Scoles, J. Chem. Phys. 115, 8748 (2001). https://doi.org/10.1063/1.1412004, Google ScholarScitation, ISI
  55. 55. Q. Wu and W. Yang, J. Chem. Phys. 116, 515 (2002). https://doi.org/10.1063/1.1424928, Google ScholarScitation, ISI
  56. 56. S. Grimme, J. Comput. Chem. 25, 1463 (2004). https://doi.org/10.1002/jcc.20078, Google ScholarCrossref, ISI
  57. 57. S. Grimme, J. Comput. Chem. 27, 1787 (2006). https://doi.org/10.1002/jcc.20495, Google ScholarCrossref, ISI
  58. 58. P. Jurečka, J. Černý, P. Hobza, and D. R. Salahub, J. Comput. Chem. 28, 555 (2007). https://doi.org/10.1002/jcc.20570, Google ScholarCrossref, ISI
  59. 59. A. D. Becke and E. R. Johnson, J. Chem. Phys. 122, 154104 (2005). https://doi.org/10.1063/1.1884601, Google ScholarScitation, ISI
  60. 60. A. D. Becke and E. R. Johnson, J. Chem. Phys. 123, 154101 (2005). https://doi.org/10.1063/1.2065267, Google ScholarScitation, ISI
  61. 61. A. Olasz, K. Vanommeslaeghe, A. Krishtal, T. Veszprémi, C. V. Alsenoy, and P. Geerlings, J. Chem. Phys. 127, 224105 (2007). https://doi.org/10.1063/1.2805391, Google ScholarScitation
  62. 62. A. Krishtal, K. Vanommeslaeghe, A. Olasz, T. Veszprémi, C. V. Alsenoy, and P. Geerlings, J. Chem. Phys. 130, 174101 (2009). https://doi.org/10.1063/1.3126248, Google ScholarScitation
  63. 63. F. O. Kannemann and A. D. Becke, J. Chem. Theory Comput. 5, 719 (2009). https://doi.org/10.1021/ct800522r, Google ScholarCrossref, ISI
  64. 64. J. Kong, Z. Gan, E. Proynov, M. Freindorf, and T. R. Furlani, Phys. Rev. A 79, 042510 (2009). https://doi.org/10.1103/PhysRevA.79.042510, Google ScholarCrossref
  65. 65. K. Rapcewicz and N. W. Ashcroft, Phys. Rev. B 44, 4032 (1991). https://doi.org/10.1103/PhysRevB.44.4032, Google ScholarCrossref
  66. 66. Y. Andersson, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 76, 102 (1996). https://doi.org/10.1103/PhysRevLett.76.102, Google ScholarCrossref, ISI
  67. 67. J. F. Dobson and B. P. Dinte, Phys. Rev. Lett. 76, 1780 (1996). https://doi.org/10.1103/PhysRevLett.76.1780, Google ScholarCrossref
  68. 68. M. Dion, H. Rydberg, E. Schröder, D. C. Langreth, and B. I. Lundqvist, Phys. Rev. Lett. 92, 246401 (2004). https://doi.org/10.1103/PhysRevLett.92.246401, Google ScholarCrossref, ISI
  69. 69. P. L. Silvestrelli, Phys. Rev. Lett. 100, 053002 (2008). https://doi.org/10.1103/PhysRevLett.100.053002, Google ScholarCrossref, ISI
  70. 70. P. L. Silvestrelli, J. Phys. Chem. A 113, 5224 (2009). https://doi.org/10.1021/jp811138n, Google ScholarCrossref
  71. 71. P. L. Silvestrelli, K. Benyahia, S. Grubisiĉ, F. Ancilotto, and F. Toigo, J. Chem. Phys. 130, 074702 (2009). https://doi.org/10.1063/1.3077288, Google ScholarScitation, ISI
  72. 72. S. D. Chakarova-Käck, E. Schröder, B. I. Lundqvist, and D. C. Langreth, Phys. Rev. Lett. 96, 146107 (2006). https://doi.org/10.1103/PhysRevLett.96.146107, Google ScholarCrossref
  73. 73. A. Puzder, M. Dion, and D. C. Langreth, J. Chem. Phys. 124, 164105 (2006). https://doi.org/10.1063/1.2189229, Google ScholarScitation, ISI
  74. 74. V. R. Cooper, T. Thonhauser, and D. C. Langreth, J. Chem. Phys. 128, 204102 (2008). https://doi.org/10.1063/1.2924133, Google ScholarScitation
  75. 75. Y. Zhang and W. Yang, Phys. Rev. Lett. 80, 890 (1998). https://doi.org/10.1103/PhysRevLett.80.890, Google ScholarCrossref, ISI
  76. 76. O. A. Vydrov, Q. Wu, and T. van Voorhis, J. Chem. Phys. 129, 014106 (2008). https://doi.org/10.1063/1.2948400, Google ScholarScitation, ISI
  77. 77. O. A. Vydrov and G. E. Scuseria, J. Chem. Phys. 125, 234109 (2006). https://doi.org/10.1063/1.2409292, Google ScholarScitation, ISI
  78. 78. A. Gulans, M. J. Puska, and R. M. Nieminen, Phys. Rev. B 79, 201105 (2009). https://doi.org/10.1103/PhysRevB.79.201105, Google ScholarCrossref, ISI
  79. 79. G. Román-Pérez and J. M. Soler, Phys. Rev. Lett. 103, 096102 (2009). https://doi.org/10.1103/PhysRevLett.103.096102, Google ScholarCrossref, ISI
  80. 80. H. C. Longuet-Higgins, Discuss. Faraday Soc. 40, 7 (1965). https://doi.org/10.1039/df9654000007, Google ScholarCrossref
  81. 81. E. Zaremba and W. Kohn, Phys. Rev. B 13, 2270 (1976). https://doi.org/10.1103/PhysRevB.13.2270, Google ScholarCrossref
  82. 82. R. McWeeny, Methods of Molecular Quantum Mechanics (Academic, London, 1992). Google Scholar
  83. 83. D. C. Langreth, M. Dion, H. Rydberg, E. Schröder, P. Hyldgaard, and B. I. Lundqvist, Int. J. Quantum Chem. 101, 599 (2005). https://doi.org/10.1002/qua.20315, Google ScholarCrossref, ISI
  84. 84. A. J. Stone and C. -S. Tong, Chem. Phys. 137, 121 (1989). https://doi.org/10.1016/0301-0104(89)87098-3, Google ScholarCrossref
  85. 85. G. J. Williams and A. J. Stone, J. Chem. Phys. 119, 4620 (2003). https://doi.org/10.1063/1.1594722, Google ScholarScitation, ISI
  86. 86. C. Hättig, G. Jansen, B. A. Hess, and J. G. Ángyán, Mol. Phys. 91, 145 (1997). https://doi.org/10.1080/002689797171841, Google ScholarCrossref
  87. 87. J. G. Ángyán, J. Chem. Phys. 127, 024108 (2007). https://doi.org/10.1063/1.2749512, Google ScholarScitation
  88. 88. A. D. Becke, J. Chem. Phys. 88, 2547 (1988). https://doi.org/10.1063/1.454033, Google ScholarScitation, ISI
  89. 89. J. Hepburn and G. Scoles, Chem. Phys. Lett. 36, 451 (1975). https://doi.org/10.1016/0009-2614(75)80278-8, Google ScholarCrossref, ISI
  90. 90. R. Ahlrichs, R. Penco, and G. Scoles, Chem. Phys. 19, 119 (1977). https://doi.org/10.1016/0301-0104(77)85124-0, Google ScholarCrossref, ISI
  91. 91. K. T. Tang and J. P. Toennies, J. Chem. Phys. 66, 1496 (1977). https://doi.org/10.1063/1.434113, Google ScholarScitation, ISI
  92. 92. K. T. Tang and J. P. Toennies, J. Chem. Phys. 68, 5501 (1978). https://doi.org/10.1063/1.435678, Google ScholarScitation
  93. 93. K. T. Tang and J. P. Toennies, J. Chem. Phys. 80, 3726 (1984). https://doi.org/10.1063/1.447150, Google ScholarScitation, ISI
  94. 94. A. Tkatchenko and M. Scheffler, Phys. Rev. Lett. 102, 073005 (2009). https://doi.org/10.1103/PhysRevLett.102.073005, Google ScholarCrossref, ISI
  95. 95. T. Brinck, J. S. Murray, and P. Politzer, J. Chem. Phys. 98, 4305 (1993). https://doi.org/10.1063/1.465038, Google ScholarScitation, ISI
  96. 96. E. R. Johnson and A. D. Becke, J. Chem. Phys. 124, 174104 (2006). https://doi.org/10.1063/1.2190220, Google ScholarScitation, ISI
  97. 97. A. D. Becke and E. R. Johnson, J. Chem. Phys. 124, 014104 (2006). https://doi.org/10.1063/1.2139668, Google ScholarScitation, ISI
  98. 98. A. D. Becke and E. R. Johnson, J. Chem. Phys. 127, 154108 (2007). https://doi.org/10.1063/1.2795701, Google ScholarScitation, ISI
  99. 99. A. Heßelmann, J. Chem. Phys. 130, 084104 (2009). https://doi.org/10.1063/1.3077939, Google ScholarScitation
  100. 100. N. Marzari and D. Vanderbilt, Phys. Rev. B 56, 12847 (1997). https://doi.org/10.1103/PhysRevB.56.12847, Google ScholarCrossref, ISI
  101. 101. T. Tsuneda, T. Suzumura, and K. Hirao, J. Chem. Phys. 110, 10664 (1999). https://doi.org/10.1063/1.479012, Google ScholarScitation, ISI
  102. 102. J. -W. Song, T. Hirosawa, T. Tsuneda, and K. Hirao, J. Chem. Phys. 126, 154105 (2007). https://doi.org/10.1063/1.2721532, Google ScholarScitation, ISI
  103. 103. M. J. Frisch, G. W. Trucks, H. B. Schlegel et al., GAUSSIAN03, Revision D.02, Gaussian, Inc., Pittsburgh, PA, 2003. Google Scholar
  104. 104. P. M. W. Gill, B. G. Johnson, and J. A. Pople, Chem. Phys. Lett. 209, 506 (1993). https://doi.org/10.1016/0009-2614(93)80125-9, Google ScholarCrossref, ISI
  105. 105. V. I. Lebedev, Dokl. Math. 45, 587 (1992). Google Scholar
  106. 106. S. F. Boys and F. Bernardi, Mol. Phys. 19, 553 (1970). https://doi.org/10.1080/00268977000101561, Google ScholarCrossref, ISI
  107. 107. A. Kumar and W. J. Meath, Mol. Phys. 54, 823 (1985). https://doi.org/10.1080/00268978500103191, Google ScholarCrossref, ISI
  108. 108. A. J. Thakkar, H. Hettema, and P. E. S. Wormer, J. Chem. Phys. 97, 3252 (1992). https://doi.org/10.1063/1.463012, Google ScholarScitation, ISI
  109. 109. T. J. Giese, V. M. Audette, and D. M. York, J. Chem. Phys. 119, 2618 (2003). https://doi.org/10.1063/1.1587684, Google ScholarScitation
  110. 110. J. F. Ogilvie and F. Y. H. Wang, J. Mol. Struct. 273, 277 (1992). https://doi.org/10.1016/0022-2860(92)87094-C, Google ScholarCrossref, ISI
  111. 111. J. F. Ogilvie and F. Y. H. Wang, J. Mol. Struct. 291, 313 (1993). https://doi.org/10.1016/0022-2860(93)85053-W, Google ScholarCrossref
  112. 112. P. Jurečka, J. Šponer, J. Černý, and P. Hobza, Phys. Chem. Chem. Phys. 8, 1985 (2006). https://doi.org/10.1039/b600027d, Google ScholarCrossref, ISI
  113. 113. Y. Zhao and D. G. Truhlar, J. Chem. Theory Comput. 3, 289 (2007). https://doi.org/10.1021/ct6002719, Google ScholarCrossref, ISI
  1. © 2009 American Institute of Physics