Optimal operators for Hartree–Fock exchange from long-range corrected hybrid density functionals
Jeng-Da Chai
*, Martin Head-Gordon
Department of Chemistry, University of California and Chemical Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA
a r t i c l e i n f o
Article history:
Received 14 August 2008 In final form 23 October 2008 Available online 31 October 2008
a b s t r a c t
The long-range operators for Hartree–Fock exchange from two recently proposed long-range corrected hybrid functionals,xB97 andxB97X [J.-D. Chai, M. Head-Gordon, J. Chem. Phys. 128 (2008) 084106], are discussed. A conserved property is found in the middle-range region of the operator fromxB97X.
We argue that the fine details of the Hartree–Fock exchange mixing in this region are responsible for the accuracy of the long-range corrected hybrid functionals in thermochemistry and barrier heights.
Published by Elsevier B.V.
In the last two decades, Kohn–Sham density-functional theory (KS-DFT) [1,2] has been the most popular quantum-chemical method for calculations on large-scale ground-state systems[3–
5]. As a result, developing an accurate exchange-correlation energy functional Exc[
q
], continues being the subject of intense current interest.Semi-local density functionals based on the local spin density approximation (LSDA) and generalized gradient approximations (GGAs) (commonly denoted as DFAs for density-functional approx- imations) are computationally favorable for large systems and have been shown to be successful in diverse quantum-chemical applications[3–5]. However, in circumstances where an accurate treatment of the non-locality of the exchange-correlation (XC) hole is crucial, DFAs can and do produce qualitatively incorrect results [6–15].
On the other hand, hybrid density functionals, combining den- sity functionals with the exact Hartree–Fock (HF) exchange EHFx , provide cost-effective ways of including the non-local effects of the XC hole. The concept of a hybrid scheme can be deduced from the adiabatic connection formalism for Exc[16]. The most widely used hybrid density functionals at present are global hybrids, which can be written as linear combinations of EHFx and EDFAxc : Exc¼ cxEHFx þ ð1 $ cxÞEDFAx þ EDFAc ð1Þ where cxis a small fractional number whose optimal value depends on the systems of interest[5].
Unlike the global hybrid functionals, the long-range corrected (LC) hybrid functionals[14,15,17–32]employ 100% HF exchange for a long-range (LR) part of the interelectron repulsion operator L(r)/r, DFA exchange for the complementary short-range (SR)
operator [1 $ L(r)]/r, and DFA correlation for the entire Coulomb operator 1/r:
ELCxc¼ ELR-HFx þ ESR-DFAx þ EDFAc ð2Þ
where L(r) is a fraction of HF exchange at r.
There are, however, two important issues relevant to the overall accuracy of the LC hybrid functionals. One is to develop accurate DFAs for the SR exchange and correlation, and the other is to seek optimal LR operators for the partition. For the first issue, an accu- rate DFA for SR exchange may be the main concern, as accurate DFAs for correlation are widely available[3–5]. Although the LSDA for SR exchange ESR-LSDAx has been available[24,25], the resulting LC hybrid LSDA functional is insufficiently accurate for typical quan- tum-chemical applications [26]. To make progress, Hirao and co-workers have proposed an ansatz for the generalization of any EDFAx to ESR-DFAx [27,28]. Although their LC hybrid GGA functionals outperform the LC hybrid LSDA functional, they are still inferior to global hybrid GGA functionals for properties such as thermochemistry.
Recently, we have also proposed a simple ansatz to extend any EDFAx to ESR-DFAx , as long as the SR operator has considerable spatial extent[14]. With the use of flexible DFAs, our LC hybrid functionals outperform the corresponding global hybrid functionals, which must be attributed to the difference of the operators used for HF exchange in these functionals (i.e. the second issue). Note that the operator cx/r is effectively used for the HF exchange mixing in the global hybrid functionals. This outcome seems quite puz- zling, as the constraint of being LC was thought to be the main rea- son limiting the accuracy of the LC hybrid functionals, as the Coulomb tails have been shown to contribute insignificantly to rel- ative energies, such as atomization energies[33].
In this Letter, we try to unlock the above mystery by searching for a distinctive common feature of the optimized LR operators used for HF exchange mixing in xB97 and xB97X [14]. Upon uncovering such a feature, we then also compare to other proposed 0009-2614/$ - see front matter Published by Elsevier B.V.
doi:10.1016/j.cplett.2008.10.070
*Corresponding author.
E-mail addresses: [email protected](J.-D. Chai), [email protected] (M. Head-Gordon).
Chemical Physics Letters 467 (2008) 176–178
Contents lists available atScienceDirect
Chemical Physics Letters
j o u r n a l h o m e p a g e : w w w . e l s e v i e r . c o m / l o c a t e / c p l e t t
SR/LR partitions and comment on the prospects for future develop- ment of optimal LR operators.
By the SR/LR partition, Eq.(2)can be re-written as
ELCxc¼ ESR-DFAxc þ ðELR-HFx þ ELR-DFAc Þ: ð3Þ
A DFA XC hole is used for the SR operator, while the HF (exact) exchange hole and a DFA correlation hole are used for the LR oper- ator. Due to its semi-local nature, a DFA XC hole is known to be rea- sonably accurate for small inter-electronic separations from the reference electron, but becomes inaccurate for large inter-elec- tronic separations[34–37]. As a result, a DFA performs far better for the SR operator than does for the LR one. One thus expects only a small fraction of HF exchange is needed for the region where the DFA XC hole is already a good description of the exact XC hole.
Therefore, L(r) should be small for the SR region, and it should ap- proach unity asymptotically. The standard error function, erf(xr), provides a smooth transition between these two extremes, and is the preferred L(r) in many LC hybrid functionals, such asxB97 [14], because it facilitates efficient evaluation of the exact ex- change LR matrix elements over Gaussian basis functions[38].
Based on optimization over a large and diverse training set, the optimal parameter for the partition,x, is found to be 0.4 Bohr$1for xB97[14], in agreement with the work of Vydrov et al.[29,30]. We thus expect that erf(xr) is optimal atx= 0.4 Bohr$1for the LC hy- brid functionals with flexible and accurate GGAs.
To improve uponxB97, we have previously proposed to include a small fraction of ESR-HFx , as this should still be important for reduc- ing the remaining self-interaction-error of ESR-DFAxc [14]. Instead of using the erf operator, we effectively use the following LR operator (denoted as erfx) in Eq.(2)
erfxðr;
x
; cxÞ ¼ erfðx
rÞ þ cxerfcðx
rÞ ð4Þ where cxis a small fractional number.This erfx operator has been shown to be superior to the erf operator[14], as the root-mean-square (RMS) errors of the training set forxB97X (with erfx) are significantly reduced (by 0.5 kcal/
mol), when compared with xB97 (with erf). The RMS errors of the training set forxB97X optimized at different values ofxare plotted inFig. 1. Atx= 0.3 Bohr$1, the optimization is done self- consistently, while at other values ofx, non-self-consistent orbi- tals are used. At x= 0.0, 0.1, and 0.5 Bohr$1, the corresponding RSHXLDA orbitals [26] are used, and at x= 0.2 and 0.4 Bohr$1, thexB97X-D[15] (a re-optimizedxB97X with empirical atom–
atom dispersion corrections) and xB97 orbitals [14] are used respectively. We have previously demonstrated that these results should only change insignificantly with the self-consistent orbitals [14]. Another type of operator that effectively includes short-range Hartree–Fock exchange has been developed with similar argu- ments[39,40].
To investigate which part of the erfx operator is responsible for this improvement, the optimal erfx operators for xB97X[14] at different values ofx(0.2, 0.3, 0.4, and 0.5 Bohr$1) are plotted in Fig. 2. Surprisingly, all of the operators almost intersect at one point (r ’ 0.8 Bohr), which implies that the values of cx in this range ofxvalues are strongly dependent onx(indeed cxbecomes negativeatx= 0.5 Bohr$1). This also indicates that these operators must not change significantly in the middle-range (MR) region in order to achieve good balanced performance for thermochemistry and barrier heights.
As the RMS errors still show strongx-dependence, the fine de- tails of the MR region of the operators are clearly crucial. The diver- sity of behavior of these operators in the SR region (r [ 0.5 Bohr) and in the LR region (rJ 1.5 Bohr) implies that exact exchange in these two regions is relatively insignificant for thermochemistry and barrier heights. Similar conclusions for the importance of MR region of operators have also been presented by Henderson et al.
[41]using different arguments.
The exact exchange contribution from the SR region is expected to be unimportant for properties that do not involve changes in the core contributions to Exc, such as thermochemistry and barrier heights, where their effects on relative energies nearly cancel.
However, for properties sensitive to the core contributions, such as atomic energies and core excitation energies, the details of HF exchange mixing in the SR region could still be important.
The LR region of Coulomb operators contributes insignificantly to relative energies, like atomization energies, as demonstrated by Adamson et al.[33]. Of course, for the properties sensitive to these tail contributions, such as dissociation of cations with odd number of electrons[6–10,14,15]and long-range charge-transfer excitation energies[11–15], the details of exact exchange mixing in the LR region are crucial.
Finally, a GGA XC hole seem to be reasonably accurate even in the MR region (0.5 Bohr [r [ 1.5 Bohr), although less accurate than in the SR region. From the values of the optimal erfx operator used inxB97X, only 30% HF exchange is needed at r = 0.5 Bohr, and 56% HF exchange is needed at r = 1.5 Bohr, which reflects the decreasing appropriateness of the underlying GGA XC hole with in- ter-electronic distance. We infer that the constraint of having zero
2.6 2.8 3 3.2 3.4 3.6 3.8 4 4.2
0 0.1 0.2 0.3 0.4 0.5
RMS Error (kcal/mol)
ω (Bohr-1)
Fig. 1. The root-mean-square (RMS) errors of the training set forxB97X optimized at different values of x. At x= 0.3 Bohr$1, the optimization is done self- consistently, while at other values ofx, non-self-consistent orbitals are used (see text).
-0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 2.2 2.4 2.6 2.8 3
erfx(r)
r (Bohr)
(ω=0.2, cx= 0.238794) (ω=0.3, cx= 0.157706) (ω=0.4, cx= 0.044520) (ω=0.5, cx=-0.099103)
Fig. 2. The optimal erfx operators forxB97X[14]at different values ofx.
J.-D. Chai, M. Head-Gordon / Chemical Physics Letters 467 (2008) 176–178 177
slope in the operator used for the HF exchange mixing in a global hybrid functional results in its inferior performance to a corre- sponding LC hybrid functional.
To test the transferability of the above observations and argu- ments, and demonstrate the importance of the MR region of the operator used for exact exchange, we consider another LR operator, the erfgau operator[32,42]:
erfgauðr;
x
; k; aÞ ¼ erfðx
rÞ $ k2x
rffiffiffiffi
p e
p
$ð1=aÞx2r2 ð5ÞThe erfgau operator has a sharper separation (with positive k and a) between the LR and SR region than the erf operator. In the work of Toulouse et al.[32], k = 1 and a = 3 are chosen to ensure the value of erfgau(r)/r and its first derivative vanish at r = 0. Re- cently, Song et al. have found that this operator gives poor results for atomization energies and barrier heights[31], with an opti- mizedx(=2.6 Bohr$1). Their results are greatly improved by opti- mizing all the three parameters (x, k, and a) of the erfgau operator.
They have found that a negative k is essential for good atomization energies and barrier heights. To see the significance of their find- ings, the three optimal erfgau operators in their work are plotted in Fig. 3, and are compared with the optimal erf operator used forxB97, as well as the optimal erfx operators used forxB97X [14], and for the corresponding global hybrid (atx= 0.0 Bohr$1).
As can be seen, the original erfgau(x= 2.60 Bohr$1, k = 1.00, a= 3.00) has a very sharp transition in the MR region, while their newly proposed erfgau(x= 0.39 Bohr$1, k = $16.4, and a = 0.013) and erfgau(x= 0.42 Bohr$1, k = $18.0, and a = 0.011) resemble the optimal erf operator ofxB97[14]in that region, which should be responsible for their improvements upon the original erfgau operator. In the SR and LR regions, all operators differ significantly, which supports our arguments that these are relatively unimpor- tant regions for thermochemistry and barrier heights.
Of course, functionals modeled on other types of LR operators can also be considered [43], and we believe that the MR region of operators is the most important. Interestingly, it has been previ- ously suggested that Coulomb-attenuated functionals cannot be simultaneously optimized for good thermochemistry and charge- transfer excitation energies[44]. Since we have shown that the for- mer is controlled by the MR region, and the latter is likely con- trolled by the LR region, it seems worthwhile to revisit this question in the future.
In conclusion, we have found a conserved feature of the optimal erfx operators forxB97X at different values ofx. On this basis, we
argue that the MR region (r = 0.5 to 1.5 Bohr) of the partitioning operators seems to be crucial for thermochemistry and barrier heights. The fine details of the MR region of the operators for the HF exchange mixing in the LC hybrid functionals are responsible for their improved performance over the corresponding global hy- brid functionals. The MR regions of the optimal LR operators of xB97 andxB97X[14], erf(x= 0.4 Bohr$1) and erfx(x= 0.3 Bohr$1, cx= 0.157706) respectively, can be used an initial criterion for devising new and improved LR operators. As the LR region is found to be unimportant for thermochemistry and barrier heights, one may devise new LR operators that approach to unity faster (better for the LR properties) than the optimal erf and erfx operators, with- out reducing the accuracy in thermochemistry and barrier heights.
We are now investigating this possibility.
Acknowledgement
This work was supported by the US Department of Energy through the Theory and Modeling in Nanoscience Program. MHG is a part-owner of Q-Chem Inc.
References
[1] P. Hohenberg, W. Kohn, Phys. Rev. 136 (1964) B864.
[2] W. Kohn, L.J. Sham, Phys. Rev. 140 (1965) A1133.
[3] R.G. Parr, W. Yang, in: Density-Functional Theory of Atoms and Molecules, Oxford University Press, 1989.
[4] R.M. Dreizler, E.K.U. Gross, in: Density Functional Theory: An Approach to the Quantum Many Body Problem, Springer-Verlag, Berlin, 1990.
[5] W. Kohn, A.D. Becke, R.G. Parr, J. Phys. Chem. 100 (1996) 12974.
[6] T. Bally, G.N. Sastry, J. Phys. Chem. A 101 (1997) 7923.
[7] B. Bra, P.C. Hiberty, A. Savin, J. Phys. Chem. A 102 (1998) 7872.
[8] P. Mori-Sánchez, A.J. Cohen, W. Yang, J. Chem. Phys. 125 (2006) 201102.
[9] A. Ruzsinszky, J.P. Perdew, G.I. Csonka, O.A. Vydrov, G.E. Scuseria, J. Chem. Phys.
126 (2007) 104102.
[10] A.D. Dutoi, M. Head-Gordon, Chem. Phys. Lett. 422 (2006) 230.
[11] A. Dreuw, J.L. Weisman, M. Head-Gordon, J. Phys. Chem. 119 (2003) 2943.
[12] A. Dreuw, M. Head-Gordon, J. Am. Chem. Soc. 126 (2004) 4007.
[13] A. Dreuw, M. Head-Gordon, Chem. Rev. 105 (2005) 4009.
[14] J.-D. Chai, M. Head-Gordon, J. Chem. Phys. 128 (2008) 084106.
[15] J.-D. Chai, M. Head-Gordon, Phys. Chem. Chem. Phys. 10 (2008) 6615.
[16] A.D. Becke, J. Chem. Phys. 98 (1993) 5648.
[17] H. Stoll, A. Savin, in: R.M. Dreizler, J.D. Providencia (Eds.), Density Functional Methods in Physics, Plenum, New York, 1985, p. 177.
[18] T. Leininger, H. Stoll, H.-J. Werner, A. Savin, Chem. Phys. Lett. 275 (1997) 151.
[19] J. Toulouse, F. Colonna, A. Savin, J. Chem. Phys. 122 (2005) 014110.
[20] J.G. Ángyán, I.C. Gerber, A. Savin, J. Toulouse, Phys. Rev. A 72 (2005) 012510.
[21] E. Goll, H.-J. Werner, H. Stoll, Phys. Chem. Chem. Phys. 7 (2005) 3917.
[22] E. Goll, H.-J. Werner, H. Stoll, T. Leininger, P. Gori-Giorgi, A. Savin, Chem. Phys.
329 (2006) 276.
[23] M.A. Rohrdanz, J.M. Herbert, J. Chem. Phys. 129 (2008) 034107.
[24] A. Savin, in: J.M. Seminario (Ed.), Recent Developments and Applications of Modern Density Functional Theory, Elsevier, Amsterdam, 1996, p. 327.
[25] P.M.W. Gill, R.D. Adamson, J.A. Pople, Mol. Phys. 88 (1996) 1005.
[26] I.C. Gerber, J.G. Ángyán, Chem. Phys. Lett. 415 (2005) 100.
[27] H. Iikura, T. Tsuneda, T. Yanai, K. Hirao, J. Chem. Phys. 115 (2001) 3540.
[28] J.-W. Song, T. Hirosawa, T. Tsuneda, K. Hirao, J. Chem. Phys. 126 (2007) 154105.
[29] O.A. Vydrov, J. Heyd, A.V. Krukau, G.E. Scuseria, J. Chem. Phys. 125 (2006) 074106.
[30] O.A. Vydrov, G.E. Scuseria, J. Chem. Phys. 125 (2006) 234109.
[31] J.-W. Song, S. Tokura, T. Sato, M.A. Watson, K. Hirao, J. Chem. Phys. 127 (2007) 154109.
[32] J. Toulouse, F. Colonna, A. Savin, Phys. Rev. A 70 (2004) 062505.
[33] R.D. Adamson, J.P. Dombroski, P.M.W. Gill, Chem. Phys. Lett. 254 (1996) 329.
[34] J.P. Perdew, M. Ernzerhof, K. Burke, A. Savin, Int. J. Quantum Chem. 61 (1997) 197.
[35] K. Burke, J.P. Perdew, M. Ernzerhof, J. Chem. Phys. 109 (1998) 3760.
[36] A.C. Cancio, C.Y. Fong, J.S. Nelson, Phys. Rev. A 62 (2000) 062507.
[37] J.P. Perdew, A. Ruzsinszky, J. Tao, V.N. Staroverov, G.E. Scuseria, G.I. Csonka, J.
Chem. Phys. 123 (2005) 062201.
[38] R.D. Adamson, J.P. Dombroski, P.M.W. Gill, J. Comp. Chem. 20 (1999) 921.
[39] T. Yanai, D.P. Tew, N.C. Handy, Chem. Phys. Lett. 393 (2004) 51.
[40] A.J. Cohen, P. Mori-Sánchez, W. Yang, J. Chem. Phys. 126 (2007) 191109.
[41] T.M. Henderson, A.F. Izmaylov, G.E. Scuseria, A. Savin, J. Chem. Phys. 127 (2007) 221103.
[42] P.M.W. Gill, R.D. Adamson, Chem. Phys. Lett. 261 (1996) 105.
[43] A.D. Dutoi, M. Head-Gordon, J. Phys. Chem. A 112 (2008) 2110.
[44] M.J.G. Peach, A.J. Cohen, D.J. Tozer, Phys. Chem. Chem. Phys. 8 (2006) 4543.
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 2.2 2.4 2.6 2.8 3
L(r)
r (Bohr)
erfgau (ω=2.60, k = 1.00, a= 3.00) erfgau (ω=0.39, k = -16.4, a= 0.013) erfgau (ω=0.42, k = -18.0, a= 0.011) erf (ω=0.4) erfx (ω=0.3, cx=0.157706) erfx (ω=0.0, cx=0.303766)
Fig. 3. The three optimal erfgau operators in Ref.[31]are compared with the optimal erf operator used forxB97, as well as the optimal erfx operators used for xB97X[14], and for the corresponding global hybrid (atx= 0.0 Bohr$1).
178 J.-D. Chai, M. Head-Gordon / Chemical Physics Letters 467 (2008) 176–178
