Added comments and corrections.

This commit is contained in:
Bruno Senjean 2020-04-08 10:28:58 +02:00
parent 60c2cad1c2
commit 9199daf07c

View File

@ -133,7 +133,7 @@
Gross-Oliveira-Kohn (GOK) ensemble density-functional theory (GOK-DFT) is a time-independent formalism which allows to compute excitation energies via the derivative of the ensemble energy with respect to the weight of each excited state.
Contrary to the time-dependent version of density-functional theory (TD-DFT), double excitations can be easily computed within GOK-DFT.
However, to take full advantage of this formalism, one must have access to a \textit{weight-dependent} exchange-correlation functional in order to model the infamous derivative discontinuity contribution to the excitation energies.
In the present article, we discuss the construction of first-rung (\ie, local) weight-dependent exchange-correlation density-functional approximations for two-electron atomic and molecular systems (He and H$_2$) specifically designed for the computation of double excitations within GOK-DFT.
In the present article, we discuss the construction of first-rung (\ie, local) weight-dependent exchange-correlation density-functional approximations for two-electron atomic and molecular systems (He and H$_2$) \bruno{but it can be applied to larger systems as well right ? thanks to your shift} specifically designed for the computation of double excitations within GOK-DFT.
\end{abstract}
\maketitle
@ -208,8 +208,8 @@ where $\hH = \hT + \hWee + \hVne$ contains the kinetic, electron-electron and nu
\begin{eqnarray}
\hGam{\bw} = \sum_{I=0}^{\nEns - 1} \ew{I} \dyad*{\overline{\Psi}^{(I)}},
\end{eqnarray}
where $\lbrace \overline{\Psi}^{(I)} \rbrace_{1 \le I \le \nEns}$ is a set of $\nEns$ orthonormal trial wave functions.
The lower bound of Eq.~\eqref{eq:ens_energy} is reached when the set of wave functions correspond to the exact eigenstates of $\hH$, \ie, $\lbrace \overline{\Psi}^{(I)} \rbrace_{1 \le I \le \nEns} = \lbrace \Psi^{(I)} \rbrace_{1 \le I \le \nEns}$.
where $\lbrace \overline{\Psi}^{(I)} \rbrace_{0 \le I \le \nEns-1}$ is a set of $\nEns$ orthonormal trial wave functions.
The lower bound of Eq.~\eqref{eq:ens_energy} is reached when the set of wave functions correspond to the exact eigenstates of $\hH$, \ie, $\lbrace \overline{\Psi}^{(I)} \rbrace_{0 \le I \le \nEns-1} = \lbrace \Psi^{(I)} \rbrace_{0 \le I \le \nEns-1}$.
Multiplet degeneracies can be easily handled by assigning the same weight to the degenerate states.
One of the key feature of the GOK ensemble is that individual excitation energies are extracted from the ensemble energy via differentiation with respect to individual weights:
\begin{equation}
@ -234,7 +234,7 @@ $\Ts{\bw}[\n{}{}]$ is the noninteracting ensemble kinetic energy functional,
\begin{equation}
\hgam{\bw} = \sum_{I=0}^{\nEns-1} \ew{I} \dyad{\Det{I}{\bw}}
\end{equation}
is the density matrix operator, $\lbrace \Det{I}{\bw} \rbrace_{1 \le I \le \nEns}$ are single-determinant wave functions (or configuration state functions) built with KS orbitals $\lbrace \MO{p}{\bw}(\br{}) \rbrace_{1 \le p \le \nOrb}$, and
is the density matrix operator, $\lbrace \Det{I}{\bw} \rbrace_{0 \le I \le \nEns-1}$ are single-determinant wave functions (or configuration state functions) built with KS orbitals $\lbrace \MO{p}{\bw}(\br{}) \rbrace_{1 \le p \le \nOrb}$, and
\begin{equation}
\label{eq:exc_def}
\begin{split}
@ -258,10 +258,10 @@ From the GOK-DFT ensemble energy expression in Eq.~\eqref{eq:Ew-GOK}, we obtain
where
\begin{align}
\label{eq:nw}
\n{}{\bw}(\br{}) & = \sum_{I=0}^{\nEns-1} \ew{I} \n{}{(I)}(\br{}),
\n{}{\bw}(\br{}) & = \sum_{I=0}^{\nEns-1} \ew{I} \n{\Det{I}{\bw}}{}(\br{}),
\\
\label{eq:nI}
\n{}{(I)}(\br{}) & = \sum_{p}^{\nOrb} \ON{p}{(I)} [\MO{p}{\bw}(\br{})]^2
\n{\Det{I}{\bw}}{}(\br{}) & = \sum_{p}^{\nOrb} \ON{p}{(I)} [\MO{p}{\bw}(\br{})]^2
\end{align}
are the ensemble and individual one-electron densities, respectively,
\begin{equation}
@ -272,7 +272,7 @@ is the weight-dependent KS energy of state $I$, and $\eps{p}{\bw}$ is the KS orb
The latters are determined by solving the ensemble KS equation
\begin{equation}
\label{eq:eKS}
\qty( \hHc(\br{}) + \fdv{\E{\Hxc}{\bw}[\n{}{}]}{\n{}{}(\br{})}) \MO{p}{\bw}(\br{}) = \eps{p}{\bw} \MO{p}{\bw}(\br{}),
\qty( \hHc(\br{}) + \fdv{\E{\Hxc}{\bw}[\n{}{\bw}]}{\n{}{}(\br{})}) \MO{p}{\bw}(\br{}) = \eps{p}{\bw} \MO{p}{\bw}(\br{}),
\end{equation}
where $\hHc(\br{}) = -\nabla^2/2 + \vne(\br{})$, and
\begin{equation}
@ -286,7 +286,10 @@ where $\hHc(\br{}) = -\nabla^2/2 + \vne(\br{})$, and
\end{equation}
is the Hxc potential.
Equation \eqref{eq:dEdw} is our working equation for computing excitation energies from a practical point of view.
Note that, although we have dropped the weight-dependency in the individual densities $\n{}{(I)}(\br{})$ defined in Eq.~\eqref{eq:nI}, these do not match the \textit{exact} individual-state densities as the non-interacting KS ensemble is expected to reproduce the true interacting ensemble density $\n{}{\bw}(\br{})$ defined in Eq.~\eqref{eq:nw}, and not each individual density.
Note that the individual densities $\n{\Det{I}{\bw}}{}(\br{})$ defined in Eq.~\eqref{eq:nI} do not match the \textit{exact} individual-state densities as the non-interacting KS ensemble is expected to reproduce the true interacting ensemble density $\n{}{\bw}(\br{})$ defined in Eq.~\eqref{eq:nw}, and not each individual density.
Nevertheless,
these densities can still be extracted in principle exactly
from the KS ensemble as shown by Fromager~\cite{Fromager_2020}.
In the following, we adopt the usual decomposition, and write down the weight-dependent xc functional as
\begin{equation}
@ -305,9 +308,9 @@ For more details about the self-consistent implementation of GOK-DFT, we refer t
For all calculations, we use the aug-cc-pVXZ (X = D, T, and Q) Dunning's family of atomic basis sets.
Numerical quadratures are performed with the \texttt{numgrid} library using 194 angular points (Lebedev grid) and a radial precision of $10^{-6}$. \cite{Becke_1988,Lindh_2001}
This study deals only with spin-unpolarised systems, \ie, $\n{\uparrow}{} = \n{\downarrow}{} = \n{}{}/2$ (where $\n{\uparrow}{}$ and $\n{\downarrow}{}$ are the spin-up and spin-down electron densities).
Moreover, we restrict our study to the case of a two-state ensemble (\ie, $\nEns = 2$) where both the ground state ($I=0$) and the first doubly-excited state ($I=1$) are considered.
Although we will sometimes ``violate'' this variational constraint, we should have $0 \le \ew{} \le 1/2$ to ensure the GOK variational principle.
However, the limit $\ew{} = 1$ is of particular interest as it corresponds to a genuine saddle point of the KS equations, and match perfectly the results obtained wit the maximum overlap method developed by Gilbert, Gill and coworkers. \cite{Gilbert_2008,Barca_2018a,Barca_2018b}
Moreover, we restrict our study to the case of a two-state ensemble (\ie, $\nEns = 2$) where both the ground state ($I=0$ with weight $1 - \ew{}$) and the first doubly-excited state ($I=1$ with weight $\ew{}$) are considered.
Although we should have $0 \le \ew{} \le 1/2$ to ensure the GOK variational principle, we will sometimes ``violate'' this variational constraint.
Indeed, the limit $\ew{} = 1$ is of particular interest as it corresponds to a genuine saddle point of the KS equations, and match perfectly the results obtained with the maximum overlap method developed by Gilbert, Gill and coworkers. \cite{Gilbert_2008,Barca_2018a,Barca_2018b}
Moreover, the limits $\ew{} = 0$ and $\ew{} = 1$ are the only two weights for which there is no ghost-interaction error.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -315,7 +318,7 @@ Moreover, the limits $\ew{} = 0$ and $\ew{} = 1$ are the only two weights for wh
\label{sec:res}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
First, we compute the ensemble energy of the \ce{H2} molecule using the aug-cc-pVTZ basis set and the weight-independent Slater local exchange functional, \cite{Dirac_1930, Slater_1951} which is explicitly given by
First, we compute the ensemble energy of the \ce{H2} molecule using the aug-cc-pVTZ basis set and the weight-independent Slater local exchange functional, \cite{Dirac_1930, Slater_1951} \bruno{notée S51 sur les figures?} which is explicitly given by
\begin{align}
\e{\ex}{\text{S}}(\n{}{}) & = \Cx{} \n{}{1/3},
&
@ -323,7 +326,7 @@ First, we compute the ensemble energy of the \ce{H2} molecule using the aug-cc-p
\end{align}
The ensemble energy $\E{}{w}$ is depicted in Fig.~\ref{fig:Ew_H2} as a function of the weight $0 \le \ew{} \le 1$.
Because this exchange functional does not depend on the ensemble weight, there is no contribution from the ensemble derivative term [last term in Eq.~\eqref{eq:dEdw}].
As anticipated, $\E{}{w}$ is far from being linear, which means that the excitation energy obtained via the derivative of the local energy varies greatly with the weight of the double excitation (see Fig.~\ref{fig:Om_H2}).
As anticipated, $\E{}{w}$ is far from being linear, which means that the excitation energy obtained via the derivative of the local energy varies significantly with the weight of the double excitation (see Fig.~\ref{fig:Om_H2}).
Note that the exact xc correlation ensemble functional would yield a perfectly linear energy and the same value of the excitation energy independently of $\ew{}$.
\begin{figure}
@ -360,7 +363,7 @@ and
&
\gamma & = - 0.367\,189,
\end{align}
makes the ensemble much more linear (see Fig.~\ref{fig:Ew_H2}), and the excitation energy much more stable and closer to the full configuration interaction (FCI) reference of $28.75$ eV \cite{Barca_2018a} (see Fig.~\ref{fig:Om_H2})
makes the ensemble much more linear (see Fig.~\ref{fig:Ew_H2})\bruno{C'est celle notée ``GIC'' sur la figure ? Pourquoi pas MSFL ? A clarifier pour le lecteur}, and the excitation energy much more stable and closer to the full configuration interaction (FCI) reference of $28.75$ eV \cite{Barca_2018a} (see Fig.~\ref{fig:Om_H2})
As readily seen from Eq.~\eqref{eq:Cxw}, $\Cx{\ew{}}$ reduces to $\Cx{}$ for $\ew{} = 0$.
Note that we are not only using data from $\ew{} = 0$ to $\ew{} = 1/2$, but we also consider $1/2 < \ew{} \le 1$.
We ensure that the weight-dependent functional does not affect the two ghost-interaction-free limit at $\ew{} = 0$ and $1$.
@ -540,6 +543,7 @@ Hence, we employ a simple embedding scheme where the two-electron FUEG (the impu
The weight-dependence of the xc functional is then carried exclusively by the impurity [\ie, the functionals defined in Eqs.~\eqref{eq:exw} and \eqref{eq:ecw}], while the remaining effects are produced by the bath (\ie, the usual LDA xc functional).
Consistently with such a strategy, Eqs.~\eqref{eq:exw} and \eqref{eq:ecw} are ``centred'' on their corresponding jellium reference
\bruno{you commented the exchange part, why ?}
\begin{equation}
\label{eq:becw}
\be{\xc}{\ew{}}(\n{}{}) = (1-\ew{}) \be{\xc}{(0)}(\n{}{}) + \ew{} \be{\xc}{(1)}(\n{}{})
@ -596,19 +600,12 @@ This embedding procedure can be theoretically justified by the generalised adiab
(where $\E{\xc}{}[\n{}{}]$ is the usual ground-state xc functional) originally derived by Franck and Fromager. \cite{Franck_2014}
Within this in-principle-exact formalism, the (weight-dependent) xc energy of the ensemble is constructed from the (weight-independent) ground-state functional.
In the case of a homogeneous system (or equivalently within the LDA), substituting Eq.~\eqref{eq:dexcdw} into \eqref{eq:GACE} yields, in the case of a bi-ensemble, Eq.~\eqref{eq:eLDA}. \bruno{La formule me semble pas juste. Les exposants sont pas bons d'après moi. Pour le bi-ensemble, on devrait avoir
$\pdv{\E{\xc}{(1-\xi,\xi)}[\n{}{}]}{\xi}$ dans l'intégrale, non ? Can it be :
\begin{equation}
\label{eq:GACE2}
\E{\xc}{\bw}[\n{}{}]
= \E{\xc}{}[\n{}{}]
+ \sum_{I=0}^{\nEns-1} \int_0^{\ew{I}} \pdv{\E{\xc}{(0,\ldots,0,\xi,\ew{I+1},\ldots,\ew{\nEns-1})}[\n{}{}]}{\xi} d\xi
\end{equation}
?}
$\pdv{\E{\xc}{(1-\xi,\xi)}[\n{}{}]}{\xi}$ dans l'intégrale, non ?}
%%% TABLE I %%%
\begin{table*}
\caption{
Excitation energies (in eV) of \ce{H2} with $\RHH = 1.4$ bohr for various methods and basis sets.
Excitation energies (in eV) of \ce{H2} with $\RHH = 1.4$ bohr for various methods and basis sets. \bruno{Why you don't report results from the eLDA functional which is not system-specific like MSFL ? What is MSFL for the correlation part ? Is it what you referred to eLDA in the text ?}
\label{tab:Energies}
}
\begin{ruledtabular}