This commit is contained in:
Pierre-Francois Loos 2019-11-13 11:32:09 +01:00
parent 2794738908
commit 52d2970a29
3 changed files with 1820 additions and 142 deletions

1718
FarDFT/FarDFT.nb Normal file

File diff suppressed because it is too large Load Diff

View File

@ -1,7 +1,7 @@
%% This BibTeX bibliography file was created using BibDesk.
%% http://bibdesk.sourceforge.net/
%% Created for Pierre-Francois Loos at 2019-09-09 09:45:30 +0200
%% Created for Pierre-Francois Loos at 2019-11-12 22:50:55 +0100
%% Saved with string encoding Unicode (UTF-8)

View File

@ -110,7 +110,7 @@
\begin{document}
\title{Weight-dependent local-density approximation for molecules}
\title{Weight-dependent exchange-correlation functionals for molecules: the local-density approximation}
\author{Clotilde \surname{Marut}}
\affiliation{\LCPQ}
@ -119,7 +119,11 @@
\affiliation{\LCPQ}
\begin{abstract}
Here is a nice abstract.
We report a first generation of local, weight-dependent exchange-correlation density-functional approximations (DFAs) for molecules.
These density-functional approximations for ensembles (eDFAs) incorporate information about both ground and excited states in the context of density-functional theory for ensembles (eDFT).
They are specially designed for the computation of double excitations within eDFT, and can be seen as a natural extension of the ubiquitous local-density approximation (LDA) to ensembles.
The resulting eDFAs, dubbed eLDA, which are based on both finite and infinite uniform electron gas models, automatically incorporate the infamous derivative discontinuity contributions to the excitation energies through their explicit ensemble weight dependence.
Their accuracy is illustrated by computing on the prototypical H$_2$ molecule.
\end{abstract}
\maketitle
@ -130,6 +134,17 @@ Here is a nice abstract.
\section{Introduction}
Here is a nice introduction.
%Most of the standard local and semi-local DFAs rely on the infinite uniform electron gas (UEG) model (also known as jellium). \cite{ParrBook, Loos_2016}
%One major drawback of the jellium paradigm, when it comes to develop eDFAs, is that the ground and excited states cannot be easily identified like in a molecule. \cite{Gill_2012, Loos_2012, Loos_2014a, Loos_2014b, Agboola_2015, Loos_2017a}
%Moreover, because the infinite UEG model is a metal, it is gapless, which means that both the fundamental and optical gaps are zero.
%From this point of view, using finite UEGs \cite{Loos_2011b, Gill_2012} (which have, like an atom, discrete energy levels) to construct eDFAs can be seen as more relevant. \cite{Loos_2014a, Loos_2014b, Loos_2017a}
%
%As a finite uniform electron gas, we consider the glomium model in which electrons move on the surface of a four-dimensional sphere (also known as a glome). \cite{Loos_2011b, Agboola_2015}
%The most appealing feature of glomium (regarding the development of functionals in the context of eDFT) is the fact that both the ground- and excited-state densities have uniform densities.
%As a result, the ensemble density will remain constant (and uniform) as the ensemble weights vary.
%This is a necessary condition for being able to model derivative discontinuities.
%%%%%%%%%%%%%%%%%%%%
%%% THEORY %%%
%%%%%%%%%%%%%%%%%%%%
@ -141,56 +156,24 @@ Here is the theory.
%%%%%%%%%%%%%%%%%%
\section{Functional}
We decompose the weight-dependent functional as
We adopt the usual decomposition, and write down the weight-dependent exchange-correlation functional as
\begin{equation}
\be{Hxc}{\bw}(\n{}{}) = \be{Hx}{\bw}(\n{}{}) + \be{c}{\bw}(\n{}{}),
\e{xc}{\ew{}}(\n{}{}) = \e{x}{\ew{}}(\n{}{}) + \e{c}{\ew{}}(\n{}{}),
\end{equation}
where $\be{Hx}{\bw}(\n{}{})$ is a weight-dependent Hartree-exchange functional designed to correct the ghost interaction \cite{Gidopoulos_2002, Pastorczak_2014, Alam_2016, Alam_2017, Gould_2017} [see Subsec.~\ref{sec:GIC}] and $\be{c}{\bw}(\n{}{})$ is a weight-dependent correlation functional [see Subsec.~\ref{sec:Ec}].
where $\e{x}{\ew{}}(\n{}{})$ and $\e{c}{\ew{}}(\n{}{})$ are the weight-dependent exchange and correlation functionals, respectively.
The construction of these two functionals is described below.
Note that, because we consider strict 1D systems, one cannot decompose further the Hartree-exchange contribution as each component diverges independently but their sum is finite. \cite{Astrakharchik_2011, Lee_2011a, Loos_2012, Loos_2013, Loos_2013a}
Here, we restrict our study to spin-unpolarized systems, i.e., closed-shell systems.
Most of the standard local and semi-local DFAs rely on the infinite uniform electron gas (UEG) model (also known as jellium). \cite{ParrBook, Loos_2016}
One major drawback of the jellium paradigm, when it comes to develop eDFAs, is that the ground and excited states cannot be easily identified like in a molecule. \cite{Gill_2012, Loos_2012, Loos_2014a, Loos_2014b, Agboola_2015, Loos_2017a}
Moreover, because the infinite UEG model is a metal, it is gapless, which means that both the fundamental and optical gaps are zero.
From this point of view, using finite UEGs \cite{Loos_2011b, Gill_2012} (which have, like an atom, discrete energy levels) to construct eDFAs can be seen as more relevant. \cite{Loos_2014a, Loos_2014b, Loos_2017a}
Here, we propose to construct a weight-dependent eDFA for the calculations of excited states in 1D systems.
As a finite uniform electron gas, we consider the ringium model in which electrons move on a perfect ring (i.e., a circle). \cite{Loos_2012, Loos_2013a, Loos_2014b}
The most appealing feature of ringium (regarding the development of functionals in the context of eDFT) is the fact that both ground- and excited-state densities are uniform.
As a result, the ensemble density will remain constant (and uniform) as the ensemble weights vary.
This is a necessary condition for being able to model derivative discontinuities.
The present weight-dependent eDFA is specifically designed for the calculation of excitation energies within eDFT.
As mentioned previously, we consider a three-state ensemble including the ground state ($I=0$), the first singly-excited state ($I=1$), and the first doubly-excited state ($I=2$) of the (spin-polarized) two-electron ringium system.
All these states have the same (uniform) density $\n{}{} = 2/(2\pi R)$ where $R$ is the radius of the ring where the electrons are confined.
We refer the interested reader to Refs.~\onlinecite{Loos_2012, Loos_2013a, Loos_2014b} for more details about this paradigm.
The present weight-dependent eDFA is specifically designed for the calculation of double excitations within eDFT.
As mentioned previously, we consider a two-state ensemble including the ground state ($I=0$), and the first doubly-excited state ($I=1$) of the two-electron glomium system.
All these states have the same (uniform) density $\n{}{} = 2/(4\pi/3 R^3)$ where $R$ is the radius of the glome where the electrons are confined.
We refer the interested reader to Refs.~\onlinecite{Loos_2011b} for more details about this paradigm.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{Ghost-interaction correction}
\label{sec:GIC}
\subsection{Weight-dependent exchange functional}
\label{sec:Ex}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
The GIC weight-dependent Hartree-exchange functional is defined as
\begin{multline}
\be{Hx}{\bw}(\n{}{\bw}) = (1-\sum_{I>0} \ew{I}) \be{Hx}{}(\n{}{(0)}) + \sum_{I>0} \ew{I} \be{Hx}{}(\n{}{(I)})
\\
- \be{Hx}{(I)}(\n{}{\bw}),
\end{multline}
where
\begin{equation}
\be{Hx}{}(\n{}{}) = \iint \frac{\n{}{}(\br_1) \n{}{}(\br_2) - \n{}{}(\br_1,\br_2)^2}{r_{12}} d\br_1 d\br_2,
\end{equation}
and
\begin{equation}
\n{}{(I)}(\omega) = (\pi R)^{-1} \cos[(I+1) \omega/2]
\end{equation}
is the first-order density matrix with $\omega$ the interelectronic angle.
It yields
\begin{equation}
\be{Hx}{}(\n{}{}) = \n{}{} \qty[ a_1 \ew{1} (\ew{1} - 1) + a_2 \ew{1} \ew{2} + a_3 \ew{2} (\ew{2} - 1)],
\end{equation}
with $a_1 = 2 \ln 2 - 1/3$, $a_2 = 8/3$ and $a_3 = 32/15$.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{Weight-dependent correlation functional}
\label{sec:Ec}
@ -199,97 +182,23 @@ with $a_1 = 2 \ln 2 - 1/3$, $a_2 = 8/3$ and $a_3 = 32/15$.
Based on highly-accurate calculations (see below), one can write down, for each state, an accurate analytical expression of the reduced (i.e., per electron) correlation energy \cite{Loos_2013a, Loos_2014a} via the following Pad\'e approximant
\begin{equation}
\label{eq:ec}
\e{c}{(I)}(\n{}{}) = \frac{c_1^{(I)}\,\n{}{}}{\n{}{} + c_2^{(I)} \sqrt{\n{}{}} + c_3^{(I)}},
\e{xc}{(I)}(\n{}{}) = \frac{c_1^{(I)}\,\n{}{}}{\n{}{} + c_2^{(I)} \sqrt{\n{}{}} + c_3^{(I)}},
\end{equation}
where the $c_k^{(I)}$'s are state-specific fitting parameters, which are provided in Table \ref{tab:OG_func}.
The value of $c_1^{(I)}$ is obtained via the exact high-density expansion of the correlation energy. \cite{Loos_2013a, Loos_2014a}
Equation \eqref{eq:ec} provides three state-specific correlation DFAs based on a two-electron system.
Combining these, one can build a three-state weight-dependent correlation eDFA:
Equation \eqref{eq:ec} provides two state-specific correlation DFAs based on a two-electron system.
Combining these, one can build a two-state weight-dependent correlation eDFA:
\begin{equation}
\label{eq:ecw}
\e{c}{\bw}(\n{}{}) = (1-\ew{1}-\ew{2}) \e{c}{(0)}(\n{}{}) + \ew{1} \e{c}{(1)}(\n{}{}) + \ew{2} \e{c}{(2)}(\n{}{}).
\e{c}{\ew{}}(\n{}{}) = (1-\ew{}) \e{c}{(0)}(\n{}{}) + \ew{} \e{c}{(1)}(\n{}{}).
\end{equation}
%%% TABLE 1 %%%
\begin{table*}
\caption{
\label{tab:OG_func}
Parameters of the correlation DFAs defined in Eq.~\eqref{eq:ec}.}
% \begin{ruledtabular}
\begin{tabular}{lcddd}
\hline\hline
State & $I$ & \tabc{$c_1^{(I)}$} & \tabc{$c_2^{(I)}$} & \tabc{$c_3^{(I)}$} \\
\hline
Ground state & $0$ & -0.0137078 & 0.0538982 & 0.0751740 \\
Singly-excited state & $1$ & -0.0238184 & 0.00413142 & 0.0568648 \\
Doubly-excited state & $2$ & -0.00935749 & -0.0261936 & 0.0336645 \\
\hline\hline
\end{tabular}
% \end{ruledtabular}
\end{table*}
%%% %%% %%% %%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{LDA-centered functional}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
In order to make the two-electron-based eDFA defined in Eq.~\eqref{eq:ecw} more universal and to ``center'' it on the jellium reference (as commonly done in DFT), we propose to \emph{shift} it as follows:
\begin{equation}
\label{eq:becw}
\be{c}{\bw}(\n{}{}) = (1-\ew{1}-\ew{2}) \be{c}{(0)}(\n{}{}) + \ew{1} \be{c}{(1)}(\n{}{}) + \ew{2} \be{c}{(2)}(\n{}{}),
\end{equation}
where
\begin{equation}
\be{c}{(I)}(\n{}{}) = \e{c}{(I)}(\n{}{}) + \e{c}{\text{LDA}}(\n{}{}) - \e{c}{(0)}(\n{}{}).
\end{equation}
The local-density approximation (LDA) correlation functional,
\begin{equation}
\e{c}{\text{LDA}}(\n{}{}) = c_1^\text{LDA} \, F\qty[1,\frac{3}{2},c_3^\text{LDA}, \frac{c_1^\text{LDA}(1-c_3^\text{LDA})}{c_2^\text{LDA}} {\n{}{}}^{-1}],
\end{equation}
specifically designed for 1D systems in Ref.~\onlinecite{Loos_2013} as been used, where $F(a,b,c,x)$ is the Gauss hypergeometric function, \cite{NISTbook} and
\begin{align}
c_1^\text{LDA} & = - \frac{\pi^2}{360},
&
c_2^\text{LDA} & = \frac{3}{4} - \frac{\ln{2\pi}}{2},
&
c_3^\text{LDA} & = 2.408779.
\end{align}
Equation \eqref{eq:becw} can be recast
\begin{equation}
\label{eq:eLDA}
\begin{split}
\be{c}{\bw}(\n{}{})
& = \e{c}{\text{LDA}}(\n{}{})
\\
& + \ew{1} \qty[\e{c}{(1)}(\n{}{})-\e{c}{(0)}(\n{}{})] + \ew{2} \qty[\e{c}{(2)}(\n{}{})-\e{c}{(0)}(\n{}{})],
\end{split}
\end{equation}
which nicely highlights the centrality of the LDA in the present eDFA.
In particular, $\be{c}{(0,0)}(\n{}{}) = \e{c}{\text{LDA}}(\n{}{})$.
Consequently, in the following, we name this correlation functional ``eLDA'' as it is a natural extension of the LDA for ensembles.
This procedure can be theoretically justified by the generalized adiabatic connection formalism for ensembles (GACE) which was originally derived by Franck and Fromager. \cite{Franck_2014}
Within this in-principle-exact formalism, the (weight-dependent) correlation energy of the ensemble is constructed from the (weight-independent) ground-state functional (such as the LDA), yielding Eq.~\eqref{eq:eLDA}.
This is a crucial point as we intend to incorporate into standard functionals (which are ``universal'' in the sense that they do not depend on the number of electrons) information about excited states that will be extracted from finite systems (whose properties may depend on the number of electrons).
Finally, we note that, by construction,
\begin{equation}
\left. \pdv{\be{c}{\bw}[\n{}{}]}{\ew{J}}\right|_{\n{}{} = \n{}{\bw}(\br)} = \be{c}{(J)}[\n{}{\bw}(\br)] - \be{c}{(0)}[\n{}{\bw}(\br)].
\end{equation}
\alert{As shown by Gould and Pittalis, comment on density- and and state-driven errors. \cite{Gould_2019}}
The density-functional approximations designed in this manuscript are based on highly-accurate energies for the ground state ($I=0$), the first singly-excited state ($I=1$), and the first doubly-excited state ($I=2$) of the (spin-polarized) two-electron glomium model.
We refer the interested reader to Refs.~\onlinecite{} for more details about this paradigm.
The reduced (i.e.~per electron) HF energy for these three states is:
The reduced (i.e.~per electron) HF energy for these two states is:
\begin{subequations}
\begin{align}
\e{HF}{(0)}(n) & = \frac{\pi^2}{8} n^2 + n,
\e{HF}{(0)}(n) & = \ldots,
\\
\e{HF}{(1)}(n) & = \frac{\pi^2}{2}n^2 + \frac{4}{3} n,
\\
\e{HF}{(2)}(n) & = \frac{9\pi^2}{8}n^2 + \frac{23}{15} n.
\e{HF}{(1)}(n) & = \ldots.
\end{align}
\end{subequations}
All these states have the same (uniform) density $n = 2/(2\pi R)$ where $R$ is the radius of the ring on which the electrons are confined.
@ -298,35 +207,26 @@ The total energy of the ground and doubly-excited states are given by the two lo
\begin{equation}
\begin{split}
H_{ij}
& = \int_0^\pi \qty[ \frac{\psi_i(\omega)}{R} \frac{\psi_j(\omega)}{R} + \frac{\psi_i(\omega)\psi_j(\omega)}{2R\sin(\omega/2)} ] d\omega
& = \int_0^\pi \qty[ \frac{\psi_i(\omega)}{R} \frac{\psi_j(\omega)}{R} + \frac{\psi_i(\omega)\psi_j(\omega)}{2R\sin(\omega/2)} ] \sin^2 \omega \, d\omega
\\
& = \frac{\sqrt{\pi}}{2 R} \qty[ \frac{\Gamma\qty(\frac{i+j}{2})}{\Gamma\qty(\frac{i+j+1}{2})} + \frac{ij}{4R} \frac{\Gamma\qty(\frac{i+j-1}{2})}{\Gamma\qty(\frac{i+j+2}{2})} ],
& = \ldots,
\end{split}
\end{equation}
where $\omega = \theta_1 - \theta_2$ is the interelectronic angle, $\Gamma(x)$ is the Gamma function, \cite{NISTbook} and
where $\omega$ is the interelectronic angle, $\Gamma(x)$ is the Gamma function, \cite{NISTbook} and
\begin{equation}
\psi_i(\omega) = \sin(\omega/2) \sin^{i-1}(\omega/2), \quad i=1,\ldots,M
\psi_i(\omega) = \sin^{i-1}(\omega/2), \quad i=1,\ldots,M
\end{equation}
are (non-orthogonal) explicitly-correlated basis functions with overlap matrix elements
\begin{equation}
S_{ij}
= \int_0^\pi \psi_i(\omega)\psi_j(\omega) d\omega
= \sqrt{\pi} \frac{\Gamma\qty(\frac{i+j+1}{2})}{\Gamma\qty(\frac{i+j+2}{2})}.
= \int_0^\pi \psi_i(\omega)\psi_j(\omega) \sin^2 \omega \, d\omega
= \ldots.
\end{equation}
Thanks to this explicitly-correlated basis, the convergence rate of the energy is exponential with respect to $M$.
Therefore, high accuracy is reached with a very small number of basis functions.
Here, we typically use $M=10$.
For the singly-excited state, one has to modify the basis functions as
\begin{equation}
\psi_i(\omega) = \cos(\omega/2) \sin^{i-1} (\omega/2),
\end{equation}
and its energy is obtained by the lowest root of the Hamiltonian in this basis, and the matrix elements reads
\begin{align}
H_{ij} & = \frac{\sqrt{\pi}}{4 R} \qty[ \frac{\Gamma\qty(\frac{i+j}{2})}{\Gamma\qty(\frac{i+j+1}{3})} + \frac{3ij+i+j-1}{4R} \frac{\Gamma\qty(\frac{i+j-1}{2})}{\Gamma\qty(\frac{i+j+4}{2})} ],
\\
S_{ij} & = \frac{\sqrt{\pi}}{2} \frac{\Gamma\qty(\frac{i+j+1}{2})}{\Gamma\qty(\frac{i+j+4}{2})}.
\end{align}
The numerical values of the correlation energy for various $R$ are reported in Table \ref{tab:Ref} for the three states of interest.
The numerical values of the correlation energy for various $R$ are reported in Table \ref{tab:Ref} for the two states of interest.
%%% FIG 1 %%%
\begin{figure}
@ -372,6 +272,66 @@ The value of $c_1^{(I)}$ is obtained via the exact high-density expansion of the
Equation \eqref{eq:ec} is depicted in Fig.~\ref{fig:Ec} for each state alongside the data gathered in Table \ref{tab:Ref}.
%%% TABLE 1 %%%
\begin{table}
\caption{
\label{tab:OG_func}
Parameters of the correlation DFAs defined in Eq.~\eqref{eq:ec}.}
\begin{ruledtabular}
\begin{tabular}{lcddd}
State & $I$ & \tabc{$c_1^{(I)}$} & \tabc{$c_2^{(I)}$} & \tabc{$c_3^{(I)}$} \\
\hline
Ground state & $0$ & & & \\
Doubly-excited state & $1$ & & & \\
\end{tabular}
\end{ruledtabular}
\end{table}
%%% %%% %%% %%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{LDA-centered functional}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
In order to make the two-electron-based eDFA defined in Eq.~\eqref{eq:ecw} more universal and to ``center'' it on the jellium reference (as commonly done in DFT), we propose to \emph{shift} it as follows:
\begin{equation}
\label{eq:becw}
\be{xc}{\ew{}}(\n{}{}) = (1-\ew{}) \be{xc}{(0)}(\n{}{}) + \ew{} \be{c}{(1)}(\n{}{}),
\end{equation}
where
\begin{equation}
\be{xc}{(I)}(\n{}{}) = \e{xc}{(I)}(\n{}{}) + \e{xc}{\text{LDA}}(\n{}{}) - \e{xc}{(0)}(\n{}{}).
\end{equation}
The local-density approximation (LDA) exchange-correlation functional is
\begin{equation}
\e{xc}{\text{LDA}}(\n{}{}) = \e{x}{\text{LDA}}(\n{}{}) + \e{c}{\text{LDA}}(\n{}{}).
\end{equation}
Equation \eqref{eq:becw} can be recast
\begin{equation}
\label{eq:eLDA}
\be{xc}{\ew{}}(\n{}{})
= \e{xc}{\text{LDA}}(\n{}{}) + \ew{} \qty[\e{xc}{(1)}(\n{}{})-\e{xc}{(0)}(\n{}{})],
\end{equation}
which nicely highlights the centrality of the LDA in the present eDFA.
In particular, $\be{xc}{(0)}(\n{}{}) = \e{xc}{\text{LDA}}(\n{}{})$.
Consequently, in the following, we name this correlation functional ``eLDA'' as it is a natural extension of the LDA for ensembles.
This procedure can be theoretically justified by the generalized adiabatic connection formalism for ensembles (GACE) which was originally derived by Franck and Fromager. \cite{Franck_2014}
Within this in-principle-exact formalism, the (weight-dependent) correlation energy of the ensemble is constructed from the (weight-independent) ground-state functional (such as the LDA), yielding Eq.~\eqref{eq:eLDA}.
This is a crucial point as we intend to incorporate into standard functionals (which are ``universal'' in the sense that they do not depend on the number of electrons) information about excited states that will be extracted from finite systems (whose properties may depend on the number of electrons).
Finally, we note that, by construction,
\begin{equation}
\left. \pdv{\be{xc}{\ew{}}[\n{}{}]}{\ew{J}}\right|_{\n{}{} = \n{}{\ew{}}(\br)} = \be{xc}{(J)}[\n{}{\ew{}}(\br)] - \be{xc}{(0)}[\n{}{\ew{}}(\br)].
\end{equation}
\alert{As shown by Gould and Pittalis, comment on density- and and state-driven errors. \cite{Gould_2019}}
The density-functional approximations designed in this manuscript are based on highly-accurate energies for the ground state ($I=0$), the first singly-excited state ($I=1$), and the first doubly-excited state ($I=2$) of the (spin-polarized) two-electron glomium model.
We refer the interested reader to Refs.~\onlinecite{} for more details about this paradigm.
%%%%%%%%%%%%%%%%%%
%%% CONCLUSION %%%
%%%%%%%%%%%%%%%%%%