J. Chem. Phys. 146, 044102 (2017); https://doi.org/10.1063/1.4974163 146, 044102
© 2017 Author(s).
Role of exact exchange in thermally- assisted-occupation density functional
theory: A proposal of new hybrid schemes
Cite as: J. Chem. Phys. 146, 044102 (2017); https://doi.org/10.1063/1.4974163
Submitted: 06 October 2016 . Accepted: 04 January 2017 . Published Online: 23 January 2017 Jeng-Da Chai
ARTICLES YOU MAY BE INTERESTED IN
Density functional theory with fractional orbital occupations
The Journal of Chemical Physics 136, 154104 (2012); https://doi.org/10.1063/1.3703894 Perspective: Kohn-Sham density functional theory descending a staircase
The Journal of Chemical Physics 145, 130901 (2016); https://doi.org/10.1063/1.4963168 Long-range interactions from the many-pair expansion: A different avenue to dispersion in DFT
The Journal of Chemical Physics 146, 024111 (2017); https://doi.org/10.1063/1.4973728
Role of exact exchange in thermally-assisted-occupation density functional theory: A proposal of new hybrid schemes
Jeng-Da Chaia)
Department of Physics, Center for Theoretical Sciences, and Center for Quantum Science and Engineering, National Taiwan University, Taipei 10617, Taiwan
(Received 6 October 2016; accepted 4 January 2017; published online 23 January 2017)
We propose hybrid schemes incorporating exact exchange into thermally assisted-occupation-density functional theory (TAO-DFT) [J.-D. Chai, J. Chem. Phys. 136, 154104 (2012)] for an improved description of nonlocal exchange effects. With a few simple modifications, global and range-separated hybrid functionals in Kohn-Sham density functional theory (KS-DFT) can be combined seamlessly with TAO-DFT. In comparison with global hybrid functionals in KS-DFT, the resulting global hybrid functionals in TAO-DFT yield promising performance for systems with strong static correlation effects (e.g., the dissociation of H2 and N2, twisted ethylene, and electronic properties of linear acenes), while maintaining similar performance for systems without strong static correlation effects.
Besides, a reasonably accurate description of noncovalent interactions can be efficiently achieved through the inclusion of dispersion corrections in hybrid TAO-DFT. Relative to semilocal density functionals in TAO-DFT, global hybrid functionals in TAO-DFT are generally superior in performance for a wide range of applications, such as thermochemistry, kinetics, reaction energies, and optimized geometries. Published by AIP Publishing. [http://dx.doi.org/10.1063/1.4974163]
I. INTRODUCTION
Over the past two decades, Kohn-Sham density functional theory (KS-DFT)1,2 has emerged as one of the most popular electronic structure methods for the study of large ground- state systems, due to its low computational cost and reasonable accuracy.3–6Nevertheless, the essential ingredient of KS-DFT, the exact exchange-correlation (XC) energy functional Exc[ ρ]
remains unknown, and needs to be approximated. Conse- quently, density functional approximations (DFAs) for Exc[ ρ]
have been continuously developed to improve the accuracy of KS-DFT for a broad range of applications.
Functionals based on the conventional semilocal DFAs, such as the local density approximation (LDA)7,8and gener- alized gradient approximations (GGAs),9–11can yield reason- ably accurate predictions of the properties governed by short- range XC effects, and possess high computational efficiency for very large systems (for brevity, hereafter we use “DFAs”
for “the conventional semilocal DFAs”). Nonetheless, owing to the inappropriate treatment of nonlocal XC effects,12,13 KS-DFAs can perform very poorly in situations where the self-interaction error (SIE),12–15noncovalent interaction error (NCIE),16–18 or static correlation error (SCE)12,19–22 is pro- nounced. Over the years, considerable efforts have been made to resolve the qualitative failures of KS-DFAs at a reasonable computational cost.
To date, global hybrid functionals23–30 and range- separated hybrid functionals,31–33 which incorporate the Hartree-Fock (HF) exchange energy into KS-DFAs, are per- haps the most successful schemes that provide an improved
a)Electronic mail: jdchai@phys.ntu.edu.tw
description of nonlocal exchange effects. Relative to KS-DFAs, the hybrid schemes, which greatly reduce the SIE problems, are reliably accurate for a wide variety of applications, such as thermochemistry and kinetics.34,35
To properly describe noncovalent interactions, a reason- ably accurate treatment of middle- and long-range dynami- cal correlation effects is critical. Accordingly, KS-DFAs and hybrid functionals may be combined with the DFT-D (KS-DFT with empirical dispersion corrections) schemes17,36–40and the double-hybrid (mixing both the HF exchange energy and the second-order Møller-Plesset (MP2) correlation energy41into KS-DFAs) schemes,30,33,42 showing an overall satisfactory accuracy for the NCIE problems.
In spite of their computational efficiency, KS-DFAs, hybrid functionals, and double-hybrid functionals can per- form very poorly for systems with strong static correlation effects (i.e., multi-reference systems).12,19–22Within KS-DFT, fully nonlocal XC functionals, such as those based on the random phase approximation (RPA), may be adopted for a reliably accurate description of strong static correlation effects.
However, RPA-type functionals remain computationally very demanding for large systems.4,13,43,44
To reduce the SCE problems with low computational complexity, we have recently developed thermally-assisted- occupation density functional theory (TAO-DFT),20,21 an efficient electronic structure method for studying the ground- state properties of very large systems (e.g., containing up to a few thousand electrons) with strong static correlation effects.45–48 Unlike finite-temperature DFT,49 TAO-DFT is developed for ground-state systems at zero temperature. In contrast to KS-DFT, TAO-DFT is a DFT with fractional orbital occupations given by the Fermi-Dirac distribution (con- trolled by a fictitious temperature θ), wherein strong static
0021-9606/2017/146(4)/044102/20/$30.00 146, 044102-1 Published by AIP Publishing.
correlation is explicitly described by the entropy contribu- tion (e.g., see Eq. (26) of Ref.20). Interestingly, TAO-DFT is as efficient as KS-DFT for single-point energy and analyti- cal nuclear gradient calculations, and is reduced to KS-DFT in the absence of strong static correlation effects. Besides, exist- ing DFA XC functionals in KS-DFT may also be adopted in TAO-DFT. The resulting TAO-DFAs have been shown to consistently improve upon KS-DFAs for multi-reference systems. Nevertheless, TAO-DFAs perform similarly to KS- DFAs for single-reference systems (i.e., systems without strong static correlation). In addition, the SIEs and NCIEs of TAO-DFAs may remain enormous in situations where these failures occur.
In this work, we aim to improve the accuracy of TAO- DFAs for a wide variety of single-reference systems. Specif- ically, we develop hybrid schemes that incorporate exact exchange into TAO-DFAs for an improved description of non- local exchange effects. Hybrid functionals (e.g., global and range-separated hybrids) in KS-DFT can be easily modified, and seamlessly combined with TAO-DFT. The rest of the paper is organized as follows. A brief review of the essen- tials of TAO-DFT is provided in SectionII. In SectionIII, the exact exchange in TAO-DFT is defined, and the corresponding global and range-separated hybrid schemes are proposed. In SectionIV, the optimal θ values for global hybrid function- als in TAO-DFT are defined, and the performance of global hybrid functionals in TAO-DFT (with the optimal θ values) is examined for various single- and multi-reference systems.
Our conclusions are given in SectionV.
II. TAO-DFT
A. Rationale for fractional orbital occupations
Consider an interacting N-electron Hamiltonian for an external potential v(r) at zero temperature, the exact ground- state density ρ(r) is interacting v-representable, as it can be obtained from the ground-state wavefunction calculated using the full configuration interaction (FCI) method at the complete basis set limit50
ρ(r) =
∞
X
i=1
ni|χi(r)|2, (1) which can be expressed in terms of the natural orbitals (NOs) {χi(r)} and natural orbital occupation numbers (NOONs) {ni} (i.e., the eigenfunctions and eigenvalues, respectively, of one-electron reduced density matrix (1-RDM)).51 Here, the NOONs {ni}, obeying the following two conditions:
∞
X
i=1
ni= N, 0 ≤ ni≤1, (2) are related to the variationally determined coefficients of the FCI expansion. As shown in Eq.(1), the exact ground-state density ρ(r) can be represented by orbitals and their occupation numbers, showing the significance of an ensemble represen- tation (via fractional orbital occupations) of the ground-state density.
By contrast, in KS-DFT, the ground-state density ρ(r) is assumed to be noninteracting pure-state vs-representable,
as it belongs to a one-determinant ground-state wavefunc- tion of a noninteracting N-electron Hamiltonian for some local potential vs(r) at zero temperature.52–54Accordingly, the Kohn-Sham (KS) orbital occupation numbers should be either 0 or 1. Due to the search over a restricted domain of densi- ties, some ground-state densities cannot be obtained within the framework of KS-DFT (i.e., even with the exact Exc[ ρ]).53–57 Baerends and co-workers55argued that the ground-state den- sity ρ(r) of a system with strong static correlation effects may not be noninteracting pure-state vs-representable, wherein an ensemble representation of the ground-state density is essen- tial. Arguments supporting this are also available from other studies.56,58
To rectify the above situation, KS-DFT has been extended to ensemble DFT,59,60wherein ρ(r) is assumed to be nonin- teracting ensemble vs-representable, as it is associated with an ensemble of pure determinantal states of the noninteract- ing KS system at zero temperature. Accordingly, the orbital occupation numbers in ensemble DFT are 0, 1, and fractional (between 0 and 1) for the orbitals above, below, and at the Fermi level, respectively. Within the framework of ensem- ble DFT, the development of DFT fractional-occupation- number (DFT-FON) method,58,61–64spin-restricted ensemble- referenced KS (REKS) method,65,66and fractional-spin DFT (FS-DFT) method12,19 has yielded great success for some systems with strong static correlation effects. Nevertheless, the practical implementation of DFT-FON and related meth- ods has been hindered by several factors, such as a possible double-counting of correlation effects and the sharp increase of computational cost for large systems.
On the other hand, the inclusion of fractional occu- pation numbers (FONs) in electronic structure calculations has a long history.12,19–21,49,58–83 In particular, the Fermi- Dirac distribution, which appears in finite-temperature DFT49 and finite-temperature HF schemes,72,76–83 has been a pop- ular distribution function for the FON-related schemes. For example, finite-temperature techniques have been developed for improving self-consistent field (SCF) convergence.70The grand canonical orbitals have been used for subsequent com- plete active space configuration interaction (CASCI) calcula- tions.71–73Recently, a fractional occupation number weighted electron density has been adopted for a real-space measure and visualization of static correlation effects.77
In TAO-DFT,20,21 the representation of the ground-state density from the exact theory (see Eq. (1)) has been high- lighted. In contrast to the orbital occupation numbers in KS-DFT and ensemble DFT, the NOONs can be fractional (between 0 and 1) for all the NOs. While the exact NOONs are intractable for large systems (due to the exponential complex- ity), the distribution of NOONs (the microcanonical averaging of NOONs) can, however, be approximately described by the Fermi-Dirac distribution with renormalized parameters (i.e., orbital energies, chemical potential, and temperature) based on the statistical arguments of Flambaum et al.84Accordingly, in TAO-DFT, the ground-state density ρ(r) of a system of N inter- acting electrons moving in an external potential vext(r) at zero temperature is assumed to be noninteracting thermal ensemble vs-representable, as it is expressed as the thermal equilibrium density of an auxiliary system of N noninteracting electrons
moving in some local potential vs(r) at a fictitious temperature θ. Consequently, ρ(r) can be represented by
ρ(r) =
∞
X
i=1
fi|ψi(r)|2, (3) where the orbital occupation number fi is the Fermi-Dirac distribution
fi= {1 + exp[(i−µ)/θ]}−1, (4) which satisfies the following two conditions:
∞
X
i=1
fi= N, 0 ≤ fi ≤1, (5)
iis the orbital energy of the ith orbital ψi(r), and µ is the chem- ical potential determined by the conservation of the number of electrons N.
As discussed in Ref.20, for a given fictitious tempera- ture θ, the Hohenberg-Kohn theorems1and the Mermin theo- rems49can be employed for the physical and auxiliary systems, respectively, to derive a set of self-consistent equations in TAO-DFT for determining the remaining “renormalized parameters” (i.e., the orbital energies {i}and chemical poten- tial µ) of the orbital occupation numbers { fi} and the orbitals {ψi(r)}, which can then be used to represent the ground-state density ρ(r), and evaluate the ground-state energy of the phys- ical system at zero temperature. In addition, due to the simi- larity of Eqs.(1)and(3), when the fictitious temperature θ in TAO-DFT is so chosen that the NOONs {ni} are approxi- mately described by the orbital occupation numbers { fi} (in the sense of statistical average, as mentioned above), the NOs {χi(r)} will be approximately described by the orbitals {ψi(r)}.
This implies that the exact ρ(r) is likely to be noninteracting thermal ensemble vs-representable at this θ value (plus some range of possible other values around it). In addition, as dis- cussed in Ref.20, strong static correlation has been shown to be properly described by the entropy contribution (e.g., see Eq. (26) of Ref.20) in TAO-DFT at this θ value (plus some range of possible other values around it).
While also adopting the Fermi-Dirac distribution, TAO- DFT is developed for the ground-state density and ground-state energy of a physical system at zero temperature, which is different from the aforementioned finite-temperature FON- related schemes (which mostly focus on the SCF conver- gence, the adoption of grand canonical orbitals and density for different purposes, and the thermodynamic properties of a physical system at finite temperature). On the other hand, while KS-DFT, ensemble DFT, and TAO-DFT all belong to zero-temperature DFT, the representations of the ground-state density are, however, different in these methods (as mentioned above). While the entropy contribution in TAO-DFT plays an important role in simulating strong static correlation (even though at the price of adding an extra θ parameter that is related to the distribution of NOONs), this term is, however, absent in KS-DFT and ensemble DFT.
B. Self-consistent equations
Consider a system of Nα up-spin and Nβ down-spin electrons moving in an external potential vext(r) at zero
(physical) temperature. In spin-polarized (spin-unrestricted) TAO-DFT,20,21 two noninteracting reference systems at the same fictitious (reference) temperature θ (measured in energy units) are employed: one described by the spin function α and the other described by the spin function β, with the corre- sponding thermal equilibrium density distributions ρs,α(r) and ρs,β(r) exactly equal to the up-spin density ρα(r) and down- spin density ρβ(r), respectively, in the original interacting system at zero temperature. The resulting self-consistent equa- tions for the σ-spin electrons (σ = α or β) can be expressed as (i runs for the orbital index)
−1
2∇2 + vs,σ(r)
ψiσ(r)= iσψiσ(r), (6) where
vs,σ(r)= vext(r)+
ρ(r0)
|r − r0|dr0+δExc[ ρα, ρβ]
δ ρσ(r) +δEθ[ ρα, ρβ] δ ρσ(r)
(7) is the effective potential (atomic units, i.e., ~ = me= e = 4πε0
= 1, are adopted throughout this work). Here, Exc[ ρα, ρβ]
≡ Ex[ ρα, ρβ]+ Ec[ ρα, ρβ] is the XC energy (i.e., the sum of the exchange energy Ex[ ρα, ρβ] and correlation energy Ec[ ρα, ρβ]) defined in spin-polarized KS-DFT,58,85 and Eθ[ ρα, ρβ] ≡ Aθ=0s [ ρα, ρβ] − Asθ[ ρα, ρβ] is the difference between the noninteracting kinetic free energy at zero tem- perature and that at the fictitious temperature θ. The σ-spin density
ρσ(r)=
∞
X
i=1
fiσ|ψiσ(r)|2 (8) is expressed in terms of the thermally-assisted-occupation (TAO) orbitals {ψiσ(r)} and their occupation numbers { fiσ},
fiσ= {1 + exp[(iσ−µσ)/θ]}−1, (9) which are given by the Fermi-Dirac distribution. Here, the chemical potential µσ is determined by the conservation of the number of σ-spin electrons Nσ,
∞
X
i=1
{1+ exp[(iσ−µσ)/θ]}−1= Nσ. (10) The two sets (one for each spin function) of self-consistent equations, Equations (6)–(10), for ρα(r) and ρβ(r), respec- tively, are coupled with the ground-state density
ρ(r) =
α,β
X
σ
ρσ(r). (11)
The self-consistent procedure described in Ref.20may be employed to obtain ρσ(r) and ρ(r). After self-consistency is achieved, the noninteracting kinetic free energy
Aθs[{fiα, ψiα}, {fiβ, ψiβ}]= Tsθ[{fiα, ψiα}, {fiβ, ψiβ}] + ESθ[{fiα}, {fiβ}] (12) can be computed, in an exact manner, as the sum of the kinetic energy
Tsθ[{fiα, ψiα}, {fiβ, ψiβ}]= −1 2
α,β
X
σ
∞
X
i=1
fiσ
ψ∗iσ(r)∇2ψiσ(r)dr (13)
and entropy contribution ESθ[{fiα}, {fiβ}]= θ
α,β
X
σ
∞
X
i=1
fiσln(fiσ)+ (1 − fiσ) ln(1 − fiσ)
(14) of noninteracting electrons at the fictitious temperature θ. The ground-state energy of the original interacting system at zero temperature is given by
E[ ρα, ρβ]= Aθs[{fiα, ψiα}, {fiβ, ψiβ}] +
ρ(r)vext(r)dr+ EH[ ρ]+ Exc[ ρα, ρβ]
+ Eθ[ ρα, ρβ], (15)
where EH[ ρ] ≡ 12 ∫ ∫ ρ(r)ρ(r
0)
|r − r0|drdr0 is the Hartree energy.
Spin-unpolarized (spin-restricted) TAO-DFT can be formu- lated by imposing the constraints of ψiα(r) = ψiβ(r) and fiα
= fiβto spin-polarized TAO-DFT.
C. Density functional approximations
As the exact Exc[ ρα, ρβ] and Eθ[ ρα, ρβ] (i.e., the essen- tial ingredients of spin-polarized TAO-DFT) remain unknown, DFAs for both of them (denoted as TAO-DFAs) are necessary for practical applications. Consequently, the performance of TAO-DFAs depends on the accuracy of DFAs and the choice of the fictitious temperature θ. Note that EDFAxc [ ρα, ρβ] can be readily obtained from that of KS-DFA, and EθDFA[ ρα, ρβ] can be obtained with the knowledge of ADFA,θs [ ρα, ρβ] as follows
EθDFA[ ρα, ρβ] ≡ ADFA,θs =0[ ρα, ρβ] − ADFA,θs [ ρα, ρβ]
= 1
2(ADFA,θs =0[2 ρα]+ ADFA,θs =0[2 ρβ])
−1
2(ADFA,θs [2 ρα]+ ADFA,θs [2 ρβ]), (16) where ADFA,θs [ ρα, ρβ] is expressed in terms of ADFA,θs [ ρ] (in its spin-unpolarized form) based on the spin-scaling relation of Aθs[ ρα, ρβ].86 Note that EDFAθ=0[ ρα, ρβ] = 0 (i.e., an exact property of Eθ[ ρα, ρβ]) is ensured by Eq.(16). Accordingly, TAO-DFAs at θ= 0 reduce to KS-DFAs.
D. Strong static correlation from TAO-DFAs
In 2012, we developed TAO-LDA,20employing the LDA XC functional ExcLDA[ ρα, ρβ]7,8 and EθLDA[ ρα, ρβ] (given by Eq.(16)with ALDA,θs [ ρ], the LDA for Aθs[ ρ] (see Appendix A of Ref. 87 and Eq. (37) of Ref.20)) in TAO-DFT. Even at the simplest LDA level, TAO-LDA was shown to provide a reasonably accurate treatment of static correlation via the entropy contribution ESθ[{fiα}, {fiβ}] (see Eq.(14)), when the distribution of TAO orbital occupation numbers (TOONs) {fiσ} (related to the chosen θ) is close to the distribution of NOONs. However, this implies that a θ related to the distribu- tion of NOONs should be employed to properly describe strong static correlation effects. For simplicity, an optimal value of θ = 7 mhartree was previously defined for TAO-LDA, based on physical arguments and numerical investigations. TAO-LDA (with θ = 7 mhartree) was shown to consistently outperform
KS-LDA for multi-reference systems (due to the appropriate treatment of static correlation), while performing comparably to KS-LDA for single-reference systems (i.e., in the absence of strong static correlation effects).
To improve the accuracy of TAO-LDA for single- reference systems, in 2014, we developed TAO-GGAs,21 adopting the GGA XC functionals EGGAxc [ ρα, ρβ] and EθGEA[ ρα, ρβ] (given by Eq.(16)with AGEA,θs [ ρ], the gradient expansion approximation (GEA) for Asθ[ ρ] (see Appendices A and B of Ref. 87)) in TAO-DFT. As TAO-GGAs should improve upon TAO-LDA mainly for the properties governed by short-range XC effects (due to the more accurate treat- ment of on-top hole density),3–6,13 and the orbital energy gaps of TAO-LDA and TAO-GGAs should be similar,34,35the optimal θ values for TAO-LDA and TAO-GGAs should remain similar (when the same physical arguments and numerical investigations are adopted to define the optimal θ values).
Therefore, we adopted an optimal value of θ = 7 mhartree for both TAO-LDA and TAO-GGAs. While EθGEA[ ρα, ρβ] should be more accurate than ELDAθ [ ρα, ρβ] for the nearly uniform electron gas, for a small value of θ (i.e., 7 mhartree), their difference was found to be much smaller than the difference between two different XC energy functionals. Unsurprisingly, since Eθ=0LDA = EGEAθ=0 = 0, the difference between EθLDA and EθGEA should remain small for a small value of θ (i.e., 7 mhartree). Accordingly, ELDAθ [ ρα, ρβ] may also be adopted for TAO-GGAs.
While TAO-DFAs (i.e., TAO-LDA and TAO-GGAs) out- perform KS-DFAs for multi-reference systems, they perform similarly to KS-DFAs for single-reference systems. As men- tioned previously, hybrid functionals in KS-DFT, which pro- vide an improved description of nonlocal exchange effects, are generally superior to KS-DFAs in performance for a broad range of applications.34,35Therefore, a possible hybrid func- tional in TAO-DFT is expected to outperform TAO-DFAs for a wide variety of single-reference systems. In Sec. III, we define the exact exchange in TAO-DFT, and propose the corresponding global and range-separated hybrid schemes in TAO-DFT.
III. HYBRID SCHEMES IN TAO-DFT A. Exact exchange
In KS-DFT, the exact exchange Ex[ ρα, ρβ] is defined as the HF exchange energy of the occupied KS orbitals {φiσ(r)},3–6
Ex[ ρα, ρβ] ≡ ExHF[{φiα}, {φiβ}]= −1 2
α,β
X
σ Nσ
X
i,j=1
× φ∗iσ(r1)φ∗jσ(r2)φjσ(r1)φiσ(r2) r12 dr1dr2
= −1 2
α,β
X
σ
|γσKS(r1, r2)|2
r12 dr1dr2, (17) where r12 = |r1 − r2| is the interelectronic distance. Here, γσKS(r1, r2) = PNi=1σφ∗iσ(r1)φiσ(r2) is the σ-spin 1-RDM in KS-DFT, and its diagonal element γKSσ (r, r)= PNi=1σ|φiσ(r)|2
= ρσ(r) is the σ-spin density in KS-DFT.
In TAO-DFT, the exact exchange Fxθ[ ρα, ρβ] can be defined as the HF exchange free energy of the TAO orbitals {ψiσ(r)} and their occupation numbers {fiσ}at the fictitious temperature θ,
Fxθ[ ρα, ρβ] ≡ FxHF,θ[{fiα, ψiα}, {fiβ, ψiβ}]= −1 2
α,β
X
σ
∞
X
i,j=1
fiσfjσ
×
ψ∗iσ(r1)ψ∗jσ(r2)ψjσ(r1)ψiσ(r2) r12 dr1dr2
= −1 2
α,β
X
σ
|γTAOσ (r1, r2)|2
r12 dr1dr2. (18) Here,
γTAOσ (r1, r2)=
∞
X
i=1
fiσψiσ∗ (r1)ψiσ(r2) (19) is the σ-spin 1-RDM in TAO-DFT, and its diagonal element
γTAOσ (r, r)=
∞
X
i=1
fiσ|ψiσ(r)|2=
∞
X
i=1
ρiσ(r)= ρσ(r) (20)
is the σ-spin density in TAO-DFT, where ρiσ(r) ≡ fiσ|ψiσ(r)|2 is the ith σ-spin orbital density in TAO-DFT. Note that the TAO orbitals {ψiσ(r)} and their occupation numbers {fiσ}are the eigenfunctions and eigenvalues, respectively, of γσTAO(r1, r2),
γTAOσ (r1, r2)ψiσ(r1)dr1=
∞
X
j=1
fjσψjσ(r2)
×
ψjσ∗ (r1)ψiσ(r1)dr1
=
∞
X
j=1
fjσψjσ(r2)δij= fiσψiσ(r2), (21) where δijis the Kronecker delta function. At θ= 0, TAO-DFT is the same as KS-DFT, and hence, Eq. (18)is reduced to Eq.(17).
To justify the use of Fxθ[ ρα, ρβ] ≡ FxHF,θ[{fiα, ψiα}, {fiβ, ψiβ}] (given by Eq. (18)) as the definition of exact exchange in TAO-DFT, here, we comment on the self- interaction energy associated with the exact exchange in TAO-DFT. On the basis of Equations(8)and(11), the Hartree energy can be expressed as
EH[ ρ] ≡ 1 2
ρ(r1) ρ(r2) r12 dr1dr2
=1 2
α,β
X
σ α,β
X
σ0
ρσ(r1) ρσ0(r2) r12 dr1dr2
=1 2
α,β
X
σ α,β
X
σ0
∞
X
i,j=1
fiσfjσ0
×
|ψiσ(r1)|2|ψjσ0(r2)|2
r12 dr1dr2. (22) Accordingly, the self-Hartree energy,14
self-Hartree energy ≡
α,β
X
σ
∞
X
i=1
EH[ ρiσ]=1 2
α,β
X
σ
∞
X
i=1
fiσ2
×
|ψiσ(r1)|2|ψiσ(r2)|2 r12 dr1dr2,
(23) which is the sum of the (σ = σ0and i = j) terms in Eq.(22), can be exactly cancelled by the self-exchange energy,14
self-exchange energy ≡
α,β
X
σ
∞
X
i=1
Fxθ[ ρiσ, 0]= −1 2
α,β
X
σ
∞
X
i=1
fiσ2
×
|ψiσ(r1)|2|ψiσ(r2)|2 r12 dr1dr2
= −
α,β
X
σ
∞
X
i=1
EH[ ρiσ], (24) which is the sum of the (i = j) terms in Eq. (18), on an orbital-by-orbital basis (i.e., each term in Eq. (23) can be exactly cancelled by a term in Eq. (24)). Therefore, com- plete cancellation of the self-interaction in the Hartree energy would require the full exact exchange (given by Eq.(18)) in TAO-DFT. By contrast, such perfect cancellation may not be achieved by the HF exchange (given by Eq. (17)) with the KS orbitals being replaced by the TAO orbitals. Besides, the self-Hartree energy in TAO-DFT is unlikely to be exactly can- celled by the self-XC energy, i.e., Pα,β
σ P∞
i=1ExcDFA[ ρiσ, 0], associated with the DFA XC functional ExcDFA[ ρα, ρβ], implying that the SIEs associated with TAO-DFAs may remain pronounced for both single- and multi-reference systems!
From Eq. (17), EHFx [{φiα}, {φiβ}] (i.e., the exact exchange in KS-DFT) can be expressed as
ExHF[{φiα}, {φiβ}]= Ex[ ρα, ρβ]
=Fxθ[ ρα, ρβ]+ (Ex[ ρα, ρβ] − Fxθ[ ρα, ρβ])
=FxHF,θ[{fiα, ψiα}, {fiβ, ψiβ}]+ Ex,θ[ ρα, ρβ], (25) the sum of FxHF,θ[{fiα, ψiα}, {fiβ, ψiβ}] (i.e., the exact exchange in TAO-DFT) and Ex,θ[ ρα, ρβ] ≡ Ex[ ρα, ρβ] − Fxθ[ ρα, ρβ]
= Fxθ=0[ ρα, ρβ] − Fxθ[ ρα, ρβ] (i.e., the difference between the exchange free energy at zero temperature and that at the fic- titious temperature θ). Subsequently, a DFA can be made for Ex,θ[ ρα, ρβ] as follows:
Ex,θDFA[ ρα, ρβ] ≡ FxDFA,θ=0[ ρα, ρβ] − FxDFA,θ[ ρα, ρβ], (26) where FxDFA,θ[ ρα, ρβ] is the DFA for Fxθ[ ρα, ρβ]. Note that Ex,θ=0DFA [ ρα, ρβ]= 0 (i.e., an exact property of Ex,θ[ ρα, ρβ]) can be readily achieved by Eq.(26). Besides, from the spin-scaling relation of Fxθ[ ρα, ρβ],86 Ex,θDFA[ ρα, ρβ] can be conveniently expressed in terms of FxDFA,θ[ ρ] (in its spin-unpolarized form),
Ex,θDFA[ ρα, ρβ]=1
2(FxDFA,θ=0[2 ρα]+ FxDFA,θ=0[2 ρβ])
−1
2(FxDFA,θ[2 ρα]+ FxDFA,θ[2 ρβ]). (27)
From Eqs.(25)and(26), the exact exchange in KS-DFT is approximately given by
Ex[ ρα, ρβ]= ExHF[{φiα}, {φiβ}]
≈FxHF,θ[{fiα, ψiα}, {fiβ, ψiβ}]+ Ex,θDFA[ ρα, ρβ], (28) the sum of the exact exchange in TAO-DFT and Ex,θDFA[ ρα, ρβ].
Note that the approximation becomes exact, when the exact Ex,θDFA[ ρα, ρβ] is employed.
While the exact exchange in TAO-DFT is free of the SIE, the scheme is not expected to perform satisfactorily for most systems, due to the lack of correlation energy Ec[ ρα, ρβ].
Besides, it is well known that the exact exchange is incom- patible with the DFA correlation in KS-DFT, implying that TAO-DFT with the exact exchange and DFA correlation would not perform well for single-reference systems (i.e., in the absence of strong static correlation effects). There- fore, similar to the hybrid schemes in KS-DFT, it may be useful to incorporate the exact exchange with the DFA XC functional in TAO-DFT. In Subsections III BandIII C, the global and range-separated hybrid schemes in TAO-DFT are proposed.
B. Global hybrid scheme
In KS-DFT, a global hybrid (GH) functional23–30 is generally expressed as
EKS-GHxc = axExHF[{φiα}, {φiβ}]+ (1 − ax)ExDFA[ ρα, ρβ]
+ EcDFA[ ρα, ρβ], (29)
where ExHF is the HF exchange energy (given by Eq. (17)), ExDFAis the DFA exchange energy, and EcDFAis the DFA cor- relation energy. The fraction of HF exchange ax, typically ranging from 0.2 to 0.25 for thermochemistry and from 0.4 to 0.6 for kinetics, can be determined by empirical fitting or physical arguments.
After substituting Eq.(28)into Eq.(29), the correspond- ing global hybrid functional in TAO-DFT can be defined as
ExcTAO-GH= ax
FxHF,θ[{fiα, ψiα}, {fiβ, ψiβ}]+ Ex,θDFA[ ρα, ρβ]
+ (1 − ax)ExDFA[ ρα, ρβ]+ EcDFA[ ρα, ρβ], (30) and the resulting ground-state energy is evaluated by
ETAO-GH= Aθs[{fiα, ψiα}, {fiβ, ψiβ}]+
ρ(r)vext(r)dr+ EH[ ρ]
+ ExcTAO-GH+ EθDFA[ ρα, ρβ]. (31) While an evaluation of the functional derivative of FxHF,θ[{fiα, ψiα}, {fiβ, ψiβ}] (i.e., an explicit functional of the TAO orbitals and their occupation numbers) with respect to the density ρσ (see Eq.(7)) can be achieved with the finite- temperature exact-exchange and related schemes,88the result- ing scheme can be computationally demanding. To reduce the computational complexity, in this work, the electronic energy for a global hybrid functional in TAO-DFT is minimized with respect to the 1-RDM γσTAO(as is usual in the finite-temperature
HF (FT-HF) and related schemes72,76–83). The resulting self- consistent equations for the σ-spin electrons can be expressed as
−1
2∇2+ vs,σloc(r)
ψiσ(r) − ax
∞
X
j=1
fjσ
×
ψjσ∗ (r0)ψiσ(r0)
|r − r0| ψjσ(r)dr0= iσψiσ(r),
(32) where
vs,σloc(r)= vext(r)+
ρ(r0)
|r − r0|dr0+δEDFAθ [ ρα, ρβ] δ ρσ(r) + (1 − ax)δEDFAx [ ρα, ρβ]
δ ρσ(r) +δEcDFA[ ρα, ρβ] δ ρσ(r) + ax
δEx,θDFA[ ρα, ρβ]
δ ρσ(r) (33)
is the local part of the effective potential. The two sets (one for each spin function) of self-consistent equations, Equations(8)–(10),(32), and(33), for ρα(r) and ρβ(r), respec- tively, are coupled with the ground-state density (given by Eq.(11)).
Note that ExcTAO-GH reduces to ExcDFA (i.e., the DFA XC functional) for ax = 0, and reduces to FxHF,θ + Ex,θDFA + EcDFA (i.e., the exact exchange in TAO-DFT, the DFA for Ex,θ, and the DFA correlation functional) for ax = 1. At θ
= 0, TAO-DFT with ExcTAO-GH is the same as KS-DFT with ExcKS-GH.
On the other hand, if the constraints of ax= 1 and EcDFA
= EθDFA = Ex,θDFA = 0 are imposed on the global hybrid scheme in TAO-DFT, the resulting scheme resembles the FT-HF scheme. Therefore, the computational cost of the global hybrid scheme in TAO-DFT is similar to that of the global hybrid scheme in KS-DFT or the FT-HF scheme.
C. Range-separated hybrid scheme
In KS-DFT, a range-separated hybrid (RSH) func- tional31–33is generally given by
ExcKS-RSH= EHFx (I)[{φiα}, {φiβ}]+ EDFAx (¯I)[ ρα, ρβ] + EDFAc [ ρα, ρβ], (34) where EHFx (I) is the HF exchange energy of an interelectronic repulsion operator I(r12),
ExHF(I)[{φiα}, {φiβ}]= −1 2
α,β
X
σ Nσ
X
i,j=1
I(r12)φ∗iσ(r1)φ∗jσ(r2)
×φjσ(r1)φiσ(r2)dr1dr2, (35) and ExDFA(¯I) is the DFA exchange energy of the complementary operator ¯I(r12) ≡ 1/r12−I(r12). Similar to the previous trick, we replace the Coulomb operator 1/r12 in Eq. (28) by the operator I(r12), yielding the following expression:
EHFx (I)[{φiα}, {φiβ}] ≈ FxHF,θ(I)[{fiα, ψiα}, {fiβ, ψiβ}] + Ex,θDFA(I)[ ρα, ρβ], (36) where
FxHF,θ(I)[{fiα, ψiα}, {fiβ, ψiβ}]= −1 2
α,β
X
σ
∞
X
i,j=1
fiσfjσ
×
I(r12)ψiσ∗ (r1)ψjσ∗ (r2)
× ψjσ(r1)ψiσ(r2)dr1dr2 (37) is the HF exchange free energy of the operator I(r12) at the fictitious temperature θ, and
EDFAx,θ (I)[ ρα, ρβ] ≡ FxDFA,θ=0(I)[ ρα, ρβ] − FxDFA,θ(I)[ ρα, ρβ] (38) is the difference between the DFA exchange free energy of the operator I(r12) at zero temperature and that at the fictitious tem- perature θ. Note that the approximation (see Eq.(36)) becomes exact, when the exact Ex,θDFA(I)[ ρα, ρβ] is employed.
After substituting Eq.(36)into Eq.(34), the correspond- ing range-separated hybrid functional in TAO-DFT can be defined as
ExcTAO-RSH=
FxHF,θ(I)[{fiα, ψiα}, {fiβ, ψiβ}]+ Ex,θDFA(I)[ ρα, ρβ]
+ ExDFA(¯I)[ ρα, ρβ]+ EcDFA[ ρα, ρβ]. (39) For I(r12)= ax/r12, ETAO-RSHxc reduces to ExcTAO-GH. However, for a general operator I(r12) (e.g., the erf,89 erfgau,90 or terf91 operator), while FxHF,θ(I) is defined, and EDFAx (¯I) and EcDFA are available from those of the range-separated hybrid scheme in KS-DFT, FxDFA,θ(I) (and hence, Ex,θDFA(I)) is mostly unavailable and needs to be developed for practical appli- cations. Therefore, in this work, while the range-separated hybrid scheme in TAO-DFT is proposed, our numerical results are only available for the global hybrid scheme in TAO-DFT.
IV. GLOBAL HYBRID FUNCTIONALS IN TAO-DFT A. Definition of the optimal θ values
As previously mentioned, the fictitious temperature θ in TAO-DFT should be chosen so that the distribution of TOONs is close to that of NOONs.20,21In this situation, the strong static correlation effects can be properly described by the entropy contribution. For single-reference systems, as the exact NOONs are close to either 0 or 1, the optimal θ should be sufficiently small. However, for multi-reference systems, as the distribution of NOONs can be diverse (due to the varying strength of static correlation), the optimal θ can span a wide range of values. Therefore, for a global hybrid functional in TAO-DFT, it is impossible to adopt a θ that is optimal for both single- and multi-reference systems. Nevertheless, it remains useful to define an optimal θ value for a global hybrid func- tional in TAO-DFT to provide an explicit description of orbital occupations.
FIG. 1. Mean absolute errors of the reaction energies of the 30 chemical reactions in the NHTBH38/04 and HTBH38/04 sets,92calculated using TAO- B3LYP, TAO-B3LYP-D3, TAO-PBE0, and TAO-BHHLYP (with various θ).
The θ= 0 cases correspond to KS-B3LYP, KS-B3LYP-D3, KS-PBE0, and KS-BHHLYP, respectively.
To be consistent with the previous definition of the opti- mal θ value for TAO-DFAs, in this work, the same physi- cal arguments and numerical investigations are adopted to define the optimal θ value for a global hybrid functional in TAO-DFT. Specifically, the performance of various global hybrid functionals in TAO-DFT (with θ = 0, 5, 10, 15, 20, 25, 30, 35, 40, 45, and 50 mhartree) is examined for the following single-reference systems:
FIG. 2. Mean absolute errors of the 166 bond lengths in the EXTS set,93 calculated using TAO-B3LYP, TAO-B3LYP-D3, TAO-PBE0, and TAO- BHHLYP (with various θ). The θ = 0 cases correspond to KS-B3LYP, KS-B3LYP-D3, KS-PBE0, and KS-BHHLYP, respectively.
TABLE I. Optimal fictitious temperature θ (in mhartree), given by Eq.(43), for TAO-B3LYP, TAO-B3LYP-D3, TAO-PBE0, and TAO-BHHLYP, where axis the fraction of exact exchange.
TAO-B3LYP TAO-B3LYP-D3 TAO-PBE0 TAO-BHHLYP
ax 1/5 1/5 1/4 1/2
θ 17.4 17.4 20 33
TABLE II. Statistical errors (in kcal/mol) of the reaction energies of the 30 chemical reactions in the NHTBH38/04 and HTBH38/04 sets,92calculated using TAO-B3LYP, TAO-B3LYP-D3, TAO-PBE0, and TAO-BHHLYP (with the optimal θ values given in TableI). The θ= 0 cases correspond to KS-B3LYP, KS-B3LYP-D3, KS-PBE0, and KS-BHHLYP, respectively.
KS-DFT TAO-DFT
B3LYP B3LYP-D3 PBE0 BHHLYP B3LYP B3LYP-D3 PBE0 BHHLYP
MSE −0.23 −0.27 −0.03 −1.25 −0.66 −0.70 −0.41 −1.76
MAE 2.01 1.95 2.41 3.63 2.33 2.36 2.63 3.95
rms 2.66 2.61 3.35 4.72 3.05 3.07 3.69 5.00
Max(−) −7.38 −7.41 −7.11 −14.00 −8.44 −8.46 −8.40 −14.21
Max(+) 4.46 4.13 10.20 7.63 4.34 4.01 10.52 6.55
TABLE III. Statistical errors (in Å) of the 166 bond lengths in the EXTS set,93calculated using TAO-B3LYP, TAO-B3LYP-D3, TAO-PBE0, and TAO-BHHLYP (with the optimal θ values given in TableI). The θ= 0 cases correspond to KS-B3LYP, KS-B3LYP-D3, KS-PBE0, and KS-BHHLYP, respectively.
KS-DFT TAO-DFT
B3LYP B3LYP-D3 PBE0 BHHLYP B3LYP B3LYP-D3 PBE0 BHHLYP
MSE 0.003 0.003 −0.002 −0.012 0.003 0.003 −0.002 −0.014
MAE 0.008 0.008 0.008 0.013 0.008 0.008 0.008 0.015
rms 0.013 0.013 0.012 0.017 0.013 0.014 0.013 0.019
Max(−) −0.078 −0.078 −0.082 −0.090 −0.080 −0.080 −0.085 −0.095
Max(+) 0.065 0.065 0.051 0.025 0.063 0.063 0.049 0.035
TABLE IV. Statistical errors (in kcal/mol) of the ωB97 training set,32calculated using TAO-B3LYP, TAO-B3LYP-D3, TAO-PBE0, and TAO-BHHLYP (with the optimal θ values given in TableI). The θ= 0 cases correspond to KS-B3LYP, KS-B3LYP-D3, KS-PBE0, and KS-BHHLYP, respectively.
KS-DFT TAO-DFT
System Error B3LYP B3LYP-D3 PBE0 BHHLYP B3LYP B3LYP-D3 PBE0 BHHLYP
G3/99 MSE −4.30 −1.99 3.94 −29.55 0.90 3.21 11.48 −11.32
(223) MAE 5.46 3.64 6.28 29.68 5.25 6.80 13.34 12.59
rms 7.34 5.23 8.65 34.13 6.97 8.31 17.16 16.62
IP MSE 2.18 2.17 −0.13 −1.72 0.25 0.24 −2.34 −5.66
(40) MAE 3.68 3.69 3.33 4.44 4.25 4.26 4.37 7.04
rms 4.81 4.81 3.98 5.47 5.30 5.31 5.27 8.19
EA MSE 1.71 1.71 −1.07 −4.79 −1.02 −1.02 −4.30 −9.98
(25) MAE 2.38 2.39 3.10 5.97 3.48 3.49 4.63 10.21
rms 3.27 3.29 3.53 6.84 4.50 4.52 5.42 11.38
PA MSE −0.77 −0.66 0.18 −0.12 0.14 0.26 1.25 1.75
(8) MAE 1.16 1.07 1.14 1.55 0.91 1.01 1.42 2.02
rms 1.36 1.33 1.61 1.78 1.21 1.27 2.03 2.63
NHTBH MSE −4.57 −5.09 −3.13 0.52 −4.88 −5.39 −3.53 −0.51
(38) MAE 4.69 5.19 3.63 2.21 5.08 5.56 4.18 2.75
rms 5.71 6.14 4.63 2.93 6.02 6.49 5.10 3.29
HTBH MSE −4.48 −5.12 −4.60 0.58 −5.20 −5.84 −5.55 −1.43
(38) MAE 4.56 5.14 4.60 2.48 5.20 5.84 5.55 2.40
rms 5.10 5.62 4.88 3.11 5.79 6.34 5.80 3.15
S22 MSE 3.95 −0.02 2.50 2.98 2.74 −1.22 1.10 0.18
(22) MAE 3.95 0.43 2.52 3.01 2.76 1.22 1.49 1.42
rms 5.17 0.59 3.62 4.22 3.98 1.37 2.40 1.98
Total MSE −2.77 −1.80 1.55 −16.93 −0.34 0.63 5.20 −7.76
(394) MAE 4.75 3.63 5.05 18.28 4.79 5.69 9.34 9.11
rms 6.38 5.03 7.06 25.85 6.27 7.17 13.33 13.18
