revised theory

This commit is contained in:
Pierre-Francois Loos 2020-04-24 22:21:47 +02:00
parent eff59ecfc2
commit 3b3b7e2eea

View File

@ -161,14 +161,14 @@ In the present article, we discuss the construction of first-rung (\textit{i.e.}
%%%%%%%%%%%%%%%%%%%%
\section{Introduction}
Time-dependent density-functional theory (TD-DFT) has been the dominant force in the calculation of excitation energies of molecular systems in the last two decades.\cite{Casida_1995,Ulrich_2012,Loos_2020a}
\titou{At a moderate computational cost} (at least compared to the other excited-state \textit{ab initio} methods), TD-DFT can provide accurate transition energies for low-lying excited states of organic molecules (see, for example, Ref.~\onlinecite{Dreuw_2005} and references therein).
At a moderate computational cost (at least compared to the other excited-state \textit{ab initio} methods), TD-DFT can provide accurate transition energies for low-lying excited states of organic molecules (see, for example, Ref.~\onlinecite{Dreuw_2005} and references therein).
\titou{Importantly, within the widely-used adiabatic approximation, setting up a TD-DFT calculation for a given system is an
almost pain-free process from a user perspective as the only (yet
essential) input variable is the choice of the
ground-state exchange-correlation (xc) functional.}
Similar to density-functional theory (DFT), \cite{Hohenberg_1964,Kohn_1965,ParrBook} TD-DFT is an in-principle exact theory which formal foundations relie on the Runge-Gross theorem. \cite{Runge_1984}
The Kohn-Sham (KS) \titou{formulation} of TD-DFT transfers the
The Kohn-Sham (KS) formulation of TD-DFT transfers the
complexity of the many-body problem to the xc functional thanks to a
judicious mapping between a time-dependent non-interacting reference
system and its interacting analog \titou{which have both
@ -235,23 +235,19 @@ Unless otherwise stated, atomic units are used throughout.
\section{Theory}
\label{sec:theo}
Let us consider a GOK ensemble of $\nEns$ electronic states with individual energies $\E{}{(0)} \le \ldots \le \E{}{(\nEns-1)}$, and (normalised) monotonically decreasing weights $\bw = (\ew{0},\ldots,\ew{M-1})$, \ie, $\sum_{I=0}^{\nEns-1} \ew{I} = 1$, and $\ew{0} \ge \ldots \ge \ew{\nEns-1}$.
\manu{For clarity, I usually exclude $\ew{0}$ from $\bw$ so that $\bw$
only contains the weights that are allowed to vary independently. One
should write explicitly $\ew{0}=1-\sum_{I=1}^{\nEns-1} \ew{I}$ and
define $\bw$ as $\bw = (\ew{1},\ldots,\ew{M-1})$}
Let us consider a GOK ensemble of $\nEns$ electronic states with individual energies $\E{}{(0)} \le \ldots \le \E{}{(\nEns-1)}$, and (normalised) monotonically decreasing weights $\bw = (\ew{1},\ldots,\ew{M-1})$, \ie, $\ew{0}=1-\sum_{I=1}^{\nEns-1} \ew{I}$, and $\ew{0} \ge \ldots \ge \ew{\nEns-1}$.
The corresponding ensemble energy
\begin{equation}
\E{}{\bw} = \sum_{I=0}^{\nEns-1} \ew{I} \E{}{(I)}
\end{equation}
fulfils \manu{can be obtained from?} the variational principle
can be obtained from the variational principle
as follows\cite{Gross_1988a}
\begin{eqnarray}\label{eq:ens_energy}
\E{}{\bw} = \min_{\hGam{\bw}} \Tr[\hGam{\bw} \hH],
\end{eqnarray}
where $\hH = \hT + \hWee + \hVne$ contains the kinetic,
electron-electron and nuclei-electron interaction potential operators,
respectively, $\Tr$ denotes the trace and $\hGam{\bw}$ is a trial
respectively, $\Tr$ denotes the trace, and $\hGam{\bw}$ is a trial
density matrix operator of the form
\begin{eqnarray}
\hGam{\bw} = \sum_{I=0}^{\nEns - 1} \ew{I} \dyad*{\overline{\Psi}^{(I)}},
@ -259,7 +255,7 @@ density matrix operator of the form
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 \cite{Gross_1988b}.
weight to the degenerate states. \cite{Gross_1988b}
One of the key feature of the GOK ensemble is that individual excitation
energies can be extracted from the ensemble energy via differentiation with respect to individual weights:
\begin{equation}\label{eq:diff_Ew}
@ -274,54 +270,49 @@ Turning to GOK-DFT, the extension of the Hohenberg--Kohn theorem to ensembles al
where $\vne(\br{})$ is the external potential and $\F{}{\bw}[\n{}{}]$ is the universal ensemble functional
(the weight-dependent analog of the Hohenberg--Kohn universal functional for ensembles).
In the KS formulation, this functional can be decomposed as
\begin{equation}
\F{}{\bw}[\n{}{}]
= \Ts{\bw}[\n{}{}] + \E{\Hxc}{\bw}[\n{}{}]
= \Tr[ \hgam{\bw} \hT ] + \Tr[ \hgam{\bw} \hWee ],
\end{equation}
\manu{The above equation is wrong (the correlation is missing) and the
notations are ambiguous. I should also say that Tim does not like the
original separation into H and xc. I propose the following reformulation
to get everyone satisfied. I also reorganized the theory for clarity.
%\begin{equation}
% \F{}{\bw}[\n{}{}]
% = \Ts{\bw}[\n{}{}] + \E{\Hxc}{\bw}[\n{}{}]
% = \Tr[ \hgam{\bw} \hT ] + \Tr[ \hgam{\bw} \hWee ],
%\end{equation}
%\manu{The above equation is wrong (the correlation is missing) and the
%notations are ambiguous. I should also say that Tim does not like the
%original separation into H and xc. I propose the following reformulation
%to get everyone satisfied. I also reorganized the theory for clarity.
\begin{equation}\label{eq:FGOK_decomp}
\F{}{\bw}[\n{}{}]
= \Tr[ \hgamdens{\bw} \hT ]+ \E{\Ha}{}[\n{}{}]+\E{\xc}{\bw}[\n{}{}],
= \Tr{ \hgamdens{\bw} \hT }+ \E{\Ha}{}[\n{}{}]+\E{\xc}{\bw}[\n{}{}],
\end{equation}
}
where
\manuf{$\Tr[ \hgamdens{\bw} \hT ]=\Ts{\bw}[\n{}{}]$} is the noninteracting ensemble kinetic energy functional,
$\Tr{ \hgamdens{\bw} \hT } =\Ts{\bw}[\n{}{}]$ is the noninteracting ensemble kinetic energy functional,
\begin{equation}
\hgam{\bw}[n] = \sum_{I=0}^{\nEns-1} \ew{I} \dyad{\Det{I}{\bw}[n]}
\end{equation}
is the \manuf{KS density-functional} density matrix operator, and $\lbrace
is the KS density-functional density matrix operator, and $\lbrace
\Det{I}{\bw}[n] \rbrace_{0 \le I \le \nEns-1}$ are single-determinant
wave functions (or configuration state functions). \manuf{Their
dependence on the density
is determined from the ensemble density
constraint
wave functions (or configuration state functions).
Their dependence on the density is determined from the ensemble density constraint
\begin{equation}
\sum_{I=0}^{\nEns-1} \ew{I} n_{\Det{I}{\bw}[n]}(\br)=n(\br).
\sum_{I=0}^{\nEns-1} \ew{I} \n{\Det{I}{\bw}[n]}{}(\br) = \n{}{}(\br).
\end{equation}
Note that the original decomposition \cite{Gross_1988b} shown in Eq.~(\ref{eq:FGOK_decomp}), where the
Note that the original decomposition \cite{Gross_1988b} shown in Eq.~\eqref{eq:FGOK_decomp}, where the
conventional (weight-independent) Hartree functional
\beq
\E{\Ha}{}[\n{}{}]=\frac{1}{2} \iint \frac{\n{}{}(\br{}) \n{}{}(\br{}')}{\abs{\br{}-\br{}'}} d\br{} d\br{}'
\E{\Ha}{}[\n{}{}]=\frac{1}{2} \iint \frac{\n{}{}(\br{}) \n{}{}(\br{}')}{\abs{\br{}-\br{}'}} d\br{} d\br{}'
\eeq
is separated
from the (weight-dependent) exchange-correlation (xc) functional, is
formally exact. In practice, the use of such a decomposition might be
problematic as inserting an ensemble density into $\E{\Ha}{}[\n{}{}]$
causes the infamous ghost-interaction error \cite{Gidopoulos_2002,
Pastorczak_2014, Alam_2016, Alam_2017, Gould_2017}. The latter should in
principle be removed by the exchange component of the ensemble xc
functional $\E{\xc}{\bw}[\n{}{}]\equiv
\E{\ex}{\bw}[\n{}{}]+\E{\co}{\bw}[\n{}{}]$, as readily seen from the
exact expression
causes the infamous ghost-interaction error. \cite{Gidopoulos_2002,Pastorczak_2014, Alam_2016, Alam_2017, Gould_2017}
The latter should in principle be removed by the exchange component of the ensemble xc functional
$\E{\xc}{\bw}[\n{}{}] \equiv \E{\ex}{\bw}[\n{}{}] + \E{\co}{\bw}[\n{}{}]$,
as readily seen from the exact expression
\beq
\E{\ex}{\bw}[\n{}{}]=\sum_{I=0}^{\nEns-1} \ew{I}\bra{\Det{I}{\bw}[n]}\hat{W}_{\rm ee}\ket{\Det{I}{\bw}[n]}
-\E{\Ha}{}[\n{}{}].
\E{\ex}{\bw}[\n{}{}]
= \sum_{I=0}^{\nEns-1} \ew{I}\mel{\Det{I}{\bw}[\n{}{}]}{\hWee}{\Det{I}{\bw}[\n{}{}]} - \E{\Ha}{}[\n{}{}].
\eeq
The minimum in Eq.~(\ref{eq:Ew-GOK}) is reached when the density $n$
The minimum in Eq.~\eqref{eq:Ew-GOK} is reached when the density $n$
equals the exact ensemble one
\beq\label{eq:nw}
n^{\bw}(\br)=\sum_{I=0}^{\nEns-1}
@ -341,7 +332,7 @@ result, the orbitals
$\lbrace \MO{p}{\bw}(\br{}) \rbrace_{1 \le p \le
\nOrb}$ from which the KS
wavefunctions $\left\{\Det{I}{\bw}\left[n^{\bw}\right]\right\}_{0\leq
I\leq M-1}$ are constructed can be obtained by solving the following ensemble KS equation
I\leq \nEns-1}$ are constructed can be obtained by solving the following ensemble KS equation
\begin{equation}
\label{eq:eKS}
\qty{ \hHc(\br{}) + \fdv{\E{\Hxc}{\bw}[\n{}{\bw}]}{\n{}{}(\br{})}} \MO{p}{\bw}(\br{}) = \eps{p}{\bw} \MO{p}{\bw}(\br{}),
@ -355,7 +346,7 @@ where $\hHc(\br{}) = -\nabla^2/2 + \vne(\br{})$, and
+ \fdv{\E{\xc}{\bw}[\n{}{}]}{\n{}{}(\br{})}.
\end{equation}
The ensemble density can be obtained directly (and exactly, if no
approximation is made) from those orbitals:
approximation is made) from these orbitals, \ie,
\beq\label{eq:ens_KS_dens}
\n{}{\bw}(\br{})=\sum_{I=0}^{\nEns-1} \ew{I}\left(\sum_{p}^{\nOrb}
\ON{p}{(I)} [\MO{p}{\bw}(\br{})]^2\right),
@ -363,19 +354,19 @@ approximation is made) from those orbitals:
where $\ON{p}{(I)}$ denotes the occupation of $\MO{p}{\bw}(\br{})$ in
the $I$th KS wave function $\Det{I}{\bw}\left[n^{\bw}\right]$. Turning
to the excitation energies, they can be extracted from the
density-functional ensemble as follows [see Eqs. ({\ref{eq:diff_Ew}})
and ({\ref{eq:Ew-GOK}}) and Refs.
\cite{Gross_1988b,Deur_2019}]:
density-functional ensemble as follows [see Eqs.~\eqref{eq:diff_Ew}
and \eqref{eq:Ew-GOK} and Refs.~\onlinecite{Gross_1988b,Deur_2019}]:
\beq
\label{eq:dEdw}
\Omega^{(I)}= \Eps{I}{\bw} - \Eps{0}{\bw} + \left. \pdv{\E{\xc}{\bw}[\n{}{}]}{\ew{I}} \right|_{\n{}{} = \n{}{\bw}},
\label{eq:dEdw}
\Omega^{(I)}= \Eps{I}{\bw} - \Eps{0}{\bw} + \left. \pdv{\E{\xc}{\bw}[\n{}{}]}{\ew{I}} \right|_{\n{}{} = \n{}{\bw}},
\eeq
where
\begin{equation}
\label{eq:KS-energy}
\Eps{I}{\bw} = \sum_{p}^{\nOrb} \ON{p}{(I)} \eps{p}{\bw}
\end{equation}
is the energy of the $I$th KS state.\\
is the energy of the $I$th KS state.
Equation \eqref{eq:dEdw} is our working equation for computing excitation energies from a practical point of view.
Note that the individual KS densities
$\n{\Det{I}{\bw}\left[n^{\bw}\right]}{}(\br{})=\sum_{p}^{\nOrb}
@ -383,14 +374,14 @@ $\n{\Det{I}{\bw}\left[n^{\bw}\right]}{}(\br{})=\sum_{p}^{\nOrb}
not necessarily match the \textit{exact} (interacting) individual-state
densities $n_{\Psi_I}(\br)$ 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}.\\
exactly from the KS ensemble as shown by Fromager. \cite{Fromager_2020}
In the following, we will work at the (weight-dependent) LDA
level of approximation, \ie
\beq
\E{\xc}{\bw}[\n{}{}]
&\overset{\rm LDA}{\approx}&
\int \e{\xc}{\bw}(\n{}{}(\br{})) \n{}{}(\br{}) d\br{}
\int \e{\xc}{\bw}(\n{}{}(\br{})) \n{}{}(\br{}) d\br{},
\\
\fdv{\E{\xc}{\bw}[\n{}{}]}{\n{}{}(\br{})}
&\overset{\rm LDA}{\approx}&
@ -403,11 +394,8 @@ We will also adopt the usual decomposition, and write down the weight-dependent
where $\e{\ex}{\bw{}}(\n{}{})$ and $\e{\co}{\bw{}}(\n{}{})$ are the
weight-dependent density-functional exchange and correlation energies
per particle, respectively.
}%%%%%% end manuf
The explicit construction of these functionals is discussed at length in Sec.~\ref{sec:res}.
\manu{Maybe we should say a little bit more about how we will design
such approximations, or just say the design of these functionals will be
presented in the following...}
%%%%%%%%%%%%%%%%
%%%%%%% Manu: stuff that I removed from the first version %%%%%
\iffalse%%%%
@ -477,8 +465,7 @@ is the Hxc potential, with
\section{Computational details}
\label{sec:compdet}
The self-consistent GOK-DFT calculations \manuf{[see Eqs.~(\ref{eq:eKS})
and (\ref{eq:ens_KS_dens})]} have been performed in a restricted formalism with the \texttt{QuAcK} software, \cite{QuAcK} which is freely available on \texttt{github}, and where the present weight-dependent functionals have been implemented.
The self-consistent GOK-DFT calculations [see Eqs.~\eqref{eq:eKS} and \eqref{eq:ens_KS_dens}] have been performed in a restricted formalism with the \texttt{QuAcK} software, \cite{QuAcK} which is freely available on \texttt{github}, and where the present weight-dependent functionals have been implemented.
For more details about the self-consistent implementation of GOK-DFT, we refer the interested reader to Ref.~\onlinecite{Loos_2020} where additional technical details can be found.
For all calculations, we use the aug-cc-pVXZ (X = D, T, Q, and 5) Dunning family of atomic basis sets. \cite{Dunning_1989,Kendall_1992,Woon_1994}
Numerical quadratures are performed with the \texttt{numgrid} library \cite{numgrid} using 194 angular points (Lebedev grid) and a radial precision of $10^{-6}$. \cite{Becke_1988b,Lindh_2001}
@ -490,7 +477,7 @@ Indeed, the limit $\ew{} = 1$ is of particular interest as it corresponds to a g
it stands a little bit beyond the theory discussed previously. What you
are looking at in the range $1/2\leq \ew{}\leq 1$ are, indeed, other
stationary points (than the minimizing ones) of the density matrix
operator functional in Eq.~(\ref{eq:min_KS_DM}). I would say that we
operator functional in Eq.~\eqref{eq:min_KS_DM}. I would say that we
look at these solutions for analysis purposes. I personally never looked
(formally) at these solutions and their physical meaning. One should clearly
mention that applying GOK-DFT in this range of weights would simply
@ -541,6 +528,10 @@ As anticipated, $\E{}{\ew{}}$ is far from being linear, which means that the exc
Taking as a reference the full configuration interaction (FCI) value of $28.75$ eV obtained with the aug-mcc-pV8Z basis set, \cite{Barca_2018a} one can see that the excitation energy varies by more than $8$ eV from $\ew{} = 0$ to $1/2$.
Note that the exact xc ensemble functional would yield a perfectly linear energy and, hence, the same value of the excitation energy independently of $\ew{}$.
\titou{The eDFT purist is going to surprised to see that we left out the singly-excited states $\sigma_g \sigma_u$ from the ensemble as this state is lower in energy than the doubly-excited state of configuration $1\sigma_u^2$.
As we wish to stick with a restricted formalism, the single excitation is naturally left out of the ensemble.
However, as a sanity check, we have tried to introduce the single excitations as well.}
\begin{figure}
\includegraphics[width=\linewidth]{Ew_H2}
\caption{
@ -762,7 +753,7 @@ As shown in Fig.~\ref{fig:Ew_H2}, the GIC-SeVWN5 is slightly less concave than i
For a more qualitative picture, Table \ref{tab:BigTab_H2} reports excitation energies for various methods and basis sets.
In particular, we report the excitation energies obtained with GOK-DFT
in the zero-weight limit (\ie, $\ew{} = 0$) and for the equi-ensemble
(\ie, $\ew{} = 1/2$). \manu{Maybe we should refer to Eq.~(\ref{eq:dEdw}) for clarity.}
(\ie, $\ew{} = 1/2$). \manu{Maybe we should refer to Eq.~\eqref{eq:dEdw} for clarity.}
For comparison purposes, we also report the linear interpolation method (LIM) excitation energies, \cite{Senjean_2015,Senjean_2016}
a pragmatic way of getting weight-independent
excitation energies defined as