Range-separated multideterminant density-functional theory with a short-range correlation functional of the on-top pair density

TL;DR

This paper presents a new local approximation for short-range correlation energy in range-separated density-functional theory, improving the accuracy of dissociation energy calculations for small molecules by interpolating between known limits.

Contribution

It introduces a novel local functional of density, gradient, and on-top pair density that interpolates between established correlation functionals at different range-separation parameters.

Findings

Accurate dissociation energy curves for H2, Li2, Be2 molecules.

Effective in capturing static correlation effects.

Promising approach for improved range-separated DFT calculations.

Abstract

We introduce an approximation to the short-range correlation energy functional with multide-terminantal reference involved in a variant of range-separated density-functional theory. This approximation is a local functional of the density, the density gradient, and the on-top pair density, which locally interpolates between the standard Perdew-Burke-Ernzerhof correlation functional at vanishing range-separation parameter and the known exact asymptotic expansion at large range-separation parameter. When combined with (selected) configuration-interaction calculations for the long-range wave function, this approximation gives accurate dissociation energy curves of the H2, Li2, and Be2 molecules, and thus appears as a promising way to accurately account for static correlation in range-separated density-functional theory.

Equations

$$E_0=\underset{n}{\min }\left\{\mathcal{F}[n]+\int v_{\mathrm{ne}}(\text{r})n(\text{r})\mathrm{d}\text{r}\right\},$$

$$\mathcal{F}[n]=\underset{\Psi \to n}{\min }⟨\Psi |\hat{T}+{\hat{W}}_{\mathrm{ee}}|\Psi ⟩,$$

$$\mathcal{F}[n]={\mathcal{F}}^{\mathrm{lr},\mu }[n]+{\bar{E}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu }[n],$$

$${\mathcal{F}}^{\mathrm{lr},\mu }[n]=\underset{\Psi \to n}{\min }⟨\Psi |\hat{T}+{\hat{W}}_{\mathrm{ee}}^{\mathrm{lr},\mu }|\Psi ⟩,$$

$${\hat{W}}_{\mathrm{ee}}^{\mathrm{lr},\mu }=\frac{1}{2}\iint w_{\mathrm{ee}}^{\mathrm{lr},\mu }(r_{12}){\hat{n}}_2({\text{r}}_1,{\text{r}}_2)\mathrm{d}{\text{r}}_1\mathrm{d}{\text{r}}_2,$$

$$\displaystyle E_{0}=\min_{\Psi}\Big{\{}\langle\Psi|\hat{T}+\hat{W}_{\mathrm{{ee}}}^{\mathrm{lr},\mu}+\hat{V}_{\mathrm{ne}}|\Psi\rangle+\bar{E}^{\mathrm{sr},\mu}_{\mathrm{Hxc}}[n_{\Psi}]\Big{\}},$$

$${\hat{H}}^{\mu }[n_{{\Psi }^{\mu }}]|{\Psi }^{\mu }⟩={\mathcal{E}}^{\mu }|{\Psi }^{\mu }⟩,$$

$${\hat{H}}^{\mu }[n_{{\Psi }^{\mu }}]=\hat{T}+{\hat{W}}_{\mathrm{ee}}^{\mathrm{lr},\mu }+{\hat{V}}_{\mathrm{ne}}+{\hat{\bar{V}}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu }[n_{{\Psi }^{\mu }}],$$

$$E_0=⟨{\Psi }^{\mu }⟩\hat{T}+{\hat{W}}_{\mathrm{ee}}^{\mathrm{lr},\mu }+{\hat{V}}_{\mathrm{ne}}{\Psi }^{\mu }+{\bar{E}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu }[n_{{\Psi }^{\mu }}].$$

$${\bar{E}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu }[n]=E_{\mathrm{H}}^{\mathrm{sr},\mu }[n]+E_{\mathrm{x}}^{\mathrm{sr},\mu }[n]+{\bar{E}}_{\mathrm{c}}^{\mathrm{sr},\mu }[n],$$

$$E_{\mathrm{H}}^{\mathrm{sr},\mu }[n]=\frac{1}{2}\iint w_{\mathrm{ee}}^{\mathrm{sr},\mu }(r_{12})n({\text{r}}_1)n({\text{r}}_2)\mathrm{d}{\text{r}}_1\mathrm{d}{\text{r}}_2,$$

$$E_{\mathrm{x}}^{\mathrm{sr},\mu }[n]=⟨{\Phi }^{\text{textsc}ks}[n]⟩{\hat{W}}_{\mathrm{ee}}^{\mathrm{sr},\mu }{\Phi }^{\text{textsc}ks}[n]-E_{\mathrm{H}}^{\mathrm{sr},\mu }[n],$$

$${\bar{E}}_{\mathrm{c}}^{\mathrm{sr},\mu }[n]={\bar{E}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu }[n]-⟨{\Phi }^{\text{textsc}ks}[n]⟩{\hat{W}}_{\mathrm{ee}}^{\mathrm{sr},\mu }{\Phi }^{\text{textsc}ks}[n],$$

$${\hat{W}}_{\mathrm{ee}}^{\mathrm{sr},\mu }=\frac{1}{2}\iint w_{\mathrm{ee}}^{\mathrm{sr},\mu }(r_{12}){\hat{n}}_2({\text{r}}_1,{\text{r}}_2)\mathrm{d}{\text{r}}_1\mathrm{d}{\text{r}}_2.$$

$${\bar{E}}_{x/c}^{\mathrm{sr},\mu ,\text{textsc}pbe}[n]=\int {\bar{e}}_{x/c}^{\mathrm{sr},\mu ,\text{textsc}\mathrm{pbe}}(n(\text{r}),\nabla n(\text{r}))\mathrm{d}\text{r}.$$

$${\bar{E}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu }[n]={\bar{E}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu }[n]-⟨{\Psi }^{\mu }[n]⟩{\hat{W}}_{\mathrm{ee}}^{\mathrm{sr},\mu }{\Psi }^{\mu }[n].$$

$$E_0=⟨{\Psi }^{\mu }⟩\hat{H}{\Psi }^{\mu }+{\bar{E}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu }[n_{{\Psi }^{\mu }}],$$

$${\bar{E}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu }[n]={\bar{E}}_{\mathrm{c}}^{\mathrm{sr},\mu }[n]+{\Delta }^{\mathrm{lr}\text{textsc}-\mathrm{sr},\mu }[n],$$

$${\Delta }^{\mathrm{lr}\text{textsc}-\mathrm{sr},\mu }[n]=$$

$$-⟨{\Psi }^{\mu }[n]⟩{\hat{W}}_{\mathrm{ee}}^{\mathrm{sr},\mu }{\Psi }^{\mu }[n].$$

$${\bar{E}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu =0}[n]={\bar{E}}_{\mathrm{c}}^{\mathrm{sr},\mu =0}[n]=E_{\mathrm{c}}^{\text{textsc}ks}[n].$$

$${\bar{E}}_{\mathrm{c}}^{\mathrm{sr},\mu }[n]$$

$${\Delta }^{\mathrm{lr}\text{textsc}-\mathrm{sr},\mu }[n]=-\frac{\pi }{2{\mu }^2}\int n_{2,\mathrm{c}}(\text{r},\text{r})\mathrm{d}\text{r}$$

$$-\frac{2\sqrt{\pi }(2\sqrt{2}-1)}{3{\mu }^3}\int n_2(\text{r},\text{r})\mathrm{d}\text{r}+O\left(\frac{1}{{\mu }^4}\right).$$

$${\bar{E}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu }[n]=\frac{2\sqrt{\pi }(1-\sqrt{2})}{3{\mu }^3}\int n_2(\text{r},\text{r})\mathrm{d}\text{r}+O\left(\frac{1}{{\mu }^4}\right).$$

$${\bar{E}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu ,\text{textsc}\mathrm{lda}}[n]=\int {\bar{e}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu ,\text{textsc}\mathrm{lda}}(n(\text{r}))\mathrm{d}\text{r},$$

$${\bar{E}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu ,\text{textsc}pbe\text{ontop}}[n]=$$

$$\int {\bar{e}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu ,\text{textsc}pbe\text{ontop}}(n(\text{r}),\nabla n(\text{r}),n_2(\text{r},\text{r}))$$

$${\bar{e}}_{\mathrm{c},\mathrm{md}}^{\mathrm{sr},\mu ,\text{textsc}pbe\text{ontop}}(n,\nabla n,n_2)=\frac{e_{\mathrm{c}}^{\text{textsc}pbe}(n,\nabla n)}{1+\beta (n,\nabla n,n_2){\mu }^3},$$

$$\beta (n,\nabla n,n_2)=\frac{3e_{\mathrm{c}}^{\text{textsc}pbe}(n,\nabla n)}{2\sqrt{\pi }(1-\sqrt{2})n_2}.$$

$$n_2^{\mu }(\text{r},\text{r})=n_2(\text{r},\text{r})\left(1+\frac{2}{\sqrt{\pi }\mu }\right)+O\left(\frac{1}{{\mu }^2}\right),$$

$$n_2(\text{r},\text{r})\approx n_2^{\mu }(\text{r},\text{r}){\left(1+\frac{2}{\sqrt{\pi }\mu }\right)}^{-1}.$$

Full text

I Introduction

Range-separated density-functional theory (RS-DFT) (see, e.g., Refs. Savin, 1996; Toulouse et al., 2004a) is an alternative to Kohn-Sham density-functional theory (KS-DFT) Kohn and Sham (1965) for electronic-structure calculations of atoms, molecules, and solids. It consists in rigorously combining a wave-function-type calculation for the long-range part of the Coulomb electron-electron interaction with a density functional for the complementary short-range part of the interaction. This permits to describe long-range electron correlation accurately and short-range electron correlation compactly with a fast basis-set convergence Franck et al. (2015). In particular, it has been shown that explicit static correlation effects can be effectively taken into account in the long-range part of the calculation by using methods such as configuration interaction (CI) Leininger et al. (1997); Pollet et al. (2002); Casanova (2018), multiconfiguration self-consistent field Fromager et al. (2007, 2009); Hedegård et al. (2018), multireference perturbation theory Fromager et al. (2010), density-matrix functional theory Pernal (2010); Rohr et al. (2010); Rohr and Pernal (2011), density-matrix renormalization group Hedegård et al. (2015), or pair coupled-cluster doubles Garza et al. (2015a).

A major limitation to the accuracy of RS-DFT are the semilocal density-functional approximations used for the short-range exchange-correlation energy Toulouse et al. (2004b, a, 2005a); Goll et al. (2005); Paziani et al. (2006); Goll et al. (2006), which still suffer from self-interaction (or fractional-charge) errors and static-correlation (or fractional-spin) errors Mussard and Toulouse (2017). An attractive remedy to this problem is to calculate exactly a large portion of the short-range exchange-correlation energy using the multideterminant wave function naturally available in RS-DFT, leaving only a residual short-range correlation energy functional to be approximated Toulouse (2005); Toulouse et al. (2005b). In particular, this permits to drastically reduce self-interaction errors since the exchange energy is now calculated with a wave function and not with an approximate exchange density functional. This strategy has been pursued in various RS-DFT approaches Stoyanova et al. (2013); Cornaton et al. (2013); Cornaton and Fromager (2014); Rebolini et al. (2015); Senjean et al. (2016); Alam et al. (2016, 2017); Rebolini et al. (2018). However, only a local-density approximation (LDA) for this short-range correlation energy functional with multideterminantal reference was available so far Toulouse et al. (2005b); Paziani et al. (2006), which tends to substantially overcorrelate.

In the present work, we develop an approximation for this short-range multideterminant correlation energy functional which uses the density, the density gradient, and the on-top pair density. The development of correlation functionals depending on the on-top pair density extracted from a multideterminant wave function has started long ago Colle and Salvetti (1975, 1979) and is still an active area of research (see, e.g., Ref. Hollett and Pegoretti, 2018). An important motivation for using the on-top pair density of a multideterminant wave function is that it clearly contains information about bond dissociation (see, e.g., Ref. Carlson et al., 2017), without having to artificially break spin symmetry. In practice, most of the works in this domain introduce the on-top pair density via effective spin densities which are fed into standard spin-dependent exchange and/or correlation density functionals Moscardó and San-Fabián (1991); Becke et al. (1995); Miehlich et al. (1997); Takeda et al. (2002, 2004); Gräfenstein and Cremer (2005); Tsuchimochi et al. (2010); G. Li Manni, R. K. Carlson, S. Luo, D. Ma, J. Olsen, D. G. Truhlar, and L. Gagliardi (2014); Garza et al. (2015b, a); Carlson et al. (2015); L. Gagliardi, D. G. Truhlar, G. Li Manni, R. K. Carlson, C. E. Hoyer, and J. Lucas Bao (2017). This is justified by the alternative interpretation of spin-density-functional theory Perdew et al. (1995, 1997), in which the spin densities are viewed as mere intermediate quantities for reproducing the total density and the on-top pair density. Here, instead of using effective spin densities, we introduce the dependence on the on-top pair density by exploiting the known exact asymptotic behavior of the short-range multideterminant correlation energy functional in the limit where the electron-electron interaction is very short ranged Toulouse et al. (2005b); Paziani et al. (2006).

The paper is organized as follows. In Sec. II, we briefly review RS-DFT, including the approach involving the short-range correlation functional with multideterminant reference, and develop a new approximation for this functional. After giving computational details in Sec. III, in particular on the selected CI method that we use for the long-range wave function, we discuss in Sec. IV the results concerning the convergence with respect to the number of determinants on the Ne atom and the Be2 molecule, the dependence on the range-separation parameter on the He and C atoms and the H2 molecule near dissociation, and the dissociation energy curves of the H2, Li2, and Be2 molecules. Finally, Sec. V summarizes our conclusions.

II Theory

II.1 Range-separated density-functional theory

The exact ground-state energy of a $N$-electron system with nuclei-electron potential $v_{\mathrm{ne}}(\textbf{r})$ can be expressed by the following minimization over $N$-representable densities $n$ Levy (1979); Lieb (1983)

[TABLE]

with the standard constrained-search universal density functional

[TABLE]

where $\hat{T}$ and $\hat{W}_{\mathrm{ee}}$ are the kinetic-energy and Coulomb electron-electron interaction operators, respectively. The minimizing multideterminant wave function in Eq. (2) will be denoted by $\Psi[n]$.

In RS-DFT, the universal density functional is decomposed as Savin (1996); Toulouse et al. (2004a)

[TABLE]

where $\mathcal{F}^{\mathrm{lr},\mu}[n]$ is a long-range (lr) universal density functional

[TABLE]

and $\bar{E}_{\mathrm{Hxc}}^{\,\mathrm{sr,}\mu}[n]$ is the complementary short-range (sr) Hartree-exchange-correlation (Hxc) density functional. In Eq. (4), $\hat{W}_{\mathrm{ee}}^{\mathrm{lr}}$ is the long-range electron-electron interaction defined as

[TABLE]

with the error-function potential $w_{\mathrm{ee}}^{\mathrm{lr},\mu}(r_{12})=\mathrm{erf}(\mu\,r_{12})/r_{12}$ (expressed with the interelectronic distance $r_{12}=||\textbf{r}_{1}-\textbf{r}_{2}||$) and the pair-density operator $\hat{n}_{2}(\textbf{r}_{1},\textbf{r}_{2})=\hat{n}(\textbf{r}_{1})\hat{n}(\textbf{r}_{2})-\delta(\textbf{r}_{1}-\textbf{r}_{2})\hat{n}(\textbf{r}_{1})$ where $\hat{n}(\textbf{r})$ is the density operator. The range-separation parameter $\mu$ corresponds to an inverse distance controlling the range of the separation. For a given density, we will denote by $\Psi^{\mu}[n]$ the minimizing multideterminant wave function in Eq. (4). Inserting the decomposition of Eq. (3) into Eq. (1), and recomposing the two-step minimization into a single one, leads to the following expression for the exact ground-state electronic energy

[TABLE]

where the minimization is done over normalized $N$-electron multideterminant wave functions, $\hat{V}_{\mathrm{ne}}=\int v_{\mathrm{ne}}(\textbf{r})\hat{n}(\textbf{r})\mathrm{d}\textbf{r}$, and $n_{\Psi}$ refers to the density of $\Psi$, i.e. $n_{\Psi}(\textbf{r})=\langle\Psi|\hat{n}(\textbf{r})|\Psi\rangle$. The minimizing multideterminant wave function $\Psi^{\mu}$ in Eq. (II.1) can be determined by the self-consistent eigenvalue equation

[TABLE]

with the long-range interacting Hamiltonian

[TABLE]

where $\hat{\bar{V}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu}[n]=\int\delta\bar{E}^{\mathrm{sr},\mu}_{\mathrm{Hxc}}[n]/\delta n(\textbf{r})\,\hat{n}(\textbf{r})\mathrm{d}\textbf{r}$ is the complementary short-range Hartree-exchange-correlation potential operator. Note that $\Psi^{\mu}$ is not the exact physical ground-state wave function but only an effective wave function. However, its density $n_{\Psi^{\mu}}$ is the exact physical ground-state density. Once $\Psi^{\mu}$ has been calculated, the exact electronic ground-state energy is obtained by

[TABLE]

Note that, for $\mu=0$, the long-range interaction vanishes, $w_{\mathrm{ee}}^{\mathrm{lr},\mu=0}(r_{12})=0$, and thus RS-DFT reduces to standard KS-DFT. For $\mu\to\infty$, the long-range interaction becomes the standard Coulomb interaction, $w_{\mathrm{ee}}^{\mathrm{lr},\mu\to\infty}(r_{12})=1/r_{12}$, and thus RS-DFT reduces to standard wave-function theory (WFT).

In principle, Eq. (7) should be solved at the full-configuration-interaction (FCI) level in a complete one-electron basis set. In practice, however, for typical values of the range-separation parameter used (around $\mu=0.5$ bohr*-1*) Gerber and Ángyán (2005a); Fromager et al. (2007), $\hat{H}^{\mu}[n_{\Psi^{\mu}}]$ contains only a non-diverging soft long-range electron-electron interaction, implying that the wave function $\Psi^{\mu}$ does not have an electron-electron cusp Gori-Giorgi and Savin (2006) and has a fast convergence with respect to the number of determinants or with respect to the size of the one-electron basis Franck et al. (2015). One can then accurately solve Eq. (7) using efficient truncated or selected CI approaches, such as the configuration interaction perturbatively selected iteratively (CIPSI) method Bender and Davidson (1969); Huron et al. (1973); Buenker and Peyerimholf (1974); Buenker et al. (1981); Evangelisti et al. (1983); Rubio et al. (1986); Cimiraglia and Persico (1987); Harrison (1991); Angeli and Persico (1997); Angeli et al. (2000); Angeli and Cimiraglia (2001); Caffarel et al. (2014); Giner and Angeli (2015) (see Sec. III), with relatively small basis sets. The resulting compact wave function $\Psi^{\mu}$ will accurately include the long-range electron correlation effects.

As regards the short-range density functional, it is usually decomposed into three contributions

[TABLE]

where $E^{\mathrm{sr},\mu}_{\mathrm{H}}[n]$ is the short-range Hartree energy functional

[TABLE]

with the short-range electron-electron interaction $w_{\mathrm{ee}}^{\mathrm{sr},\mu}(r_{12})=1/r_{12}-w_{\mathrm{ee}}^{\mathrm{lr},\mu}(r_{12})$, and $E^{\mathrm{sr},\mu}_{\mathrm{x}}[n]$ and $\bar{E}_{\mathrm{c}}^{\mathrm{sr,}\mu}[n]$ are the short-range exchange and correlation energy functionals

[TABLE]

[TABLE]

defined with the Kohn-Sham (KS) single-determinant wave function $\Phi^{\textsc{ks}}[n]=\Psi^{\mu=0}[n]$ and the short-range electron-electron interaction operator

[TABLE]

Whereas $E^{\mathrm{sr},\mu}_{\mathrm{H}}[n]$ is calculated exactly, approximations need to be used for $E^{\mathrm{sr},\mu}_{\mathrm{x}}[n]$ and $\bar{E}_{\mathrm{c}}^{\mathrm{sr,}\mu}[n]$. In this work, we use the short-range version of the Perdew-Burke-Ernzerhof (PBE) Perdew et al. (1996) exchange and correlation functionals of Ref. Goll et al., 2006 (see also Refs. Toulouse et al., 2005a; Goll et al., 2005) which takes the form

[TABLE]

It has been shown that such semi-local density-functional approximations become more accurate as the range of the electron-electron interaction is reduced Toulouse et al. (2004a). Nevertheless, for the values of the range-separation parameter commonly used, the short-range PBE exchange and correlation density functionals still contain substantial self-interaction and static-correlation errors Mussard and Toulouse (2017).

II.2 Short-range correlation energy functional with multideterminant reference

The definition of the short-range correlation energy functional in Eq. (13) is based on the KS single-determinant wave function $\Phi^{\textsc{ks}}[n]$. In RS-DFT, it is in fact more natural to define another short-range correlation energy functional based on the multideterminant (md) wave function $\Psi^{\mu}[n]$ Toulouse (2005); Toulouse et al. (2005b)

[TABLE]

In lieu of the standard expression of the ground-state energy in the context of RS-DFT using only the long-range electron-electron interaction in the expectation value over the wave function $\Psi^{\mu}$ as described by Eq. (9), we can now easily include the full-range interaction in the expectation value by writing the exact ground-state electronic energy as

[TABLE]

where $\hat{H}=\hat{T}+\hat{W}_{\mathrm{{ee}}}+\hat{V}_{\mathrm{ne}}$ is the complete electronic Hamiltonian. This allows one to extract as much information as possible from the wave function $\Psi^{\mu}$ by calculating exactly the short-range Hartree and “exchange” energies related to it, i.e. the term $\braket{\Psi^{\mu}}{\hat{W}^{\mathrm{sr},\mu}_{\mathrm{ee}}}{\Psi^{\mu}}$. Since the wave function $\Psi^{\mu}$ is obtained without considering the short-range component of the electron-electron interaction, some short-range correlation is still missing in $\braket{\Psi^{\mu}}{\hat{H}}{\Psi^{\mu}}$ and must be recovered by the complementary multideterminant short-range correlation energy functional ${{\bar{E}_{\mathrm{{c,md}}}^{\mathrm{sr,}{\mu}}[n]}}$. Obviously, in practice, this functional must be approximated, but calculating the energy via Eq. (17) instead of Eq. (9) reduces the demand put on density-functional approximations. In particular, since in Eq. (17) the whole exchange energy is calculated with a wave function and not with an approximate exchange density functional, we expect to eliminate most of the self-interaction error. We note that, contrary to the expression in Eq. (9), the energy expression in Eq. (17) is not variational with respect to $\Psi^{\mu}$. Even though it is possible to formulate a self-consistent version of Eq. (17) via a multideterminant extension of the optimized-effective-potential (OEP) approach Toulouse et al. (2005b); Stoyanova et al. (2013), we do not consider this possibility in this work.

In order to construct an approximation for ${{\bar{E}_{\mathrm{{c,md}}}^{\mathrm{sr,}{\mu}}[n]}}$, we now study two exact conditions on this functional. For this, it is convenient to express the functional $\bar{E}_{\mathrm{c,md}}^{\mathrm{sr},\mu}[n]$ in terms of the original functional $\bar{E}_{\mathrm{c}}^{\mathrm{sr},\mu}[n]$, using Eqs. (13) and (16),

[TABLE]

where $\Delta^{\mathrm{lr}\textsc{-}\mathrm{sr},\mu}[n]$ is a mixed long-range/short-range quantity

[TABLE]

We expect for most systems that $\Delta^{\mathrm{lr}\textsc{-}\mathrm{sr},\mu}\geq 0$, i.e. $|\bar{E}_{\mathrm{c,md}}^{\mathrm{sr},\mu}|\leq|\bar{E}^{\mathrm{sr},\mu}_{\mathrm{c}}|$.

The first condition is for $\mu=0$. In this case, since the RS-DFT wave function reduces to the KS wave function, $\Psi^{\mu=0}[n]=\Phi^{\textsc{ks}}[n]$, the short-range multideterminant correlation functional reduces to the usual KS correlation functional

[TABLE]

The second condition is for $\mu\to\infty$. In this limit, the asymptotic expansion of $\bar{E}^{\mathrm{sr},\mu}_{\mathrm{c}}[n]$ is known to be Toulouse et al. (2004a); Gori-Giorgi and Savin (2006)

[TABLE]

where $n_{2}(\textbf{r},\textbf{r})=\braket{\Psi[n]}{\hat{n}_{2}(\textbf{r},\textbf{r})}{\Psi[n]}$ is the Coulombic on-top pair density (i.e., the on-top pair density associated with the full-range wave function $\Psi[n]$) and $n_{2,\mathrm{c}}(\textbf{r},\textbf{r})=n_{2}(\textbf{r},\textbf{r})-n_{2,\textsc{ks}}(\textbf{r},\textbf{r})$ is its correlation contribution defined with respect to the KS on-top pair density $n_{2,\textsc{ks}}(\textbf{r},\textbf{r})=\braket{\Phi^{\textsc{ks}}[n]}{\hat{n}_{2}(\textbf{r},\textbf{r})}{\Phi^{\textsc{ks}}[n]}$. The asymptotic expansion of $\Delta^{\mathrm{lr}\textsc{-}\mathrm{sr},\mu}[n]$ for $\mu\to\infty$ can be obtained by generalizing the expansion given in the case of the homogeneous electron gas in Ref. Paziani et al., 2006, leading to

[TABLE]

The terms in $1/\mu^{2}$ in Eqs. (21) and (II.2) cancel each other, and we get the asymptotic expansion of $\bar{E}_{\mathrm{{c,md}}}^{\mathrm{{sr},\mu}}[n]$ for $\mu\to\infty$

[TABLE]

The short-range multideterminant correlation functional $\bar{E}_{\mathrm{{c,md}}}^{\mathrm{{sr},\mu}}[n]$ goes to zero as $1/\mu^{3}$ when $\mu\to\infty$, i.e. faster that the original short-range correlation functional $\bar{E}^{\mathrm{sr},\mu}_{\mathrm{c}}[n]$ of RS-DFT. This is not a surprise since $\bar{E}_{\mathrm{{c,md}}}^{\mathrm{{sr},\mu}}[n]$ accounts for a smaller part of the correlation energy than $\bar{E}^{\mathrm{sr},\mu}_{\mathrm{c}}[n]$. We thus see that, because of the local nature of the short-range interaction for a large value of $\mu$, the on-top pair density $n_{2}(\textbf{r},\textbf{r})$ appears naturally as a key ingredient in the short-range multideterminant correlation functional $\bar{E}_{\mathrm{{c,md}}}^{\mathrm{{sr},\mu}}[n]$.

II.3 Approximations for the short-range multideterminant correlation functional $\bar{E}_{\mathrm{{c,md}}}^{\mathrm{{sr},\mu}}[n]$

Until now, the only approximation available for the functional $\bar{E}_{\mathrm{{c,md}}}^{\mathrm{{sr},\mu}}[n]$ was the short-range LDA (srLDA) approximation Toulouse et al. (2005b); Paziani et al. (2006)

[TABLE]

where $\bar{e}_{\mathrm{{c,md}}}^{\mathrm{{sr},\mu,\textsc{lda}}}(n)$ is the energy density extracted from the homogeneous electron gas for which a parametrization is given in Ref. Paziani et al., 2006. Unfortunately, this srLDA approximation tends to give substantially too negative correlation energies for small values of $\mu$ (and in particular for the values commonly used, i.e. around $\mu=0.5$ bohr*-1*) Toulouse et al. (2005b); Stoyanova et al. (2013).

Here, we construct a new approximation for the functional $\bar{E}_{\mathrm{{c,md}}}^{\mathrm{{sr},\mu}}[n]$ based on the two exact conditions in Eqs. (20) and (23). We propose a local interpolation between the standard PBE correlation functional at $\mu=0$ (of course, any other generalized-gradient approximation to the KS correlation functional could be used) and the leading term of the asymptotic expansion of $\bar{E}_{\mathrm{{c,md}}}^{\mathrm{{sr},\mu}}[n]$ for $\mu\to\infty$. The resulting approximation, referred to as “srPBEontop”, is a local functional of the density, the density gradient, and the on-top pair density

[TABLE]

where the energy density is taken as

[TABLE]

which reduces to the standard PBE correlation energy density $e_{\mathrm{c}}^{\textsc{pbe}}(n,\nabla n)$ for $\mu=0$. In order to recover the correct large-$\mu$ behavior in Eq. (23), $\beta(n,\nabla n,n_{2})$ is taken as

[TABLE]

However, there is one difficulty with using the approximation in Eq. (25): the Coulombic on-top pair density $n_{2}(\textbf{r},\textbf{r})$ is not available in RS-DFT. Instead, what is available is the on-top pair density of the wave function $\Psi^{\mu}$ obtained with a long-range electron-electron interaction: $n_{2}^{\mu}(\textbf{r},\textbf{r})=\braket{\Psi^{\mu}}{\hat{n}_{2}(\textbf{r},\textbf{r})}{\Psi^{\mu}}$. Fortunately, the Coulombic on-top pair density $n_{2}(\textbf{r},\textbf{r})$ can be extrapolated from the long-range on-top pair density $n_{2}^{\mu}(\textbf{r},\textbf{r})$, as shown in Ref. Gori-Giorgi and Savin, 2006. The extrapolation is based on the asymptotic expansion of $n_{2}^{\mu}(\textbf{r},\textbf{r})$ for $\mu\to\infty$ Gori-Giorgi and Savin (2006)

[TABLE]

which, after inversion, gives the following estimation of the Coulombic on-top pair density

[TABLE]

Obviously, the Coulombic on-top pair density $n_{2}(\textbf{r},\textbf{r})$ is smaller than the long-range one $n_{2}^{\mu}(\textbf{r},\textbf{r})$. Since the latter is obtained with a reduced electron-electron repulsion, the probability of finding two electrons at the same point of space is larger. Note that, in the limit $\mu=0$, the extrapolation formula in Eq. (29) just unphysically gives $n_{2}(\textbf{r},\textbf{r})=0$ for all systems. However, this is not a problem since in the srPBEontop functional of Eq. (25) the on-top pair density $n_{2}(\textbf{r},\textbf{r})$ has an effect only for not too small values of $\mu$. In Ref. Gori-Giorgi and Savin, 2006, another extrapolation method based on the pair-distribution function of the homogeneous electron gas was also proposed. We do not consider this latter extrapolation in the present work since we have found that the simple one in Eq. (29) gives satisfying results.

For one-electron systems, the on-top pair density $n_{2}(\textbf{r},\textbf{r})$ or $n_{2}^{\mu}(\textbf{r},\textbf{r})$ vanishes, and consequently, from Eqs. (26) and (27), the srPBEontop correlation energy correctly vanishes as well. In other words, the srPBEontop correlation functional is self-interaction free for one-electron systems. In many-electron systems, we expect the same behavior in spatial regions of one-electron character.

III Computational details

We have implemented the RS-DFT approach, including the short-range multideterminant correlation functionals, in the software Quantum Package Scemama et al. (2015).

In practice, we first perform calculations with the self-consistent RS-DFT approach of Eq. (II.1) using the srPBE approximation of Ref. Goll et al., 2006 for the short-range exchange-correlation functional $\bar{E}_{\mathrm{{xc}}}^{\mathrm{sr,}{\mu}}[n]$. We calculate the multideterminant wave function $\Psi^{\mu}$ by solving Eq. (7) at the FCI or CIPSI level (see below) using Hartree-Fock orbitals. After a FCI or CIPSI calculation, the density $n_{\Psi^{\mu}}$ entering the short-range Hartree-exchange-correlation potential $\hat{\bar{V}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu}[n_{\Psi^{\mu}}]$ in Eq. (8) is updated and the procedure is iterated to achieve convergence with respect to the density (with an energy threshold of 10*-4* hartree). Depending on the type of long-range CI calculation used, we will refer to this method as “lrFCI+srPBE” or “lrCIPSI+srPBE”, or generically as “lrCI+srPBE”.

We then perform calculations according to Eq. (17) with the previously calculated wave function $\Psi^{\mu}$ and using either the srLDA approximation of Eq. (24) or the srPBEontop approximation of Eq. (25) for the short-range multideterminant correlation functional $\bar{E}_{\mathrm{{c,md}}}^{\mathrm{sr,}{\mu}}[n]$. We will refer to these calculations as “CI+Ec,md(srLDA)” and “CI+Ec,md(srPBEontop)”, where again CI can stand for either FCI or CIPSI.

We now briefly describe the CIPSI method as used here. The CIPSI method Bender and Davidson (1969); Huron et al. (1973); Buenker and Peyerimholf (1974); Buenker et al. (1981); Evangelisti et al. (1983); Rubio et al. (1986); Cimiraglia and Persico (1987); Harrison (1991); Angeli and Persico (1997); Angeli et al. (2000); Angeli and Cimiraglia (2001); Caffarel et al. (2014); Giner and Angeli (2015) is a selected CI which allows one to perform wave-function calculations at the near FCI level by keeping only the most important Slater determinants in a given FCI space. Starting from an initial guess for the wave function, $\Ket{\Psi^{\mu,(0)}}=\sum_{\mathrm{I}\in\mathcal{R}}c^{\mu}_{\textsc{i}}\Ket{\mathrm{I}}$ where $\Ket{\mathrm{I}}$ are Slater determinants in the reference variational space $\mathcal{R}$, the importance of a given Slater determinant $\Ket{\mathrm{K}}$ outside $\mathcal{R}$ is estimated using Epstein-Nesbet multireference perturbation theory. The second-order correction on the eigenvalue associated with the reference wave function $\ \mathcal{E}^{\mu,(0)}$ arising from the Slater determinant $\Ket{\mathrm{K}}$ is given by

[TABLE]

The variational space $\mathcal{R}$ is then enlarged by including the determinants associated with the largest perturbative corrections, and the procedure is iterated. In practice, the size of the variational space is doubled at each iteration until the magnitude of the total second-order Epstein-Nesbet correction on the eigenvalue, $\mathcal{E}^{\mu,(2)}=\sum_{\mathrm{K}}{\cal E}^{\mu,(2)}_{\text{K}}$, is smaller than a given threshold (10*-5* hartree or smaller). At a given iteration of the loop over the density $n_{\Psi^{\mu}}$ (entering the short-range potential $\hat{\bar{V}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu}[n_{\Psi^{\mu}}]$), we use the wave function obtained at the previous iteration as the starting guess for the CIPSI calculation. Thus, the variational space considered at the $i$th iteration is included in the variational space considered at the $(i+1)^{\text{th}}$ iteration, $\mathcal{R}^{(i)}\subset\mathcal{R}^{(i+1)}$. Note that in order to fully couple the RS-DFT calculation with the CIPSI method we should then add a perturbative correction to the total energy in Eqs. (9) or (17). However, in the present study, since all CIPSI calculations were iterated until we obtained a very small $\mathcal{E}^{\mu,(2)}$, we can neglect this perturbative correction to the total energy in comparison to the threshold used for converging the density.

All calculations were performed using correlation-consistent Dunning basis sets Dunning (1989); Woon and Dunning (1995); B. P. Prascher, D. E. Woon, K. A. Peterson, T. H. Dunning, Jr., and A. K. Wilson (2011) specified later.

IV Results and discussion

IV.1 Convergence with respect to the number of determinants

We first report in Fig. 1 a comparison of the convergence of the standard CIPSI and lrCIPSI+srPBE total variational energies as a function of the number of selected determinants for the Ne atom and for the Be2 molecule using the aug-cc-pCVQZ and aug-cc-pCVTZ basis sets, respectively, and correlating all the electrons in the CI calculations. For lrCIPSI+srPBE we use a range-separation parameter of $\mu=0.5\ \mathrm{bohr}^{-1}$. This figure clearly illustrates that the cuspless long-range wave function $\Psi^{\mu}$ of RS-DFT is much more compact than its Coulombic counterpart of standard WFT. Indeed, one sees that with a mere hundreds or thousands of determinants the lrCIPSI+srPBE total energy is already converged as much as the standard CIPSI total energy is with hundreds of thousands of determinants. This shows that the coupling of RS-DFT with a selected CI procedure such as the CIPSI method allows one to reduce by several orders of magnitude the dimension of the variational space of the wave function required to obtain a given accuracy.

IV.2 Total energies as a function of the range-separation parameter

We now discuss the accuracy of the total energy obtained with the different approximate RS-DFT schemes as a function of the range-separation parameter $\mu$. Figure 2 reports the results for the He and C atoms, and for the H2 molecule near dissociation as an example of a strongly correlated system. For He and H2, the calculations were performed at the FCI level using the cc-pVTZ basis set. For C, the calculations were performed at the CIPSI level using the cc-pCVTZ basis set and allowing core excitations.

IV.2.1 lrCI+srPBE total energy

As previously noted, RS-DFT reduces to KS-DFT for $\mu=0$ and to standard WFT for $\mu\to\infty$. The behavior of the lrCI+srPBE total energy as a function of $\mu$ is in agreement with these limits. Indeed, for $\mu\to 0$, the lrCI+srPBE energy is substantially above the CI energy and goes toward the KS-DFT energy obtained with the PBE exchange-correlation functional (KS-PBE, not shown). For $\mu\to\infty$, the lrCI+srPBE energy converges asymptotically to the standard CI energy.

For an optimal intermediate value of $\mu$, which is dependent on the system, the lrCI+srPBE total energy is comparable to the full-range CI total energy, or even more accurate in the case of He at the energy minimum. We must stress, however, that since the lrCI+srPBE total energy is not necessarily an upper bound of the exact energy, the value of $\mu$ minimizing the total energy cannot generally be considered as the optimal value of $\mu$. For the cases of C and H2, the lrCI+srPBE total energy is not more accurate (or only marginally) than the CI total energy, but the use of lrCI+srPBE allows one to obtain near FCI quality results with not too large a value of $\mu$ leading to a more compact wave function, as discussed in Sec. 1.

We note that the optimal value of $\mu$ required to obtain an accurate total energy is substantially larger for C than for He and H2. This is due to the contribution of the core spatial region of C which is associated to high densities and thus to small interelectronic distances. In order to have a part of the exchange-correlation energy of the core electrons of C treated via the long-range CI wave function, the long-range electron-electron interaction $w_{\mathrm{ee}}^{\mathrm{lr},\mu}(r_{12})$ must include the interaction between electrons at sufficiently small distances, i.e. $\mu$ must be sufficiently large.

IV.2.2 CI+Ec,md(srLDA) total energy

For $\mu\to\infty$, the CI+Ec,md(srLDA) total energy converges to the standard CI total energy, as was the case for lrCI+srPBE. However, the $\mu=0$ limit is different. Since the whole exchange energy is extracted from the wave function in the approach using the short-range multideterminantal correlation functional [Eq. (17)], the $\mu=0$ limit corresponds to a KS-DFT calculation with exact exchange. More precisely, at $\mu=0$, the CI+Ec,md(srLDA) energy reduces to $\braket{\Phi^{\textsc{pbe}}}{\hat{H}}{\Phi^{\textsc{pbe}}}+E_{\mathrm{c}}^{\textsc{lda}}[n_{\Phi^{\textsc{pbe}}}]$ where $\Phi^{\textsc{pbe}}$ is the KS single determinant obtained by solving the KS equation with the PBE exchange-correlation functional, and $E_{\mathrm{c}}^{\textsc{lda}}[n]$ is the standard LDA correlation functional. This explains why, for small values of $\mu$, CI+Ec,md(srLDA) is inaccurate. For He and C, it gives far too negative total energies because the well-known overestimation (in absolute value) of the correlation energy by the LDA functional by about a factor of 2 is not compensated by an approximate exchange functional as in standard KS-DFT. For H2 near dissociation, the missing static correlation effects makes the CI+Ec,md(srLDA) total energy too high for small values of $\mu$.

Due to its very inaccurate $\mu=0$ limit, the CI+Ec,md(srLDA) total energies tend to have large variations with respect to $\mu$ and become reasonably accurate only for values of $\mu$ similar to those required for lrCI+srPBE. Therefore, CI+Ec,md(srLDA) cannot be considered as an improvement over lrCI+srPBE. A better approximation must be used for the short-range multideterminantal correlation functional ${{\bar{E}_{\mathrm{{c,md}}}^{\mathrm{sr,}{\mu}}[n]}}$.

IV.2.3 CI+Ec,md(srPBEontop) total energy

As was the case for CI+Ec,md(srLDA), the CI+Ec,md(srPBEontop) total energy goes to the standard CI total energy for $\mu\to\infty$. For $\mu=0$, the CI+Ec,md(srPBEontop) energy reduces to the KS-PBE energy with exact exchange, i.e. $\braket{\Phi^{\textsc{pbe}}}{\hat{H}}{\Phi^{\textsc{pbe}}}+E_{\mathrm{c}}^{\textsc{pbe}}[n_{\Phi^{\textsc{pbe}}}]$ where $E_{\mathrm{c}}^{\textsc{pbe}}$ is the standard PBE correlation functional. By contrast, we remind that lrCI+srPBE reduces to standard KS-PBE at $\mu=0$. One must have this in mind when comparing CI+Ec,md(srPBEontop) and lrCI+srPBE at small $\mu$.

For He and C, it turns out that KS-PBE with exact exchange is more accurate than standard KS-PBE, which makes CI+Ec,md(srPBEontop) more accurate than lrCI+srPBE at small and intermediate $\mu$. Also, the fact that the PBE approximation to the KS correlation functional is a better approximation than the LDA approximation makes CI+Ec,md(srPBEontop) give much more accurate total energies at small and intermediate $\mu$ in comparison to CI+Ec,md(srLDA). For these weakly correlated systems, CI+Ec,md(srPBEontop) gives quite accurate total energies over the whole range of $\mu$.

For H2 near dissociation, because of the presence of static correlation effects, the total energy given by CI+Ec,md(srPBEontop) is much too high for small values of $\mu$. In particular, at $\mu=0$, we recover the known fact that KS-PBE with exact exchange gives a larger error than KS-PBE for strongly correlated systems. However, the error rapidly decreases with $\mu$, and the CI+Ec,md(srPBEontop) total energy converges to the accurate CI total energy significantly faster than both lrCI+srPBE and CI+Ec,md(srLDA). This must be due to the use of the on-top pair density which imposes the correct asymptotic behavior for $\mu\to\infty$.

We thus conclude that the srPBEontop approximation to the short-range multideterminant correlation functional ${{\bar{E}_{\mathrm{{c,md}}}^{\mathrm{sr,}{\mu}}[n]}}$ constitutes overall a large improvement over the srLDA approximation.

IV.3 Dissociation energy curves of the H2, Li2, and Be2 molecules

The dissociation energy curves of the homonuclear diatomic molecules H2, Li2, and Be2 are reported in Figs. 3, 4, and 5, respectively. These molecules cover different types of bonding and correlation effects. The RS-DFT calculations were performed at the frozen-core FCI level using the cc-pVTZ basis set for H2 and Li2 and at the frozen-core CIPSI level using the aug-cc-pVTZ basis set for Be2. We did not attempt to find an optimal value for the range-separation parameter and we simply used the common value of $\mu=0.5$ bohr*-1* Gerber and Ángyán (2005a). We did not try to remove the basis-set superposition error (BSSE) in the dissociation energy curve of the weakly bound Be2 molecule since the BSSE is known to be small for this system for frozen-core calculations with triple-zeta basis sets, and even more so with range separation Gerber and Ángyán (2005b); Reinhardt et al. (2018).

IV.3.1 H2 molecule

The electronic ground-state of the H2 molecule is one of the standard toy model of quantum chemistry owing to the range of correlation effects that it presents, from dynamic correlation at the equilibrium internuclear distance to static correlation at dissociation.

The KS-PBE total energy curve showed in Fig. 3 is a good example of the success and failure of standard KS-DFT with semilocal density approximations. KS-PBE gives an accurate energy near the equilibrium which illustrates the fact that the PBE approximation correctly describes dynamic correlation effects. By contrast, the KS-PBE results are far from being satisfying near the dissociation, which shows the incapacity of the PBE approximation to deal with static correlation effects.

The lrCI+srPBE method provides a way to partly correct the description of static correlation. Indeed, the long-range wave function accounts for part of the static correlation, while the srPBE functional accounts for the dynamic correlation. Thus, near the equilibrium, the lrCI+srPBE energy curve is essentially as accurate as the KS-PBE one (in fact slightly more accurate), and in the dissociation limit lrCI+srPBE greatly improves upon KS-PBE by giving an energy which correctly saturates. However, for the value of the range-separation parameter used ($\mu=0.5$ bohr*-1*), a substantial part of the electron-electron interaction is still taken into account via the short-range Hartree-exchange-correlation functional which leads to the important remaining error at dissociation.

The CI+Ec,md(srLDA) total energy curve is below the exact energy curve. This is of course due to the overestimation (in absolute value) of the short-range multideterminant correlation energy ${{\bar{E}_{\mathrm{{c,md}}}^{\mathrm{sr,}{\mu}}}}$ by the srLDA correlation functional. Still, the CI+Ec,md(srLDA) approach constitutes for this system a substantial improvement over KS-PBE and lrCI+srPBE, especially in terms of the relative shape of the dissociation curve.

We see that CI+Ec,md(srPBEontop) provides by far the most accurate total energy curve, either in terms of absolute energy or relative shape. The H2 molecule is simple enough to easily understand why our new srPBEontop functional gives accurate results. At dissociation, the two electrons are so far away from each other that the electron-electron interaction becomes negligible. Therefore, the exact long-range interacting Hamiltonian of RS-DFT $\hat{H}^{\mu}$ [Eq. (8)] becomes equivalent to the physical Hamiltonian $\hat{H}$, and consequently the exact long-range wave function $\Psi^{\mu}$ reduces to the exact ground-state wave function $\Psi$ of the system. Hence, at dissociation, the term $\braket{\Psi^{\mu}}{\hat{H}}{\Psi^{\mu}}$ in Eq. (17) should be equal to the exact energy, and the short-range multideterminant correlation energy ${{\bar{E}_{\mathrm{{c,md}}}^{\mathrm{sr,}{\mu}}[n]}}$ should vanish. This exact behavior is correctly recovered thanks to the on-top pair density. Indeed, at dissociation, the on-top pair density, $n_{2}(\textbf{r},\textbf{r})$ or $n_{2}^{\mu}(\textbf{r},\textbf{r})$, goes to zero since the two electrons are far away from each other, and it is easy to check from Eqs. (26) and (27) that this makes the srPBEontop correlation energy vanish. We thus see that the dependence on the on-top pair density is the key to obtain the correct dissociation limit. The reason why in practice CI+Ec,md(srPBEontop) still gives a small error at dissociation is that the short-range potential $\hat{\bar{V}}_{\mathrm{Hxc}}^{\mathrm{sr},\mu}[n]$ in Eq. (8) does not exactly vanish at dissociation due to the use of the PBE approximation. Consequently, the long-range wave function $\Psi^{\mu}$ does not exactly reduce to the exact wave function in the dissociation limit but is a good approximation to it.

IV.3.2 Li2 and Be2 molecules

We now consider the electronic ground-state energy curves of the Li2 and Be2 molecules. Although they are still relatively small systems, they raise more difficulties than H2, not only because of the increasing number of electrons but especially because of the more subtle mix between dynamic and static correlations that has to be described. Also, these two molecules are characterized by two different types of bond: Li2 is a strongly bonded molecule, while Be2 is a weakly bonded molecule with a very shallow well of only a few millihartree.

Here the distinction between dynamic correlation at the equilibrium geometry and static correlation at dissociation that exists in H2 is no longer valid. Static correlation effects are present at all internuclear distances for these molecules, which limits the accuracy of lrCI+srPBE total energies. Indeed, as shown in the left panels of Figs. 4 and 5, the lrCI+srPBE total energy is well above the exact energy for all internuclear distances. Using the short-range multideterminant correlation approach, we obtain results following the same trends observed for H2. For both Li2 and Be2, the CI+Ec,md(srLDA) total energy curve is far below the exact one, while the CI+Ec,md(srPBEontop) total energy curve is quite close to the exact one.

We note that, in the case of Li2, the FCI+Ec,md(srPBEontop) total energy curve, calculated with the cc-pVTZ basis set without core excitations, is much more accurate than the frozen-core FCI total energy curve calculated with the same basis set, and even slightly closer to the exact energy curve than the CIPSI total energy curve calculated with the cc-pCVTZ basis set and allowing core excitations. This is so because, in the FCI+Ec,md(srPBEontop) method, core correlation being a short-range effect is included in the srPBEontop functional. Thus, RS-DFT allows one to drop core excitations in the expansion of the wave function without losing accuracy, which is another important advantage in terms of computational cost.

We now discuss the relative dissociation energy curves, $E(R)-E(R\to\infty)$ where $R$ is the internuclear distance, shown in the right panels of Figs. 4 and 5. For Li2, all the methods tested here give almost the same relative dissociation energy curve and is very close to the exact relative energy curve. Thus, even though these methods give very different total energies, they all provide an accurate estimation of both the equilibrium distance and the dissociation energy.

For Be2, the different methods give more diverse relative dissociation energy curves. This is due to the fact that we are looking at a much smaller energy scale in comparison to Li2, and also to the fact that the Be2 bond involves a complex mix of correlation effects. The full-range frozen-core CIPSI calculation using the aug-cc-pVTZ basis set gives a substantially underestimated dissociation energy and a slightly overestimated equilibrium distance, i.e. it favors too much the separated atoms over the more correlated bonded molecule. This is due to the incompleteness of the basis set and possibly also to the missing of core excitations. On the contrary, lrCIPSI+srPBE largely overestimates the dissociation energy and slightly underestimates the equilibrium distance, i.e. it favors too much the bonded molecule over the dissociated atoms. In this case, the main source of error comes from the srPBE exchange-correlation functional (fractional-charge and/or fractional spin errors). The short-range multideterminant correlation approach gives quite good relative dissociation energy curves. The CI+Ec,md(srLDA) relative dissociation energy curve is a bit too high, in particular at long distances. The CI+Ec,md(srPBEontop) relative dissociation energy curve is the closest to the exact one, showing the srPBEontop functional properly accounts for differential correlation effects.

V Conclusion

In this work, we have developed a new approximation to the short-range multideterminant correlation functional ${{\bar{E}_{\mathrm{{c,md}}}^{\mathrm{sr,}{\mu}}[n]}}$ involved in the variant of RS-DFT given by Eq. (17). This approximation, named srPBEontop, is a local functional of the density, the density gradient, and the on-top pair density, which locally interpolates between the standard PBE correlation functional at vanishing range-separation parameter $\mu$ and the known exact asymptotic expansion of the functional at large $\mu$. By combining this srPBEontop correlation functional with (selected) CI calculations for the long-range wave function, one expects to obtain a multideterminant RS-DFT method which is essentially free from self-interaction errors and appropriately accounts for both short-range dynamic correlation and static correlation. This is supported by the accurate dissociation energy curves of the small but diversely correlated molecules H2, Li2, and Be2 that we have obtained with the multideterminant RS-DFT approach with the srPBEontop approximation.

Besides assessing the present method on more systems, possible future developments include adding the second-order CIPSI perturbative correction, performing self-consistent calculations with the srPBEontop approximation, combining this approximation with the recent local-$\mu$ approach of Ref. Giner et al., 2018, and calculating excited states for example using perturbation theory along the ground-state range-separated adiabatic connection Rebolini et al. (2015, 2018) or using ghost-interaction-corrected ensemble RS-DFT Senjean et al. (2016); Alam et al. (2016, 2017).

Acknowledgements.

We thank Pierre-François Loos and Cyrus Umrigar for comments on the manuscript.