loos/BSE-PES
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity
bb0d93e42f
BSE-PES/BSE-PES.tex

542 lines
31 KiB
TeX
Raw Blame History

\documentclass[aps,prb,reprint,noshowkeys,superscriptaddress]{revtex4-1}
\usepackage{graphicx,dcolumn,bm,xcolor,microtype,multirow,amscd,amsmath,amssymb,amsfonts,physics,longtable,wrapfig,txfonts}
\usepackage[version=4]{mhchem}
\usepackage[utf8]{inputenc}
\usepackage[T1]{fontenc}
\usepackage{txfonts}
\usepackage[
colorlinks=true,
citecolor=blue,
breaklinks=true
]{hyperref}
\urlstyle{same}
\newcommand{\ie}{\textit{i.e}}
\newcommand{\eg}{\textit{e.g}}
\newcommand{\alert}[1]{\textcolor{red}{#1}}
\definecolor{darkgreen}{HTML}{009900}
\usepackage[normalem]{ulem}
\newcommand{\titou}[1]{\textcolor{red}{#1}}
\newcommand{\denis}[1]{\textcolor{purple}{#1}}
\newcommand{\xavier}[1]{\textcolor{darkgreen}{#1}}
\newcommand{\trashPFL}[1]{\textcolor{red}{\sout{#1}}}
\newcommand{\trashDJ}[1]{\textcolor{purple}{\sout{#1}}}
\newcommand{\trashXB}[1]{\textcolor{darkgreen}{\sout{#1}}}
\newcommand{\PFL}[1]{\titou{(\underline{\bf PFL}: #1)}}
\renewcommand{\DJ}[1]{\denis{(\underline{\bf DJ}: #1)}}
\newcommand{\XB}[1]{\xavier{(\underline{\bf XB}: #1)}}
\newcommand{\mc}{\multicolumn}
\newcommand{\fnm}{\footnotemark}
\newcommand{\fnt}{\footnotetext}
\newcommand{\tabc}[1]{\multicolumn{1}{c}{#1}}
\newcommand{\SI}{\textcolor{blue}{supporting information}}
\newcommand{\QP}{\textsc{quantum package}}
\newcommand{\T}[1]{#1^{\intercal}}
% coordinates
\newcommand{\br}[1]{\mathbf{r}_{#1}}
\newcommand{\dbr}[1]{d\br{#1}}
% methods
\newcommand{\GW}{$GW$}
\newcommand{\evGW}{ev\GW}
\newcommand{\qsGW}{qs\GW}
\newcommand{\GOWO}{$G_0W_0$}
\newcommand{\Hxc}{\text{Hxc}}
\newcommand{\xc}{\text{xc}}
\newcommand{\Ha}{\text{H}}
\newcommand{\co}{\text{x}}
%
\newcommand{\Norb}{N}
\newcommand{\Nocc}{O}
\newcommand{\Nvir}{V}
\newcommand{\IS}{\lambda}
% operators
\newcommand{\hH}{\Hat{H}}
% energies
\newcommand{\Enuc}{E^\text{nuc}}
\newcommand{\Ec}{E_\text{c}}
\newcommand{\EHF}{E^\text{HF}}
\newcommand{\EBSE}{E^\text{BSE}}
\newcommand{\EcRPA}{E_\text{c}^\text{dRPA}}
\newcommand{\EcBSE}{E_\text{c}^\text{BSE}}
\newcommand{\EcsBSE}{{}^1\EcBSE}
\newcommand{\EctBSE}{{}^3\EcBSE}
\newcommand{\IP}{\text{IP}}
\newcommand{\EA}{\text{EA}}
% orbital energies
\newcommand{\e}[1]{\epsilon_{#1}}
\newcommand{\eHF}[1]{\epsilon^\text{HF}_{#1}}
\newcommand{\eKS}[1]{\epsilon^\text{KS}_{#1}}
\newcommand{\eQP}[1]{\epsilon^\text{QP}_{#1}}
\newcommand{\eGOWO}[1]{\epsilon^\text{\GOWO}_{#1}}
\newcommand{\eGW}[1]{\epsilon^\text{\GW}_{#1}}
\newcommand{\eevGW}[1]{\epsilon^\text{\evGW}_{#1}}
\newcommand{\eGnWn}[2]{\epsilon^\text{\GnWn{#2}}_{#1}}
\newcommand{\Om}[2]{\Omega_{#1}^{#2}}
% Matrix elements
\newcommand{\A}[2]{A_{#1}^{#2}}
\newcommand{\tA}[2]{\Tilde{A}_{#1}^{#2}}
\newcommand{\B}[2]{B_{#1}^{#2}}
\renewcommand{\S}[1]{S_{#1}}
\newcommand{\ABSE}[2]{A_{#1}^{#2,\text{BSE}}}
\newcommand{\BBSE}[2]{B_{#1}^{#2,\text{BSE}}}
\newcommand{\ARPA}[2]{A_{#1}^{#2,\text{RPA}}}
\newcommand{\BRPA}[2]{B_{#1}^{#2,\text{RPA}}}
\newcommand{\G}[1]{G_{#1}}
\newcommand{\LBSE}[1]{L_{#1}}
\newcommand{\XiBSE}[1]{\Xi_{#1}}
\newcommand{\Po}[1]{P_{#1}}
\newcommand{\W}[2]{W_{#1}^{#2}}
\newcommand{\Wc}[1]{W^\text{c}_{#1}}
\newcommand{\vc}[1]{v_{#1}}
\newcommand{\Sig}[1]{\Sigma_{#1}}
\newcommand{\SigGW}[1]{\Sigma^\text{\GW}_{#1}}
\newcommand{\Z}[1]{Z_{#1}}
\newcommand{\MO}[1]{\phi_{#1}}
\newcommand{\ERI}[2]{(#1|#2)}
\newcommand{\sERI}[2]{[#1|#2]}
% excitation energies
\newcommand{\OmRPA}[2]{\Omega_{#1}^{#2,\text{RPA}}}
\newcommand{\OmRPAx}[2]{\Omega_{#1}^{#2,\text{RPAx}}}
\newcommand{\OmBSE}[2]{\Omega_{#1}^{#2,\text{BSE}}}
\newcommand{\spinup}{\downarrow}
\newcommand{\spindw}{\uparrow}
\newcommand{\singlet}{\uparrow\downarrow}
\newcommand{\triplet}{\uparrow\uparrow}
% Matrices
\newcommand{\bO}{\mathbf{0}}
\newcommand{\bI}{\mathbf{1}}
\newcommand{\bvc}{\mathbf{v}}
\newcommand{\bSig}{\mathbf{\Sigma}}
\newcommand{\bSigX}{\mathbf{\Sigma}^\text{x}}
\newcommand{\bSigC}{\mathbf{\Sigma}^\text{c}}
\newcommand{\bSigGW}{\mathbf{\Sigma}^\text{\GW}}
\newcommand{\be}{\mathbf{\epsilon}}
\newcommand{\beGW}{\mathbf{\epsilon}^\text{\GW}}
\newcommand{\beGnWn}[1]{\mathbf{\epsilon}^\text{\GnWn{#1}}}
\newcommand{\bde}{\mathbf{\Delta\epsilon}}
\newcommand{\bdeHF}{\mathbf{\Delta\epsilon}^\text{HF}}
\newcommand{\bdeGW}{\mathbf{\Delta\epsilon}^{GW}}
\newcommand{\bOm}[1]{\mathbf{\Omega}^{#1}}
\newcommand{\bA}[1]{\mathbf{A}^{#1}}
\newcommand{\btA}[1]{\Tilde{\mathbf{A}}^{#1}}
\newcommand{\bB}[1]{\mathbf{B}^{#1}}
\newcommand{\bX}[1]{\mathbf{X}^{#1}}
\newcommand{\bY}[1]{\mathbf{Y}^{#1}}
\newcommand{\bZ}[1]{\mathbf{Z}^{#1}}
\newcommand{\bK}{\mathbf{K}}
\newcommand{\bP}[1]{\mathbf{P}^{#1}}
% units
\newcommand{\IneV}[1]{#1 eV}
\newcommand{\InAU}[1]{#1 a.u.}
\newcommand{\InAA}[1]{#1 \AA}
\newcommand{\kcal}{kcal/mol}
\newcommand{\NEEL}{Universit\'e Grenoble Alpes, CNRS, Institut NEEL, F-38042 Grenoble, France}
\newcommand{\CEISAM}{Laboratoire CEISAM - UMR CNRS 6230, Universit\'e de Nantes, 2 Rue de la Houssini\`ere, BP 92208, 44322 Nantes Cedex 3, France}
\newcommand{\LCPQ}{Laboratoire de Chimie et Physique Quantiques (UMR 5626), Universit\'e de Toulouse, CNRS, UPS, France}
\begin{document}
\title{Ground-State Potential Energy Surfaces Within the Bethe-Salpeter Formalism: Pros and Cons}
\author{Xavier \surname{Blase}}
\email{xavier.blase@neel.cnrs.fr }
\affiliation{\NEEL}
\author{Denis \surname{Jacquemin}}
\email{denis.jacquemin@univ-nantes.fr}
\affiliation{\CEISAM}
\author{Pierre-Fran\c{c}ois \surname{Loos}}
\email{loos@irsamc.ups-tlse.fr}
\affiliation{\LCPQ}
\begin{abstract}
%\begin{wrapfigure}[12]{o}[-1.25cm]{0.4\linewidth}
% \centering
% \includegraphics[width=\linewidth]{TOC}
%\end{wrapfigure}
The many-body Green's function Bethe-Salpeter equation (BSE) formalism performed on top of a $GW$ calculation has shown to be a promising alternative to time-dependent density-functional theory (TD-DFT) in order to compute vertical transition energies for medium and large molecular systems.
Although less popular, the BSE formalism can also be employed to compute ground-state energies thanks to the adiabatic-connection fluctuation-dissipation theorem (ACFDT).
Here, we study the topological features of the ground-state potential energy surfaces (PES) of several diatomic molecules.
Our aim is to know whether or not the BSE formalism is able to reproduce faithfully the main features of these PES near equilibrium, and, in particular, the location of the minima on the ground-state potential energy surfaces.
Thanks to comparison with both similar and state-of-art computational approaches, we show that the present ACFDT-based approach is surprisingly accurate, and can even compete with coupled cluster methods.
However, we also observe, in some cases, unphysical irregularities on the ground-state PES, which are due to the appearance of a satellite resonance with a weight similar to that of the $GW$ quasiparticle peak.
\end{abstract}
\maketitle
%%%%%%%%%%%%%%%%%%%%%%%%
%\section{Introduction}
%\label{sec:intro}
%%%%%%%%%%%%%%%%%%%%%%%%
The features of ground- and excited-state potential energy surfaces (PES) are critical for the faithful description and a deeper understanding of photochemical and photophysical processes. \cite{Bernardi_1996,Olivucci_2010,Robb_2007}
For example, chemoluminescence, fluorescence and other related processes are associated with geometric relaxation of excited states, and structural changes upon electronic excitation. \cite{Navizet_2011}
Reliable predictions of these mechanisms which have attracted much experimental and theoretical interest lately require exploring the ground- and excited-state PES.
From a theoretical point of view, the accurate prediction of excited electronic states remains a challenge, \cite{Gonzales_2012, Loos_2020a} especially for large systems where state-of-the-art computational techniques (such as multiconfigurational methods \cite{Andersson_1990,Andersson_1992,Roos_1996,Angeli_2001}) cannot be afforded.
For such systems, one has to rely on more approximate, yet computationally cheaper approaches. \cite{Grimme_2004a,Ghosh_2018}
For the last two decades, time-dependent density-functional theory (TD-DFT) \cite{Casida} has been the go-to method to compute absorption and emission spectra in large molecular systems.
At a relatively low computational cost, TD-DFT can provide accurate vertical and adiabatic transition energies for low-lying excited states of organic molecules with a typical error of $0.2$--$0.4$ eV. \cite{Loos_2019b}
Similar to density-functional theory (DFT), \cite{Hohenberg_1964,Kohn_1965,ParrBook} TD-DFT is an in-principle exact theory which formal foundation relies on the Runge-Gross theorem. \cite{Runge_1984}
The Kohn-Sham (KS) formalism of TD-DFT transfers the complexity of the many-body problem to the exchange-correlation (xc) functional thanks to a judicious mapping between a time-dependent non-interacting reference system and its interacting analog which have both the exact same electronic density.
In TD-DFT, the excitation energies are obtained as the poles of the ground-state frequency-dependent linear response function, which is obtained by solving a Dyson equation in which the exchange and correlation effects are cast in the xc kernel.
This xc kernel plays for the excited states the same role as the xc functional for the ground state.
Hence, the PES for the excited states can be easily and efficiently obtained as a function of the molecular geometry by simply adding the ground-state DFT energy to the excitation energy of the selected state.
One of the strongest assets of TD-DFT is the availability of first- and second-order analytic nuclear gradients (\ie, the first derivatives of the excited-state energy with respect to the nuclear displacements), which enables the exploration of excited-state PES.\cite{Furche_2002}
Although the exact spectrum of the interacting system can be obtained in principle, both the ground-state xc potential and the xc kernel (which enters the definition of the response function) have to be approximated in practice.
This has major consequences, and it is well documented that, TD-DFT fails badly to model, for examples, charge transfer,\cite{Tozer_1999,Dreuw_2003,Sobolewski_2003,Dreuw_2004,Maitra_2017} Rydberg \cite{Tozer_1998,Tozer_2000,Casida_1998,Casida_2000,Tozer_2003} and doubly excited states. \cite{Levine_2006,Tozer_2000,Elliott_2011}
An important point from a practical point of view is that the TD-DFT xc kernel is usually considered as static instead of being frequency dependent. \cite{Romaniello_2009a,Sangalli_2011}
One key consequence of this so-called adiabatic approximation (based on the assumption that the density varies slowly with time) is that double excitations are completely absent from the TD-DFT spectra. \cite{Levine_2006,Tozer_2000,Huix-Rotllant_2010,Elliott_2011}
Yet another problem is the choice of the xc functionals as the quality of excitation energies are substantially dependent on this choice.
With a similar computational cost, the many-body Green's function Bethe-Salpeter equation (BSE) formalism \cite{Salpeter_1951,Strinati_1988} is a very valuable alternative to TD-DFT with early \textit{ab initio} calculations in condensed matter physics appearing at the end of the 90's. \cite{Albrecht_1998,Rohlfing_1998,Benedict_1998,vanderHorst_1999}
In the past few years, BSE has gained momentum for the study of molecular systems \cite{Ma_2009,Pushchnig_2002,Tiago_2003,Tiago_2008,Sai_2008,Palumno_2009,Rocca_2010,Sharifzadeh_2012,Cudazzo_2012,Boulanger_2014,Ljungberg_2015,Hirose_2015,Cocchi_2015,Ziaei_2017,Abramson_2017} and is now a serious candidate as a computationally inexpensive method that can effectively model excited states with a typical error of $0.1$--$0.3$ eV according to large and systematic benchmark calculations. \cite{Jacquemin_2015,Bruneval_2015,Blase_2016,Jacquemin_2016,Hung_2016,Hung_2017,Krause_2017,Jacquemin_2017,Blase_2018}
One of the main advantage of BSE compared to TD-DFT is that it allows a faithful description of charge-transfer states. \cite{Lastra_2011,Blase_2011b,Baumeier_2012,Duchemin_2012,Cudazzo_2013,Ziaei_2016}
Moreover, when performed on top of a (partially) self-consistently {\evGW} calculation, \cite{Hybertsen_1986, Shishkin_2007, Blase_2011, Faber_2011,Rangel_2016,Kaplan_2016,Gui_2018} BSE@{\evGW} has been shown to be weakly dependent on its starting point (\ie, on the xc functional selected for the underlying DFT calculation). \cite{Jacquemin_2016,Gui_2018}
However, similar to TD-DFT, the static version of BSE cannot describe multiple excitations. \cite{Romaniello_2009a,Sangalli_2011}
A significant limitation of the BSE formalism as compared to TD-DFT lies in the lack of analytic forces for the excited states, preventing efficient applications to the study of photoluminescence, photocatalysis, etc. While calculations of the {\GW} quasiparticle energies ionic gradients is becoming very popular,
\cite{Lazzeri_2008,Faber_2011b,Yin_2013,Faber_2015,Montserrat_2016,Zhenglu_2019} only one pioneering study of the excited-state BSE gradients has been published so far. \cite{Beigi_2003} In this study devoted to small molecules (\ce{CO} and \ce{NH3}), only the BSE excitation energy gradients were calculated, while computing the KS-DFT (LDA) forces as ground-state analogue.
Contrary to TD-DFT, the ground-state correlation energy calculated at the BSE level remains in its infancy with very few available studies for atomic and molecular systems. \cite{Olsen_2014,Holzer_2018,Li_2019,Li_2020}
As a matter of fact, in the largest recent available benchmark study of the 26 small molecules forming the HEAT test set, \cite{Holzer_2018} the BSE correlation energy, as evaluated within the adiabatic-connection fluctuation-dissipation (ACFDT) framework, was mostly discarded from the set of tested techniques due to instabilities (negative frequency modes in the BSE polarization propagator) and replaced by an approximate (RPAsX) approach where the screened-Coulomb potential matrix elements was removed from the resonant electron-hole contribution. \cite{Maggio_2016,Holzer_2018}
Such a modified BSE polarization propagator was inspired by a previous study on the homogeneous electron gas. \cite{Maggio_2016}
With such an approximation, amounting to neglect excitonic effects in the electron-hole propagator, the question of using either KS-DFT or {\GW} eigenvalues in the construction of the propagator becomes further relevant, increasing accordingly the number of possible definitions for the ground-state correlation energy.
Finally, renormalizing or not the Coulomb interaction by the coupling parameter $\lambda$ in the Dyson equation for the interacting polarizability leads to two different versions of the BSE correlation energy.
Here, in analogy to the random-phase approximation (RPA)-type formalismes, \cite{Furche_2008, Angyan_2011, Holzer_2018} the ground-state BSE energy is calculated in the adiabatic-connection fluctuation-dissipation theorem framework.
Embracing this definition, the purpose of the present study is to investigate the quality of ground--state PES near equilibrium obtained within the BSE approach.
The location of the minima on the ground-state PES is of particular interest.
This study is a first preliminary step towards the development of analytical nuclear gradients within the BSE@{\GW} formalism.
%The paper is organized as follows.
%In Sec.~\ref{sec:theo}, we introduce the equations behind the BSE formalism.
%In particular, the general equations are reported in Subsec.~\ref{sec:BSE}, and the corresponding equations obtained in a finite basis in Subsec.~\ref{sec:BSE_basis}.
%Subsection \ref{sec:BSE_energy} defines the BSE total energy, and various other quantities of interest for this study.
%Computational details needed to reproduce the results of the present work are reported in Sec.~\ref{sec:comp_details}.
%Section \ref{sec:PES} reports ground-state PES for various diatomic molecules.
%Finally, we draw our conclusion in Sec.~\ref{sec:conclusion}.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\section{Theory}
%\label{sec:theo}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\subsection{The Bethe-Salpeter equation}
%\label{sec:BSE}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
In order to compute the neutral (optical) excitations of the system and their associated oscillator strengths, the BSE expresses the two-body propagator \cite{Strinati_1988, Martin_2016}
\begin{multline}
\label{eq:BSE}
\LBSE{}(1,2,1',2') = \LBSE{0}(1,2,1',2')
\\
+ \int d3 d4 d5 d6 \LBSE{0}(1,4,1',3) \XiBSE{}(3,5,4,6) \LBSE{}(6,2,5,2')
\end{multline}
as the linear response of the one-body Green's function $\G{}$ with respect to a general non-local external potential
\begin{equation}
\XiBSE{}(3,5,4,6) = i \fdv{[\vc{\Ha}(3) \delta(3,4) + \Sig{\xc}(3,4)]}{\G{}(6,5)},
\end{equation}
which takes into account the self-consistent variation of the Hartree potential
\begin{equation}
\vc{\Ha}(1) = - i \int d2 \vc{}(2) \G{}(2,2^+),
\end{equation}
(where $\vc{}$ is the bare Coulomb operator) and the exchange-correlation self-energy $\Sig{\xc}$.
In Eq.~\eqref{eq:BSE}, $\LBSE{0}(1,2,1',2') = - i \G{}(1,2')\G{}(2,1')$, and $(1)=(\br{1},\sigma_1,t_1)$ is a composite index gathering space, spin and time variables.
In the {\GW} approximation, \cite{Hedin_1965,Aryasetiawan_1998,Onida_2002,Martin_2016,Reining_2017} we have
\begin{equation}
\SigGW{\xc}(1,2) = i \G{}(1,2) \W{}{}(1^+,2),
\end{equation}
where $\W{}{}$ is the screened Coulomb operator, and hence the BSE reduces to
\begin{equation}
\XiBSE{}(3,5,4,6) = \delta(3,4) \delta(5,6) \vc{}(3,6) - \delta(3,6) \delta(4,5) \W{}{}(3,4),
\end{equation}
where, as commonly done, we have neglected the term $\delta \W{}{}/\delta \G{}$ in the functional derivative of the self-energy. \cite{Hanke_1980,Strinati_1984,Strinati_1982}
Finally, the static approximation is enforced, \ie, $\W{}{}(1,2) = \W{}{}(\{\br{1}, \sigma_1, t_1\},\{\br{2}, \sigma_2, t_2\}) \delta(t_1-t_2)$, which corresponds to restricting $\W{}{}$ to its static limit, \ie, $\W{}{}(1,2) = \W{}{}(\{\br{1},\sigma_1\},\{\br{2},\sigma_2\}; \omega=0)$.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\subsection{BSE in a finite basis}
%\label{sec:BSE_basis}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
For a closed-shell system, in order to compute the singlet BSE excitation energies (within the static approximation) of the physical system (\ie, $\IS = 1$) in a finite basis, one must solve the following linear response problem \cite{Casida,Dreuw_2005,Martin_2016}
\begin{equation}
\label{eq:LR}
\begin{pmatrix}
\bA{\IS} & \bB{\IS} \\
-\bB{\IS} & -\bA{\IS} \\
\end{pmatrix}
\begin{pmatrix}
\bX{\IS}_m \\
\bY{\IS}_m \\
\end{pmatrix}
=
\Om{m}{\IS}
\begin{pmatrix}
\bX{\IS}_m \\
\bY{\IS}_m \\
\end{pmatrix},
\end{equation}
where $\Om{m}{\IS}$ is the $m$th excitation energy with eigenvector $\T{(\bX{\IS}_m \, \bY{\IS}_m)}$ at interpolation strength $\IS$, and we have assumed real-valued spatial orbitals $\{\MO{p}(\br{})\}_{1 \le p \le \Norb}$.
The matrices $\bA{\IS}$, $\bB{\IS}$, $\bX{\IS}$, and $\bY{\IS}$ are all of size $\Nocc \Nvir \times \Nocc \Nvir$ where $\Nocc$ and $\Nvir$ are the number of occupied and virtual orbitals (\ie, $\Norb = \Nocc + \Nvir$), respectively.
In the following, the index $m$ labels the $\Nocc \Nvir$ single excitations, $i$ and $j$ are occupied orbitals, $a$ and $b$ are unoccupied orbitals, while $p$, $q$, $r$, and $s$ indicate arbitrary orbitals.
In the absence of instabilities (\ie, $\bA{\IS} - \bB{\IS}$ is positive-definite), \cite{Dreuw_2005} Eq.~\eqref{eq:LR} is usually transformed into an Hermitian eigenvalue problem of smaller dimension
\begin{equation}
\label{eq:small-LR}
(\bA{\IS} - \bB{\IS})^{1/2} (\bA{\IS} + \bB{\IS}) (\bA{\IS} - \bB{\IS})^{1/2} \bZ{\IS} = (\bOm{\IS})^2 \bZ{\IS},
\end{equation}
where the excitation amplitudes are
\begin{subequations}
\begin{align}
\bX{\IS} + \bY{\IS} = (\bOm{\IS})^{-1/2} (\bA{\IS} - \bB{\IS})^{1/2} \bZ{\IS},
\\
\bX{\IS} - \bY{\IS} = (\bOm{\IS})^{1/2} (\bA{\IS} - \bB{\IS})^{-1/2} \bZ{\IS}.
\end{align}
\end{subequations}
In the case of BSE, the specific expression of the matrix elements are
\begin{subequations}
\begin{align}
\label{eq:LR_BSE}
\ABSE{ia,jb}{\IS} & = \delta_{ij} \delta_{ab} (\eGW{a} - \eGW{i}) + \IS \W{ia,bj}{\IS}(\omega = 0) - \lambda \ERI{ia}{jb},
\\
\BBSE{ia,jb}{\IS} & = \W{ia,jb}{\IS}(\omega = 0) - \lambda \ERI{ib}{ja} ,
\end{align}
\end{subequations}
where $\eGW{p}$ are the {\GW} quasiparticle energies,
\begin{multline}
\label{eq:W}
\W{ia,jb}{\IS}(\omega) = 2 \ERI{ia}{jb}
\\
+ \sum_m^{\Nocc \Nvir} \sERI{ia}{m} \sERI{jb}{m} \qty(\frac{1}{\omega - \OmRPA{m}{\IS} + i \eta} + \frac{1}{\omega + \OmRPA{m}{\IS} - i \eta})
\end{multline}
are the elements of the screened Coulomb operator $\W{}{\IS}$,
\begin{equation}
\sERI{pq}{m} = \sum_i^{\Nocc} \sum_a^{\Nvir} \ERI{pq}{ia} (\bX{\IS}_m + \bY{\IS}_m)_{ia}
\end{equation}
are the screened two-electron integrals,
\begin{equation}
\ERI{pq}{rs} = \iint \frac{\MO{p}(\br{}) \MO{q}(\br{}) \MO{r}(\br{}') \MO{s}(\br{}')}{\abs*{\br{} - \br{}'}} \dbr{} \dbr{}',
\end{equation}
are the bare two-electron integrals, and $\delta_{pq}$ is the Kronecker delta. \cite{NISTbook}
In Eq.~\eqref{eq:W}, $\eta$ is a positive infinitesimal, and $\OmRPA{m}{\IS}$ are the direct (\ie, without exchange) RPA neutral excitation energies computed by solving the linear eigenvalue problem \eqref{eq:LR} with the following matrix elements
\begin{subequations}
\begin{align}
\label{eq:LR_BSE}
\ARPA{ia,jb}{\IS} & = \delta_{ij} \delta_{ab} (\eHF{a} - \eHF{i}) + \IS \ERI{ia}{bj},
\\
\BRPA{ia,jb}{\IS} & = \IS \ERI{ia}{jb},
\end{align}
\end{subequations}
where $\eHF{p}$ are the HF orbital energies.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\subsection{Ground-state BSE energy}
%\label{sec:BSE_energy}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
The key quantity to define in the present context is the total ground-state BSE energy $\EBSE$.
Although this choice is not unique, \cite{Holzer_2018} we propose to define it as
\begin{equation}
\label{eq:EtotBSE}
\EBSE = \Enuc + \EHF + \EcBSE
\end{equation}
where $\Enuc$ and $\EHF$ are the nuclear repulsion energy and ground-state HF energy (respectively), and
\begin{equation}
\label{eq:EcBSE}
\EcBSE = \frac{1}{2} \int_0^1 \Tr(\bK \bP{\IS}) d\IS
\end{equation}
is the ground-state BSE correlation energy computed in the adiabatic connection framework, where
\begin{equation}
\bK =
\begin{pmatrix}
\btA{\IS=1} & \bB{\IS=1} \\
\bB{\IS=1} & \btA{\IS=1} \\
\end{pmatrix}
\end{equation}
is the interaction kernel \cite{Angyan_2011, Holzer_2018} [with $\tA{ia,jb}{\IS} = \IS \ERI{ia}{bj}$],
\begin{equation}
\bP{\IS} =
\begin{pmatrix}
\bY{\IS} \T{(\bY{\IS})} & \bY{\IS} \T{(\bX{\IS})} \\
\bX{\IS} \T{(\bY{\IS})} & \bX{\IS} \T{(\bX{\IS})} \\
\end{pmatrix}
-
\begin{pmatrix}
\bO & \bO \\
\bO & \bI \\
\end{pmatrix}
\end{equation}
is the correlation part of the two-electron density matrix at interaction strength $\IS$, $\Tr$ denotes the matrix trace and $\T{}$ the matrix transpose.
Note that the present definition of the BSE correlation energy [see Eq.~\eqref{eq:EcBSE}] has be labeled ``XBS'' for ``extended Bethe Salpeter'' by Holzer et al. \cite{Holzer_2018}
Several important comments are in order here.
For spin-restricted closed-shell molecular systems around their equilibrium geometry (such as the ones studied here), it is rare to encounter singlet instabilities.
However, they may appear in the presence of strong correlation (\eg, when the bond is stretch).
In such a case, this hampers the use of Eq.~\eqref{eq:EcBSE}.
Triplet instabilities are much more common.
However, contrary to the plasmon formulation (an alternative expression of the BSE correlation energy), \cite{Schuck_Book, Gell-Mann_1957, Rowe_1968, Sawada_1957b, Li_2020} the triplet excitations do not contribute in the adiabatic connection formulation, which is an indisputable advantage of the ACFDT appraoch.
Although at the RPA level, the plasmon and adiabatic connection formulations are equivalent, \cite{Sawada_1957b, Fukuta_1964, Furche_2008} this is not the case at the BSE level.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\section{Computational details}
%\label{sec:comp_details}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
All the preliminary {\GW} calculations performed to obtain the screened Coulomb operator and the quasiparticle energies have been done using a (restricted) Hartree-Fock (HF) starting point, which is a very adequate choice in the case of the (small) systems that we have considered here.
Perturbative {\GW} (or {\GOWO}) \cite{Hybertsen_1985a, Hybertsen_1986} calculations are employed as starting point to compute the BSE neutral excitations.
These will be labeled as BSE@{\GOWO}.
In the case of {\GOWO}, the quasiparticle energies have been obtained by linearizing the non-linear, frequency-dependent quasiparticle equation.
Further details about our implementation of {\GOWO} can be found in Refs.~\onlinecite{Loos_2018,Veril_2018}.
Finally, the infinitesimal $\eta$ has been set to zero for all calculations.
The numerical integration required to compute the correlation energy along the adiabatic path has been performed with a 21-point Gauss-Legendre quadrature.
As one-electron basis sets, we employ the Dunning family (cc-pVXZ) defined with cartesian gaussian functions.
Unless, otherwise stated, the frozen-core approximation has not been enforced in our calculations in order to provide a fairer comparison between methods.
However, we have found that the conclusions drawn in the present study hold within the frozen-core approximation (see {\SI}).
Because Eq.~\eqref{eq:EcBSE} requires the entire BSE excitation spectrum (both singlet and triplet), we perform a complete diagonalization of the $\Nocc \Nvir \times \Nocc \Nvir$ BSE linear response matrix [see Eq.~\eqref{eq:small-LR}], which corresponds to a $\order{\Nocc^3 \Nvir^3} = \order{\Norb^6}$ computational cost.
This step is, by far, the computational bottleneck in our current implementation.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\section{Potential energy surfaces}
%\label{sec:PES}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% TABLE I %%%
\begin{table*}
\caption{
Equilibrium distances (in bohr) of the ground state of diatomic molecules obtained at various levels of theory and basis sets.}
\label{tab:Req}
\begin{ruledtabular}
\begin{tabular}{llcccccccc}
& & \mc{8}{c}{Molecules} \\
\cline{3-10}
Method & Basis & \ce{H2} & \ce{LiH}& \ce{LiF}& \ce{N2} & \ce{CO} & \ce{BF} & \ce{F2} & \ce{HCl}\\
\hline
CC3 & cc-pVDZ & 1.438 & 3.043 & 3.012 & 2.114 & 2.166 & 2.444 & 2.740 & 2.435 \\
& cc-pVTZ & 1.403 & 3.011 & 2.961 & 2.079 & 2.143 & 2.392 & 2.669 & 2.413 \\
& cc-pVQZ & 1.402 & 3.019 & 2.963 & 2.075 & 2.136 & 2.390 & 2.663 & 2.403 \\
CCSD & cc-pVDZ & 1.438 & 3.044 & 3.006 & 2.101 & 2.149 & 2.435 & 2.695 & 2.433 \\
& cc-pVTZ & 1.403 & 3.012 & 2.954 & 2.064 & 2.126 & 2.382 & 2.629 & 2.409 \\
& cc-pVQZ & 1.402 & 3.020 & 2.953 & 2.059 & 2.118 & 2.118 & 2.621 & 2.398 \\
CC2 & cc-pVDZ & 1.426 & 3.046 & 3.026 & 2.146 & 2.187 & 2.444 & 2.710 & 2.427 \\
& cc-pVTZ & 1.393 & 3.008 & 2.995 & 2.109 & 2.163 & 2.394 & 2.664 & 2.406 \\
& cc-pVQZ & 1.391 & 2.989 & 2.982 & 2.106 & 2.156 & 2.393 & 2.665 & 2.396 \\
MP2 & cc-pVDZ & 1.426 & 3.041 & 3.010 & 2.133 & 2.166 & 2.431 & 2.681 & 2.426 \\
& cc-pVTZ & 1.393 & 3.004 & 2.968 & 2.095 & 2.144 & 2.383 & 2.636 & 2.405 \\
& cc-pVQZ & 1.391 & 3.008 & 2.970 & 2.091 & 2.137 & 2.382 & 2.634 & 2.395 \\
BSE@{\GOWO}@HF & cc-pVDZ & 1.437 & 3.042 & 3.000 & 2.107 & 2.153 & 2.407 & 2.700 & >2.440 \\
& cc-pVTZ & 1.404 & 3.023 & glitch & & & <2.420 & & <2.410 \\
& cc-pVQZ & 1.399 & & & & & & & \\
RPA@{\GOWO}@HF & cc-pVDZ & 1.426 & 3.019 & 2.994 & 2.083 & 2.144 & 2.403 & 2.691 & 2.436 \\
& cc-pVTZ & 1.388 & 3.013 & glitch & & & <2.420 & & <2.410 \\
& cc-pVQZ & 1.382 & & & & & & & \\
RPAx@HF & cc-pVDZ & 1.428 & 3.040 & 2.998 & 2.077 & 2.130 & 2.417 & NaN & 2.424 \\
& cc-pVTZ & 1.395 & 3.003 & <2.990 & & & <2.420 & & <2.410 \\
& cc-pVQZ & 1.394 & & & & & & & \\
RPA@HF & cc-pVDZ & 1.431 & 3.021 & 2.999 & 2.083 & 2.134 & & 2.623 & 2.424 \\
& cc-pVTZ & 1.388 & 2.978 & <2.990 & & & 2.416 & & <2.410 \\
& cc-pVQZ & 1.386 & & & & & <2.420 & & \\
% FROZEN CORE VERSION
% Method & Basis & \ce{H2} & \ce{LiH}& \ce{LiF}& \ce{N2} & \ce{CO} & \ce{BF} & \ce{F2} & \ce{HCl}\\
% \hline
% CC3 & cc-pVDZ & 1.438 & 3.052 & 3.014 & 2.115 & 2.167 & 2.447 & 2.741 & 2.438 \\
% & cc-pVTZ & 1.403 & 3.036 & 2.985 & 2.087 & 2.150 & 2.405 & 2.672 & 2.414 \\
% & cc-pVQZ & 1.402 & 3.037 & 2.985 & 2.080 & 2.142 & 2.398 & 2.667 & 2.413 \\
% CCSD & cc-pVDZ & 1.438 & 3.044 & 3.006 & 2.101 & 2.149 & 2.435 & 2.695 & 2.433 \\
% & cc-pVTZ & 1.403 & 3.012 & 2.954 & 2.064 & 2.126 & 2.382 & 2.629 & 2.409 \\
% & cc-pVQZ & 1.402 & 3.020 & 2.953 & 2.059 & 2.118 & 2.380 & 2.621 & 2.398 \\
% CC2 & cc-pVDZ & 1.426 & & & & & & & \\
% & cc-pVTZ & 1.393 & & & & & & & \\
% & cc-pVQZ & 1.391 & & & & & & & \\
% MP2 & cc-pVDZ & 1.426 & 3.049 & 3.012 & 2.134 & 2.167 & 2.433 & 2.681 & 2.429 \\
% & cc-pVTZ & 1.393 & 3.026 & 2.990 & 2.104 & 2.151 & 2.395 & 2.640 & 2.407 \\
% & cc-pVQZ & 1.391 & 3.026 & 2.990 & 2.098 & 2.144 & 2.389 & 2.638 & 2.405 \\
% BSE@{\GOWO}@HF & cc-pVDZ & 1.437 & & & & & & & \\
% & cc-pVTZ & 1.404 & & & & & & & \\
% & cc-pVQZ & 1.399 & & & & & & & \\
% RPA@{\GOWO}@HF & cc-pVDZ & 1.426 & & & & & & & \\
% & cc-pVTZ & 1.388 & & & & & & & \\
% & cc-pVQZ & 1.382 & & & & & & & \\
% RPAx@HF & cc-pVDZ & 1.428 & & & & & & & \\
% & cc-pVTZ & 1.395 & & & & & & & \\
% & cc-pVQZ & 1.394 & & & & & & & \\
% RPA@HF & cc-pVDZ & 1.431 & & & & & & & \\
% & cc-pVTZ & 1.388 & & & & & & & \\
% & cc-pVQZ & 1.386 & & & & & & & \\
\end{tabular}
\end{ruledtabular}
\end{table*}
%%% FIG 1 %%%
\begin{figure*}
\includegraphics[width=0.45\linewidth]{H2_GS_VTZ}
\includegraphics[width=0.45\linewidth]{LiH_GS_VTZ}
\includegraphics[width=0.45\linewidth]{LiF_GS_VTZ}
\includegraphics[width=0.45\linewidth]{N2_GS_VTZ}
\includegraphics[width=0.45\linewidth]{CO_GS_VTZ}
\includegraphics[width=0.45\linewidth]{BF_GS_VTZ}
\includegraphics[width=0.45\linewidth]{F2_GS_VTZ}
\includegraphics[width=0.45\linewidth]{HCl_GS_VTZ}
\caption{
PES of the ground state of diatomic molecules around their respective equilibrium geometry obtained at various levels of theory with the cc-pVTZ basis set.
Additional graphs for other basis sets and within the frozen-core approximation can be found in the {\SI}.
\label{fig:PES}
}
\end{figure*}
%%% %%% %%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\subsection{Hydrogen molecule}
%\label{sec:H2}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\subsection{Lithium hydride and lithium fluoride}
%\label{sec:LiX}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%
%\subsection{The isoelectronic sequence: \ce{N2}, \ce{CO}, and \ce{BF}}
%\label{sec:isoN2}
%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%
%\section{Conclusion}
%\label{sec:conclusion}
%%%%%%%%%%%%%%%%%%%%%%%%
\titou{A nice conclusion here saying that what we have done is pretty awesome.}
%%%%%%%%%%%%%%%%%%%%%%%%
\section*{Supporting Information}
%%%%%%%%%%%%%%%%%%%%%%%%
See {\SI} for \titou{bla bla bla.}
%%%%%%%%%%%%%%%%%%%%%%%%
\begin{acknowledgements}
This work was performed using HPC resources from GENCI-TGCC (Grant No.~2018-A0040801738) and CALMIP (Toulouse) under allocation 2019-18005.
Funding from the \textit{``Centre National de la Recherche Scientifique''} is acknowledged.
This work has been supported through the EUR grant NanoX ANR-17-EURE-0009 in the framework of the \textit{``Programme des Investissements d'Avenir''}.
\end{acknowledgements}
%%%%%%%%%%%%%%%%%%%%%%%%
\bibliography{BSE-PES,BSE-PES-control}
\end{document}
Reference in New Issue View Git Blame Copy Permalink