PTEROSOR/BSEdyn
10
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

git clean up

Browse Source
This commit is contained in:
Pierre-Francois Loos 2020-05-17 22:39:03 +02:00
parent 8df8f336cb
commit 2acfdcd32f
4 changed files with 13995 additions and 132 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

13135
BSEdyn.bib Normal file
View File

File diff suppressed because it is too large Load Diff

439
BSEdyn.tex Normal file
View File

@ -0,0 +1,439 @@
\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{natbib}
%\usepackage[extra]{tipa}
%\bibliographystyle{achemso}
%\AtBeginDocument{\nocite{achemso-control}}
\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{\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{RPA}}
\newcommand{\EcRPAx}{E_\text{c}^\text{RPAx}}
\newcommand{\EcBSE}{E_\text{c}^\text{BSE}}
\newcommand{\IP}{\text{IP}}
\newcommand{\EA}{\text{EA}}
\newcommand{\Req}{R_\text{eq}}
% 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^{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{\ARPAx}[2]{A_{#1}^{#2,\text{RPAx}}}
\newcommand{\BRPAx}[2]{B_{#1}^{#2,\text{RPAx}}}
\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^{GW}_{#1}}
\newcommand{\Z}[1]{Z_{#1}}
\newcommand{\MO}[1]{\phi_{#1}}
\newcommand{\ERI}[2]{(#1|#2)}
\newcommand{\sERI}[3]{[#1|#2]^{#3}}
%% bold in Table
\newcommand{\bb}[1]{\textbf{#1}}
\newcommand{\rb}[1]{\textbf{\textcolor{red}{#1}}}
\newcommand{\gb}[1]{\textbf{\textcolor{darkgreen}{#1}}}
% 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}^{GW}}
\newcommand{\be}{\mathbf{\epsilon}}
\newcommand{\beGW}{\mathbf{\epsilon}^{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}[2]{\mathbf{Z}_{#1}^{#2}}
\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}
\DeclareMathOperator*{\argmax}{argmax}
\DeclareMathOperator*{\argmin}{argmin}
\newcommand\vari{{\varepsilon}_i}
\newcommand\vara{{\varepsilon}_a}
\newcommand\varj{{\varepsilon}_j}
\newcommand\varb{{\varepsilon}_b}
\newcommand\varn{{\varepsilon}_n}
\newcommand\varm{{\varepsilon}_m}
\newcommand\Oms{{\Omega}_s}
\newcommand\hOms{\frac{{\Omega}_s}{2}}
\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}
\newcommand{\CEA}{Universit\'e Grenoble Alpes, CEA, IRIG-MEM-L Sim, 38054 Grenoble, France}
\begin{document}
\title{Pros and Cons of the Bethe-Salpeter Formalism for Ground-State Energies}
\author{Pierre-Fran\c{c}ois \surname{Loos}}
\email{loos@irsamc.ups-tlse.fr}
\affiliation{\LCPQ}
\author{Xavier \surname{Blase}}
\email{xavier.blase@neel.cnrs.fr }
\affiliation{\NEEL}
\begin{abstract}
This is the abstract
%\\
%\bigskip
%\begin{center}
% \boxed{\includegraphics[width=0.5\linewidth]{TOC}}
%\end{center}
%\bigskip
\end{abstract}
\maketitle
%%%%%%%%%%%%%%%%%%%%%%%%
\section{Introduction}
\label{sec:intro}
%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%
\section{Theory}
\label{sec:theory}
%%%%%%%%%%%%%%%%%%%%%%%%
The Fourier components with respect to time $t_1$ of $iL_0(1, 4; 1', 3) = G(1, 3)G(4, 1')$ reads, dropping the (space/spin)-variables:
\begin{align*}
[iL_0]( \omega_1 ) = \frac{ 1 }{ 2\pi } \int d \omega \; G(\omega - \frac{\omega_1}{2} ) G( {\omega} + \frac{\omega_1}{2} ) e^{ i \omega \tau_{34} } e^{ i \omega_1 t^{34} }
\end{align*}
with $\tau_{34} = t_3 - t_4$ and
$t^{34} = (t_3 + t_4)/2$. Plugging now the 1-body Green's function Lehman representation, e.g.
$$
G(x_1,x_3 ; \omega) = \sum_n \frac{ \phi_n(x_1) \phi_n^*(x_3) } { \omega - \varepsilon_n + i \eta \text{sgn}(\varepsilon_n - \mu) }
$$
and projecting on $\phi_a^*(x_1) \phi_i(x_{1'})$, one obtains the $\omega_1= \Oms$ component
\begin{align*}
\int dx_1 dx_{1'} \; & \phi_a^*(x_1) \phi_i(x_{1'}) L_0(x_1,3;x_{1'},4; \Oms) = e^{i \Oms t^{34} } \times \\
& \frac{ \phi_a^*(x_3) \phi_i(x_4) } { \Oms - ( \vara - \vari ) + i \eta }
\Big( \theta( \tau ) e^{i ( \vari + \hOms) \tau }
+ \theta( - \tau ) e^{i (\vara - \hOms \tau } \Big)
\end{align*}
with $\tau = \tau_{34}$.
We further obtain the spectral representation of
$\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) | N,s \rangle$
expanding the field operators over a complete orbital basis creation/destruction operators:
\begin{align*}
\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) & | N,s \rangle = - \Big( e^{ -i \Omega_s t^{34} } \Big) \sum_{mn} \phi_m(x_3) \phi_n^*(x_4) \langle N | {\hat a}_n^{\dagger} {\hat a}_m | N,s \rangle \times \nonumber \\
\times & \Big( \theta( \tau ) e^{- i ( \varepsilon_m - \hOms ) \tau }
+ \theta( -\tau ) e^{ - i ( \varepsilon_n + \hOms) \tau } \Big)
\end{align*}
with $\tau = \tau_{34}$ and where the $ \lbrace \varepsilon_{n/m} \rbrace$ are proper addition/removal energies such that e.g.
$$
e^{ i H \tau } {\hat a}_m^{\dagger} | N \rangle = e^{ i (E_0^N + \varepsilon_m ) \tau } {\hat a} _m^{\dagger} | N \rangle
$$
Selecting (n,m)=(j,b) yields the largest components
$A_{jb}^{s} = \langle N | {\hat a}_j^{\dagger} {\hat a}_b | N,s \rangle $, while (n,m)=(b,j) yields much weaker
$B_{jb}^{s} = \langle N | {\hat a}_b^{\dagger} {\hat a}_j | N,s \rangle $ contributions. We used chemist notations with (i,j) indexing occupied orbitals and (a,b) virtual ones. Neglecting the $B_{jb}^{s}$ leads to the Tamm Dancoff approximation (TDA). Obtaining similarly the spectral representation of $ \langle N | T {\hat \psi}(1) {\hat \psi}^{\dagger}(1') | N,s \rangle$ ($t_{1'} = t_1^{+}$) projected onto $\phi_a^*(x_1) \phi_i(x_{1'})$,
one obtains after a few tedious manipulations (see Supplemental Information) the dynamical Bethe-Salpeter equation (DBSE) :
\begin{align}
( \varepsilon_a - \varepsilon_i - \Omega_s ) A_{ia}^{s}
&+ \sum_{jb} \Big( v_{ai,bj} - \widetilde{W}_{ij,ab}(\Oms) \Big) A_{jb}^{s} \\
&+ \sum_{bj} \Big( v_{ai,jb} - \widetilde{W}_{ib,aj}(\Oms) \Big) B_{jb}^{s}
= 0
\end{align}
with an effective dynamically screened Coulomb potential (see Pina eq. 24):
\begin{align}
\widetilde{W}_{ij,ab}(\Oms) &= { i \over 2 \pi} \int d\omega \; e^{-i \omega 0^+ } W_{ij,ab}(\omega) \times \\
\hskip 1cm &\times \left[ \frac{1}{ (\Oms - \omega) - ( \varb - \vari ) +i \eta } + \frac{1}{ (\Oms + \omega) - ( \vara - \varj ) + i\eta } \right] \nonumber
\end{align}
In the present study, we use the exact spectral representation of $W(\omega)$ at the RPA level:
\begin{align*}
W_{ij,ab}(\omega) &= (ij|ab) + 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
& \times \Big( \frac{1}{ \omega-\Omega_m^{RPA} + i\eta } - \frac{1}{ \omega + \Omega_m^{RPA} - i\eta } \Big)
\end{align*}
so that
\begin{align}
\widetilde{W}_{ij,ab}( \Oms ) &= (ij|ab) + 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
& \times \left[ \frac{ 1 }{ \Omega_{ib}^{s} - \Omega_m^{RPA} + i\eta } + \frac{ 1}{ \Omega_{ja}^{s} - \Omega_m^{RPA} + i\eta }
\right] \nonumber
\end{align}
with e.g. $ \Omega_{ib}^{s} = \Oms - ( \varepsilon_b - \varepsilon_i) $. \textcolor{red}{Due to excitonic effects, the lowest BSE ${\Omega}_1$ excitation energy stands lower than the lowest $\Omega_m^{RPA}$ excitation energy, so that
e.g. $( \Omega_{ib}^{s} - \Omega_m^{RPA} )$ is strictly negative and cannot diverge. Further, $\Omega_{ib}^{s}$ and $\Omega_{ja}^{s}$ are necessarily negative for in-gap low lying BSE excitations, such that
$$
\left[ \frac{ 1 }{ \Omega_{ib}^{s} - \Omega_m^{RPA} + i\eta } + \frac{ 1}{ \Omega_{ja}^{s} - \Omega_m^{RPA} + i\eta }
\right]
<
\Big( \frac{1}{ \omega-\Omega_m^{RPA} + i\eta } - \frac{1}{ \omega + \Omega_m^{RPA} - i\eta } \Big) < 0
$$
in the limit $(\omega \rightarrow 0)$ of the standard adiabatic BSE . WELL, do we know the sign of
$[ij|m] [ab|m]$ ?? }
\titou{This is the theory section from the previous paper.}
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)$.
For a closed-shell system in a finite basis, to compute the singlet BSE excitation energies (within the static approximation) of the physical system (\ie, $\IS = 1$), 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 interaction strength $\IS$, $\T{}$ is the matrix transpose, and we assume 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$ is the total number of spatial orbitals), 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, when $\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{m}{\IS} = (\Om{m}{\IS})^2 \bZ{m}{\IS},
\end{equation}
where the excitation amplitudes are
\begin{subequations}
\begin{align}
(\bX{\IS} + \bY{\IS})_m = (\Om{m}{\IS})^{-1/2} (\bA{\IS} - \bB{\IS})^{+1/2} \bZ{m}{\IS},
\\
(\bX{\IS} - \bY{\IS})_m = (\Om{m}{\IS})^{+1/2} (\bA{\IS} - \bB{\IS})^{-1/2} \bZ{m}{\IS}.
\end{align}
\end{subequations}
Introducing the so-called Mulliken notation for the bare 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}
and the corresponding (static) screened Coulomb potential matrix elements at coupling strength $\IS$
\begin{equation}
\W{pq,rs}{\IS} = \iint \MO{p}(\br{}) \MO{q}(\br{}) \W{}{\IS}(\br{},\br{}') \MO{r}(\br{}') \MO{s}(\br{}') \dbr{} \dbr{}',
\end{equation}
the BSE matrix elements read
\begin{subequations}
\begin{align}
\label{eq:LR_BSE-A}
\ABSE{ia,jb}{\IS} & = \delta_{ij} \delta_{ab} (\eGW{a} - \eGW{i}) + \IS \qty[ 2 \ERI{ia}{jb} - \W{ij,ab}{\IS} ],
\\
\label{eq:LR_BSE-B}
\BBSE{ia,jb}{\IS} & = \IS \qty[ 2 \ERI{ia}{bj} - \W{ib,aj}{\IS} ],
\end{align}
\end{subequations}
where $\eGW{p}$ are the $GW$ quasiparticle energies.
In the standard BSE approach, $\W{}{\IS}$ is built within the direct RPA scheme, \ie,
\begin{subequations}
\label{eq:wrpa}
\begin{align}
\W{}{\IS}(\br{},\br{}')
& = \int \frac{\epsilon_{\IS}^{-1}(\br{},\br{}''; \omega=0)}{\abs*{\br{}' - \br{}''}} \dbr{}'' ,
\\
\epsilon_{\IS}(\br{},\br{}'; \omega)
& = \delta(\br{}-\br{}') - \IS \int \frac{\chi_{0}(\br{},\br{}''; \omega)}{\abs*{\br{}' - \br{}''}} \dbr{}'' ,
\end{align}
\end{subequations}
with $\epsilon_{\IS}$ the dielectric function at coupling constant $\IS$ and $\chi_{0}$ the non-interacting polarizability. In the occupied-to-virtual orbital product basis, the spectral representation of $\W{}{\IS}$ can be written as follows in the case of real spatial orbitals
\begin{multline}
\label{eq:W}
\W{ij,ab}{\IS}(\omega) = \ERI{ij}{ab} + 2 \sum_m^{\Nocc \Nvir} \sERI{ij}{m}{\IS} \sERI{ab}{m}{\IS}
\\
\times \qty(\frac{1}{\omega - \OmRPA{m}{\IS} + i \eta} - \frac{1}{\omega + \OmRPA{m}{\IS} - i \eta}),
\end{multline}
where the spectral weights at coupling strength $\IS$ read
\begin{equation}
\sERI{pq}{m}{\IS} = \sum_i^{\Nocc} \sum_a^{\Nvir} \ERI{pq}{ia} (\bX{\IS}_m + \bY{\IS}_m)_{ia}.
\end{equation}
In the case of complex orbitals, we refer the reader to Ref.~\onlinecite{Holzer_2019} for a correct use of complex conjugation in the spectral representation of $\W{}{}$.
Note that, in the case of {\GOWO}, the RPA neutral excitations in Eq.~\eqref{eq:W} are computed using the HF orbital energies.
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_RPA-A}
\ARPA{ia,jb}{\IS} & = \delta_{ij} \delta_{ab} (\eHF{a} - \eHF{i}) + 2 \IS \ERI{ia}{jb},
\\
\label{eq:LR_RPA-B}
\BRPA{ia,jb}{\IS} & = 2 \IS \ERI{ia}{bj},
\end{align}
\end{subequations}
where $\eHF{p}$ are the Hartree-Fock (HF) orbital energies.
The relationship between the BSE formalism and the well-known RPAx (\ie, RPA with exchange) approach can be obtained by switching off the screening so that $\W{}{\IS}$ reduces to the bare Coulomb potential $\vc{}$.
In this limit, the $GW$ quasiparticle energies reduce to the HF eigenvalues, and Eqs.~\eqref{eq:LR_BSE-A} and \eqref{eq:LR_BSE-B} to the RPAx equations:
\begin{subequations}
\begin{align}
\label{eq:LR_RPAx-A}
\ARPAx{ia,jb}{\IS} & = \delta_{ij} \delta_{ab} (\eHF{a} - \eHF{i}) + \IS \qty[ 2 \ERI{ia}{jb} - \ERI{ij}{ab} ],
\\
\label{eq:LR_RPAx-B}
\BRPAx{ia,jb}{\IS} & = \IS \qty[ 2 \ERI{ia}{bj} - \ERI{ib}{aj} ].
\end{align}
\end{subequations}
%%% FIG 1 %%%
%\begin{figure}
% \includegraphics[width=\linewidth]{}
%\caption{
%\label{fig:}
%}
%\end{figure}
%%% %%% %%%
%%%%%%%%%%%%%%%%%%%%%%%%
\section{Conclusion}
\label{sec:conclusion}
%%%%%%%%%%%%%%%%%%%%%%%%
This is the conclusion
%%%%%%%%%%%%%%%%%%%%%%%%
\acknowledgements{
%%%%%%%%%%%%%%%%%%%%%%%%
This work was performed using HPC resources from GENCI-TGCC (Grant No.~2019-A0060801738) and CALMIP (Toulouse) under allocation 2020-18005.
Funding from the \textit{``Centre National de la Recherche Scientifique''} is acknowledged.
This work has also been supported through the EUR grant NanoX ANR-17-EURE-0009 in the framework of the \textit{``Programme des Investissements d'Avenir''.}}
%%%%%%%%%%%%%%%%%%%%%%%%
\section*{Supporting Information}
%%%%%%%%%%%%%%%%%%%%%%%%
See {\SI} for plenty of stuff
\bibliography{BSEdyn}
\end{document}

421
SI/bsedyn-SI.tex Normal file
View File

@ -0,0 +1,421 @@
\documentclass[10pt]{book}
\usepackage{graphicx}
\usepackage{amsmath} %% pour mettre du texte en mode math
\usepackage{tcolorbox}
\newcommand\uom{\underline{\omega}}
\newcommand\Oms{{\Omega}_s}
\newcommand\hOms{\frac{{\Omega}_s}{2}}
\begin{document}
\chapter{Strinati for the big kids}
We wish to solve the self-consistent Bethe-Salpeter equation:
\begin{equation*}
L(1,2;1',2') = L_0(1,2;1',2') + \int d3456 \; L_0(1,4;1',3) \Xi[3,5;4,6] L(6,2;5,2')
\end{equation*}
where e.g. $1=(x_1,t_1)$ is a space/spin and time position and:
\begin{align*}
iL(1,2;1',2') &= -G_2(1,2 ; 1',2') + G(1,1')G(2,2') \\
iL_0(1,2;1',2') &= G(1,2')G(2,1') \\
iG(1,2') &= \langle N | T {\hat \psi}(1) {\hat \psi}^{\dagger}(2') | N \rangle \\
i^2 G_2(1,2 ; 1',2') &= \langle N | T \Big[ {\hat \psi}(1) {\hat \psi}(2) {\hat \psi}^{\dagger}(2') {\hat \psi}^{\dagger}(1') \Big] | N \rangle
\end{align*}
For optical properties, with instantaneous creation/destruction of excitations, we are interested in the collapse of the time variables starting with imposing $t_2'=t_2^{+}$ in the right-hand side $L(1,2;1',2')$ and left hand-side $L(6,2;5,2')$ response functions. Since $t_2'=t_2^{+}=(t_2 + 0^+)$ indicates that the (2,2') operators cannot be separated, we see that the ordering operator selects two possible channels depending on the times ordering, namely:
\begin{align*}
i^2 G_2^{I}(1,2 ; 1',2') &= \theta( t_{m}^{11'} -t_2 )\langle N | T [{\hat \psi}(1) {\hat \psi}^{\dagger}(1') ] T [ {\hat \psi}(2) {\hat \psi}^{\dagger}(2') ] | N \rangle \\
i^2 G_2^{II}(1,2 ; 1',2') &= \theta( t_2 - t_{M}^{11'} )\langle N | T [ {\hat \psi}(2) {\hat \psi}^{\dagger}(2') ]
T [{\hat \psi}(1) {\hat \psi}^{\dagger}(1') ] | N \rangle
\end{align*}
with $t_{M}^{11'} = \max(t_1,t_{1'}) = | \tau_{11'} |/2 + t^{11'}$ and $t_{m}^{11'} = \min(t_1,t_{1'}) = - | \tau_{11'} |/2 + t^{11'}$ where $\tau_{11'} = t_1 - t_{1'}$ and $t^{11'} = (t_1 + t_{1'})/2$.
Introducing the complete set of eigenstates of the N-electron systems $\lbrace | N,s \rangle \rbrace$, where s index the eigenstates and $| N,s=0 \rangle = | N \rangle$ the ground-state, one obtains for $G^{I}$ :
\begin{align*}
i^2 G_2^{I}(1,2 ; 1',2') &= \theta( t_{m}^{11'} -t_2 ) \sum_{s>0} \langle N | T {\hat \psi}(1) {\hat \psi}^{\dagger}(1') | N,s \rangle
\langle N,s | T {\hat \psi}(2) {\hat \psi}^{\dagger}(2') | N \rangle \\
& + \theta( t_{m}^{11'} -t_2 ) \frac{ G(1,1') }{ i } \frac{ G(2,2') }{ i } \\
\end{align*}
where in the right-hand side of the first line the sum $(s>0)$ avoids the (s=0) ground-state contribution that builds the 1-body Green's functions product in the 2nd line.
This yields a first channel for $L(1,2;1',2')$ namely:
$$
iL^{I}(1,2;1',2') = \theta( t_{m}^{11'} -t_2 ) \sum_{s>0} \langle N | T {\hat \psi}(1) {\hat \psi}^{\dagger}(1') | N,s \rangle
\langle N,s | T {\hat \psi}(2) {\hat \psi}^{\dagger}(2') | N \rangle
$$
Switching from the Heisenberg representation (time-dependent operators) to the Schr{\"o}dinger one for the field operators, with e.g.
$$
{\hat \psi}(2) = {\hat \psi}(x_2 t_2) = e^{ i H t_2} {\hat \psi}(x_2) e^{ -i H t_2}
$$
one obtains
\begin{align*}
\langle N,s | T {\hat \psi}(2) {\hat \psi}^{\dagger}(2') | N \rangle &=
- \langle N,s | e^{i H t_2^+} {\hat \psi}^{\dagger}(x_2') e^{i H (t_2 - t_2^+) } {\hat \psi}(x_2) e^{-i H t_2} | N \rangle \\
&= - \langle N,s | {\hat \psi}^{\dagger}(x_{2'}) {\hat \psi}(x_2) | N \rangle \times e^{i \Omega_s t_2 } \\
& = {\tilde \chi}_s(x_2,x_{2'}) e^{i \Omega_s t_2 }
\end{align*}
with $\;{\tilde \chi}_s(x_2,x_{2'}) = \langle N,s | {\hat \psi}(x_2) {\hat \psi}^{\dagger}(x_{2'}) | N \rangle$.
The energy $\Omega_s = E_s^N - E_0^N$ is the s-th neutral excitation energy of the system, namely the energy we are looking for.
The (-) sign in the right-hand side comes from the commutation of the creation and destruction operators with $t_2^+ > t_2$.
A similar calculation for the second channel [$t_2 > \max(t_1,t_{1'})$] leads to a frequency dependence in $e^{ -i \Omega_s t_2 }$.
Putting everything together, one obtains :
\begin{align*}
iL (1,2;1',2') = & \; \theta( t_{m}^{11'} -t_2 ) \sum_{s>0} \langle N | T {\hat \psi}(1) {\hat \psi}^{\dagger}(1') | N,s \rangle
{\tilde \chi}_s(x_2,x_{2'}) \times e^{i \Omega_s t_2 } \\
- & \; \theta( t_2 - t_{M}^{11'} ) \sum_{s>0} \langle N,s | T {\hat \psi}(1) {\hat \psi}^{\dagger}(1') | N \rangle
{\chi}_s(x_2,x_{2'}) \times e^{ -i \Omega_s t_2 }
\end{align*}
with $\;{\chi}_s(x_2,x_{2'}) = \langle N | {\hat \psi}(x_2) {\hat \psi}^{\dagger}(x_{2'}) | N,s \rangle$. \\
Performing the very same analysis for $L(6,2;5,2')$ and picking up a specific $e^{ +i \Omega_s t_2 }$ contribution, dividing by
$ \langle N,s | T {\hat \psi}(2) {\hat \psi}^{\dagger}(2') | N \rangle $ both sides of the Bethe-Salpeter equation leads to: \\
\begin{tcolorbox}[ width=\textwidth, colback={white}]
\begin{align}
\langle N & | T {\hat \psi}(1) {\hat \psi}^{\dagger}(1') | N,s \rangle \theta( t_{m}^{11'} -t_2 ) = \\
&= \int d3456 \; L_0(1,4;1',3) \Xi[3,5;4,6] \langle N | T {\hat \psi}(6) {\hat \psi}^{\dagger}(5) | N,s \rangle \theta( t_{m}^{56} -t_2 ) \nonumber
\end{align}
\end{tcolorbox}
\noindent where the $L_0(1,2;1',2')$ does not contribute if \textbf{ we select some $\Omega_s$ frequency below the quasiparticle gap of the system}. This is the common situation in molecular systems where due to excitonic effects, namely electron-hole interaction, the lowest optical excitation energies are lower than the photoemission quasiparticle gap. $ L_0(1,2;1',2')$ has indeed a lowest excitation energy at the energy gap $E_g$ and cannot participate to the $e^{i \Omega_s t_2}$ response of the system if $\Omega_s < E_g$. \\
This is equation (11.3) of Strinati with e.g. $$ \langle N | T {\hat \psi}(1) {\hat \psi}^{\dagger}(1') | N,s \rangle = \chi_s(x_1,x_{1'}; \tau_1) e^{-i \Omega_s t^1}$$
(see equation G.3a of Appendix G) but keeping the Heaviside fonctions for a better definition of time integrals. \\
\noindent \textbf{Spectral representation of $G(t_1-t_3)G(t_4-t_1^+)$ } \\
We now consider the $t_1$ time. Using the same algebra as above, one obtains
$$
\langle N | T {\hat \psi}(1) {\hat \psi}^{\dagger}(1') | N,s \rangle =
\langle N | T {\hat \psi}(x_1) {\hat \psi}^{\dagger}(x_{1'}) | N,s \rangle \times e^{- i \Oms t_1 }
$$
imposing now $t_1' = t_1^+$. To extract the $e^{-i \Omega_s t_1 }$ behavior in the right hand-side of the BSE equation, one must find the Fourier components with respect to the $t_1$ variable of $L_0(1,4;1',3) = -i G(1,3) G(4,1')$ by multiplying by $ \int dt_1 e^{i \omega_1 t_1}$ both sides of the BSE equation, \textcolor{red}{ with $( \omega_1 \rightarrow \Oms + i\delta )$ where we introduce the small infinitesimal $\delta = 0^+$ for defining properly the $t_1$ time integral on the left-hand side: }
$$
\int dt_1 e^{i \omega_1 t_1} e^{ - i \Oms t_1 } \theta( t_1 -t_2 ) = e^{-\delta t_2} \int d \tau e^{- \delta \tau } \theta( \tau ) = ( \delta^{-1} )
$$
This yields an annoying but useful factor that will canceled out in the following.
Expressing again the field operators in their Schr\"{o}dinger representation, one see that the 1-body Green's functions depend on the time-differences $\tau_{13} = t_1 - t_3$ and $\tau_{41} = t_4 - t_1$, respectively, with e.g.
$$
\langle N | \psi(1) \psi^{\dagger}(3) | N \rangle = \langle N | \psi(x_1) e^{-i H \tau_{13}} \psi^{\dagger}(x_3) | N \rangle e^{+i E_N^0 \tau_{13}}
$$
Taking by convention the time difference as that of the destruction operator minus that of the creation analog and introducing our notations for the Fourier representation :
\begin{equation*}
G( \tau_{13} ) = {1 \over 2\pi} \int d \omega \; G(\omega) e^{ - i \omega \tau_{13} } \;\;\;\;\;
G(\tau_{41} ) = {1 \over 2\pi} \int d \omega' \; G(\omega') e^{ - i {\omega}' \tau_{41} }
\end{equation*}
one obtains (keeping only time/frequency variables for clarity)
$$
i L_0(1,4;1',3) = \int { d\omega \over 2\pi } { d{\omega}' \over 2\pi } \; G(\omega) G({\omega}') e^{ -i \omega \tau_{13} } e^{ - i {\omega}' \tau_{41} }
$$
The Fourier components of ($iL_0$) with respect to the ($t_1$) time are:
\begin{align*}
[iL_0]( \omega_1 ) &= \int dt_1 e^{ i \omega_1 t_1 } \int { d\omega \over 2\pi } { d{\omega}' \over 2\pi } \; G(\omega ) G( {\omega}' )
e^{- i \omega (t_1 - t_3) } e^{ -i {\omega}' (t_4 - t_1) } \\
&= {1 \over (2\pi)^2 } \int d \omega d{\omega}' \; G(\omega ) G( {\omega}' ) 2\pi \delta( -\omega + \omega_1 + {\omega}' )
e^{ i \omega t_3 } e^{ - i {\omega}' t_4 } \\
&= {1 \over 2\pi } \int d \omega \; G(\omega ) G( {\omega} - \omega_1 ) e^{ i \omega t_3 } e^{ -i ( \omega - \omega_1) t_4 } \\
&= { 1 \over 2\pi } \int d \omega \; G(\omega - {\omega_1 \over 2} ) G( {\omega} + {\omega_1 \over 2} ) e^{ i \omega \tau_{34} } e^{ i \omega_1 t^{34} }
\end{align*}
where in the last line we made the change of variable $( \omega \rightarrow \omega + {\omega_1} / 2 ) $ and introduced the average time
$t^{34} = (t_3 + t_4) /2$.
We now make the identification $( \omega_1 = {\Oms} + i\delta )$ as done to properly pick-up the left-hand side $e^{-i \Oms t_1}$ contribution and obtain :
\begin{align*}
[iL_0 ]( \Omega_s)
= { 1 \over 2\pi } \int d \omega \; G(\omega - \hOms ) G( {\omega} + \hOms ) e^{ i \omega \tau_{34} } e^{ i \Oms^{+} t^{34} }
\end{align*}
with $( \Oms^{+} = \Oms + i\delta )$ keeping the infinitesimal $\delta$ where it will be needed.
Using the Lehman representation of the 1-body Green's functions, e.g.
$$
G(x_1,x_3 ; \omega) = \sum_n \frac{ \phi_n(x_1) \phi_n^*(x_3) } { \omega - \varepsilon_n + i \eta \text{sgn}(\varepsilon_n - \mu) }
$$
where we reintroduced the (space, spin) variables, one obtains by picking the residues associated with the poles of each Green's function:
\begin{align*}
\int d \omega & \; G(x_1,x_3; \omega - \hOms ) G(x_4,x_{1'}; \omega + \hOms ) e^{ i \omega \tau_{34} } = \\
& 2 {i} \pi \theta( \tau_{34} ) \sum_{nj} \frac{ \phi_n(x_1) \phi_n^*(x_3) \phi_j(x_4) \phi_j^*(x_{1'}) } { \varepsilon_j + \Omega_s - \varepsilon_n + i \eta \times \text{sgn}(\varepsilon_n - \mu) } e^{i (\varepsilon_j + {\Omega_s \over 2} ) \tau_{34} } \\
+ & 2 {i} \pi \theta( \tau_{34} ) \sum_{jn} \frac{ \phi_j(x_1) \phi_j^*(x_3) \phi_n(x_4) \phi_n^*(x_{1'}) } { \varepsilon_j - \Omega_s - \varepsilon_n + i \eta \times \text{sgn}(\varepsilon_n - \mu) } e^{i (\varepsilon_j - {\Omega_s \over 2} ) \tau_{34} } \\
- & 2 {i} \pi \theta(- \tau_{34} ) \sum_{nb} \frac{ \phi_n(x_1) \phi_n^*(x_3) \phi_b(x_4) \phi_b^*(x_{1'}) } { \varepsilon_b + \Omega_s - \varepsilon_n + i \eta \times \text{sgn}(\varepsilon_n - \mu) } e^{i (\varepsilon_b + {\Omega_s \over 2} ) \tau_{34} } \\
- & 2 {i} \pi \theta(- \tau_{34} ) \sum_{bn} \frac{ \phi_b(x_1) \phi_b^*(x_3) \phi_n(x_4) \phi_n^*(x_{1'}) } { \varepsilon_b - \Omega_s - \varepsilon_n + i \eta \times \text{sgn}(\varepsilon_n - \mu) } e^{i (\varepsilon_b - {\Omega_s \over 2} ) \tau_{34} }
\end{align*}
where for $\tau_{34} > 0$ (resp. $\tau_{34} < 0$) the contour is close in the upper (resp. lower) half complex plane.
We adopted chemist notations with (i,j) and (a,b) indexing occupied and virtual states, respectively. Projecting now on
$ \phi_a^*(x_1) \phi_i(x_{1'}) \; $ one obtains by orthonormalization :
\begin{align*}
\int dx_1 dx_{1'} \; \phi_a^*(x_1) \phi_i(x_{1'}) \int d \omega & \; G(x_1,x_3; \omega + {\Omega_s \over 2} ) G(x_4,x_{1'}; \omega - {\Omega_s \over 2} ) e^{ i \omega \tau_{34} } = \\
2 & {i} \pi \theta( \tau_{34} ) \frac{ \phi_a^*(x_3) \phi_i(x_4) } { \varepsilon_i + \Omega_s - \varepsilon_a + i \eta } e^{i (\varepsilon_i + {\Omega_s \over 2} ) \tau_{34} } \\
-2 & {i} \pi \theta( - \tau_{34} ) \frac{ \phi_a^*(x_3) \phi_i(x_4) } { \varepsilon_a - \Omega_s - \varepsilon_i - i \eta } e^{i (\varepsilon_a - {\Omega_s \over 2} ) \tau_{34} }
\end{align*}
After multiplication by $ ( \varepsilon_a - \varepsilon_i - \Omega_s -i \eta )$ (and setting $\eta \rightarrow 0$), the Bethe-Salpeter equation becomes: \\
\begin{tcolorbox}[ width=\textwidth, colback={white}]
\begin{align*}
( \varepsilon_a - \varepsilon_i & - \Omega_s ) \int dx_1 \phi_a^*(x_1) \phi_i(x_{1'}) \langle N | T {\hat \psi}(x_1) {\hat \psi}^{\dagger}(x_{1'}) | N,s \rangle ( \delta^{-1} )
= \nonumber \\
& - \int d34 \; \phi_a^*(x_3) \phi_i(x_4)
\left[ \theta( \tau_{34} ) e^{i (\varepsilon_i + {\Omega_s \over 2} ) \tau_{34} } + \theta( - \tau_{34} ) e^{i (\varepsilon_a - {\Omega_s \over 2} ) \tau_{34} } \right] \\
& \times
e^{ i \Omega_s^{+} t^{34} } \int d56 \; \Xi[3,6;4,5] \langle N | T {\hat \psi}(5) {\hat \psi}^{\dagger}(6) | N,s \rangle \times \theta( t^{56}_m - t_2 )
\end{align*}
\end{tcolorbox}
\vskip .3cm \noindent This is equation (11.5) by Strinati with the same
$$
\theta( \tau ) e^{i (\varepsilon_i + {\Omega_s \over 2} ) \tau } + \theta( - \tau ) e^{ i (\varepsilon_a - {\Omega_s \over 2} ) \tau } =
e^{i {\Omega_s \over 2} \tau } \Big( \theta( \tau ) e^{i \varepsilon_i \tau } + \theta( - \tau ) e^{i (\varepsilon_a - \Omega_s )\tau } \Big)
$$
time-structure where $\tau = \tau_{34} = t_3 - t_4 $. \\
\noindent{\textbf{The $GW$ approximation. } Adopting the $GW$ approximation, and neglecting the $(\partial W / \partial G)$ of higher order in $W$, one obtains:
\begin{equation}
\Xi[3,6;4,5] = v(3,5) \delta(3,4) \delta(5,6) - W(3^+,4) \delta(3,5) \delta(4,6)
\end{equation}
so that the Bethe-Salpeter equation becomes:
\begin{align*}
( \varepsilon_a - \varepsilon_i & - \Omega_s ) \int dx_1 \phi_a^*(x_1) \phi_i(x_{1'}) \langle N | T {\hat \psi}(x_1) {\hat \psi}^{\dagger}(x_{1'}) | N,s \rangle ( \delta^{-1} )
= \nonumber \\
& - \int d35 \; \phi_a^*(x_3) \phi_i(x_3) e^{ i \Omega_s^{+} t_3 } v(3,5) \langle N | T [{\hat \psi}(5) {\hat \psi}^{\dagger}(5) ] | N,s \rangle \times \theta( t_5 - t_2 ) \\
& + \int d34 \; \phi_a^*(x_3) \phi_i(x_4)
\left[ \theta( \tau_{34} ) e^{i (\varepsilon_i + {\Omega_s \over 2} ) \tau_{34} } + \theta( - \tau_{34} ) e^{i (\varepsilon_a - {\Omega_s \over 2} ) \tau_{34} } \right] \times \\
& \hskip 1cm \times
e^{ i \Omega_s^{+} t^{34} } W(3^+,4) \langle N | T [{\hat \psi}(3) {\hat \psi}^{\dagger}(4) ] | N,s \rangle \times \theta( t^{34}_m - t_2 )
\end{align*}
\\
\noindent\textbf{Spectral representation of $\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) | N,s \rangle$ } \\
\noindent Starting with:
\begin{align*}
\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) | N,s \rangle &= \theta( \tau_{34} ) \langle N | {\hat \psi}(3) {\hat \psi}^{\dagger}(4) | N,s \rangle \\
&- \theta( -\tau_{34} ) \langle N | {\hat \psi}^{\dagger}(4) {\hat \psi}(3) | N,s \rangle
\end{align*}
and decomposing the static (Schr\"{o}dinger) field operators over the Hartree-Fock molecular orbitals creation/destruction operators:
$$
{\hat \psi}^{\dagger}(x_4) = \sum_n \phi_n^*(x_4) {\hat a}_n^{\dagger} \;\;\; \text{ and} \;\;\; {\hat \psi}(x_3) = \sum_m \phi_m(x_3) {\hat a}_m $$
one obtains
\begin{align*}
\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) | N,s \rangle &= \sum_{mn} \phi_m(x_3) \phi_n^*(x_4) \times \\
\times & \Big( \theta( \tau_{34} ) \langle N | {\hat a}_m e^{- i H \tau_{34} } {\hat a}^{\dagger}_n | N,s \rangle e^{-i E_s^N t_4} e^{+ i E_0^N t_3} \\
&- \theta( -\tau_{34} ) \langle N | {\hat a_n}^{\dagger} e^{ i H \tau_{34} } {\hat a}_m | N,s \rangle e^{ - i E_s^N t_3} e^{ + i E_0^N t_4} \Big)
\end{align*}
Taking the conjugate of the $ {\hat a}_m e^{- i H \tau_{34} } $ and
$ {\hat a_n}^{\dagger} e^{ i H \tau_{34} } $ operators to act on the ground-state $| N \rangle$ state, one obtains:
$$
e^{ i H \tau } {\hat a}_m^{\dagger} | N \rangle = e^{ i (E_0^N + \varepsilon_m ) \tau } {\hat a} _m^{\dagger} | N \rangle
\;\;\;\;\; \text{and} \;\;\;\;
e^{ - i H \tau } {\hat a}_n | N \rangle = e^{ - i ( E_0^N - \varepsilon_n ) \tau } {\hat a}_n | N \rangle
$$
This properties hold in the case of Hartree-Fock formalism (Koopman's theorem) but holds more generally if $\lbrace \varepsilon_m , \varepsilon_n \rbrace$ are quasiparticle energies,
namely true electronic removal or addition energies such as in the $GW$ formalism.
Forming the corresponding bra, one obtains :
\begin{align*}
\langle N | T {\hat \psi}(3) & {\hat \psi}^{\dagger}(4) | N,s \rangle = - \sum_{mn} \phi_m(x_3) \phi_n^*(x_4) \langle N | {\hat a}_n^{\dagger} {\hat a}_m | N,s \rangle \times \\
\times & \Big( \theta( \tau_{34} ) e^{- i \varepsilon_m \tau_{34} } e^{ -i \Oms t_4 }
+ \theta( -\tau_{34} ) e^{ - i \varepsilon_n \tau_{34} } e^{ - i \Oms t_3 } \Big)
\end{align*}
With $t_3 = \tau_{34}/2 + t^{34}$ and $t_4 = -\tau_{34}/2 + t^{34}$, one finally obtains :
\begin{align*}
\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) & | N,s \rangle = - \Big( e^{ -i \Omega_s t^{34} } \Big) \sum_{mn} \phi_m(x_3) \phi_n^*(x_4) \langle N | {\hat a}_n^{\dagger} {\hat a}_m | N,s \rangle \times \nonumber \\
\times & \Big( \theta( \tau_{34} ) e^{- i ( \varepsilon_m - \hOms ) \tau_{34} }
+ \theta( -\tau_{34} ) e^{ - i ( \varepsilon_n + \hOms) \tau_{34} } \Big)
\end{align*}
We must now discuss the (m,n) indexes. The most natural channel for ${\hat a}_m^{\dagger} {\hat a}_n | N \rangle$ to project on an excited state occurs
for (n,m) indexing (occupied,virtual) orbitals, respectively, in order to promote an electron from the occupied to the empty orbitals.
We will label $A_{jb}^{s} = \langle N | {\hat a}_j^{\dagger} {\hat a}_b | N,s \rangle$ these large amplitude projections. However, in the case of a generic correlated system,
the ground-state does not project exclusively on the mono-determinant built from the occupied 1-body orbitals $\lbrace \phi_i \rbrace$, but contains also contributions from
singly, doubly, etc. excited determinants in the full $\lbrace \phi_n \rbrace$ space. As such, the $\langle N | {\hat a}_b^{\dagger} {\hat a}_j | N,s \rangle$ weights are small in general but not strictly zero.
We will label $B_{jb}^{s}$ these small amplitudes. Neglecting these $B_{jb}^{s}$ coefficients leads to the Tamm-Dancoff approximation (TDA). We start with the TDA in the following for simplicity, adding the "small B" contribution in a second step. \\
\noindent In the TDA, namely with $B_{jb}^{s} = 0$, one obtains :
\begin{align*}
\langle N | T {\hat \psi}(3) & {\hat \psi}^{\dagger}(4) | N,s \rangle \cong - \Big( e^{ -i \Omega_s t^{34} } \Big) \sum_{jb} \phi_b(x_3) \phi_j^*(x_4) A_{jb}^{s} \times \nonumber \\
\times & \Big( \theta( \tau_{34} ) e^{- i ( \varepsilon_b - \hOms ) \tau_{34} }
+ \theta( -\tau_{34} ) e^{ - i ( \varepsilon_j + \hOms) \tau_{34} } \Big)
\end{align*}
As a result, the integral associated with the screened Coulomb potential reads:
\begin{align*}
\int d34 \; \phi_a^*(x_3) & \phi_i(x_4)
\left[ \theta( \tau_{34} ) e^{i (\varepsilon_i + {\Omega_s \over 2} ) \tau_{34} } + \theta( - \tau_{34} ) e^{i (\varepsilon_a - {\Omega_s \over 2} ) \tau_{34} } \right] \times \\
& \times
e^{ i \Omega_s^{+} t^{34} } W(3^+,4) \langle N | T [{\hat \psi}(3) {\hat \psi}^{\dagger}(4) ] | N,s \rangle \times \theta( t^{34}_m - t_2 ) \\
= - \int d34 \; \phi_a^*(x_3) & \phi_i(x_4) \theta( t^{34}_m - t_2 ) W(3^+,4) \sum_{jb} \phi_b(x_3) \phi_j^*(x_4) A_{jb}^{s} e^{ - \delta t^{34} } \times \\
& \times \left[ \theta( \tau_{34} ) e^{i ( \varepsilon_i - \varepsilon_b + \Oms ) \tau_{34} } + \theta( - \tau_{34} ) e^{i (\varepsilon_a - \varepsilon_j - \Oms ) \tau_{34} } \right] \\
= - \sum_{jb} A_{jb}^{s} & \int d\tau_{34} dt^{34} e^{ - \delta t^{34} } \theta( t^{34} - | \tau_{34} |/2 - t_2 ) W_{ij,ab}( \tau_{34}^{+} ) \times \\
& \times \left[ \theta( \tau_{34} ) e^{i ( \varepsilon_i - \varepsilon_b + \Oms ) \tau_{34} } + \theta( - \tau_{34} ) e^{i (\varepsilon_a - \varepsilon_j - \Oms ) \tau_{34} } \right]
\end{align*}
where we have explicited $t^{34}_m$ and made a change of variable $(t_3,t_4) \rightarrow (\tau_{34}, t^{34})$ with unit wronskian, defining the screened Coulomb matrix elements
\begin{align*}
W_{ij,ab}( \tau_{34}^{+} )&= \int d{\bf r}_3 d{\bf r}_4 \; \phi_i(x_4) \phi_j^*(x_4) W(3^+,4) \phi_a^*(x_3) \phi_b(x_3) \;\;\;\;
\end{align*}
with $(\tau_{34}^{+} = t_3^{+} - t_4)$ and adopting Mulliken notations where we group together the index associated with orbitals taken at the same space position (and putting complexe conjugate orbitals as inner indexes).
The $t^{34}$ time integral yields again $( \delta^{-1} )$ as done previously for the $t_1$-time integral, leaving only the dependence on the time difference $\tau_{34} = t_3 - t_4$.
We can understand the origin of this problematic $( \delta^{-1} )$ factor : the screened- Coulomb potential contribution only depends on the time difference $\tau_{34}$ so that the time variable $t^{34}$ becomes problematic.
To deal now with the bare Coulomb contribution, we need
$\langle N | T {\hat \psi}(5) {\hat \psi}^{\dagger}(5) | N,s \rangle$ that in the TDA reads :
$$
\langle N | T {\hat \psi}(5) {\hat \psi}^{\dagger}(5) | N,s \rangle \cong - \Big( e^{ -i \Omega_s t_5 } \Big) \sum_{jb} \phi_b(x_5) \phi_j^*(x_5) A_{jb}^{s}
$$
leading to the following contribution :
\begin{align*}
- \int d35 \; \phi_a^*(x_3) & \phi_i(x_3) e^{ i \Omega_s^{+} t_3 } v(3,5) \langle N | T [{\hat \psi}(5) {\hat \psi}^{\dagger}(5) ] | N,s \rangle \times \theta( t_5 - t_2 ) \\
& = \sum_{jb} A_{jb}^{s} v_{ia,jb} \int dt_3 dt_5 \delta( t_3 - t_5) e^{ i \Omega_s^{+} t_3 } e^{ -i \Omega_s t_5 }\theta( t_5 - t_2 )
\end{align*}
with $v(3,5) = v( {\bf r}_3 - {\bf r}_5 ) \delta( t_3 - t_5)$ a static interaction and :
$$
v_{ia,jb} = \int d{\bf r}_3 d{\bf r}_5 \; \phi_i(x_3) \phi_a^*(x_3) v( {\bf r}_3 - {\bf r}_5 ) \phi_j^*(x_5) \phi_b(x_5)
$$
The time integral reads now $\int dt_3 e^{- \delta t_3} \theta(t_3 - t_2 )$ yielding $( \delta^{-1} )$ again.
This factor in front of all terms can now be factorized out. \\
\noindent Finally, for the right-hand side of the Bethe-Salpeter equation, we develop:
\begin{equation*}
\langle N | {\hat \psi}^{\dagger}(x_{1'}) {\hat \psi}(x_1) | N,s \rangle = \sum_{jb} \phi_b(x_1) \phi_j^*(x_1) A_{jb}^{s}
\end{equation*}
so that by orthonormalization:
\begin{equation*}
\int dx_1 dx_{1'} \; \phi_a^*(x_1) \phi_i(x_{1'}) \langle N | {\hat \psi}^{\dagger}(x_{1'}) {\hat \psi}(x_1) | N,s \rangle =
\sum_{jb} A_{jb}^{s} \delta_{ij} \delta_{ab} = A_{ia}^{s}
\end{equation*}
Note that the projection on $\phi_a^*(x_1) \phi_i(x_{1'})$ will always cancel the small $B_{jb}^{s}$ contribution so that this last expression is always valid.
We obtain thus a much simplified reformulation of the dynamical TDA Bethe-Salpeter equation: \\
\begin{tcolorbox}[ width=\textwidth, colback={white}]
\begin{align*}
( \varepsilon_a - \varepsilon_i - \Omega_s ) & A_{ia}^{s}
+ \sum_{jb} A_{jb}^{s} \times v_{ai,bj} \stackrel{\text{TDA}}{=} \sum_{jb} A_{jb}^{s} \times \int d\tau \; W_{ij,ab}( {\tau}^{+} ) \times \\
& \times \left[ \theta( \tau ) e^{i (\varepsilon_i - \varepsilon_b+ \Oms ) \tau } + \theta( - \tau ) e^{i (\varepsilon_a - \varepsilon_j - \Oms ) \tau } \right]
\end{align*}
\end{tcolorbox}
\noindent where we dropped the index of ${\tau}_{34}$. Taking finally the Fourier representation of $W$, the time integration reads:
\begin{align*}
\widetilde{W}_{ij,ab}(\Oms) =& \int d\tau \; W_{ij,ab}( \tau^{+} )
\times \; \left[ \theta( \tau ) e^{i (\varepsilon_i - \varepsilon_b+ \Oms ) \tau } + \theta( - \tau ) e^{i (\varepsilon_a - \varepsilon_j - \Oms ) \tau } \right] \\
=& \int { d\omega \over 2\pi } W_{ij,ab}( \omega ) \int d\tau \; e^{-i \omega ( \tau + 0^+) }
\left[ \theta( \tau ) e^{i (\varepsilon_i - \varepsilon_b+ \Oms ) \tau } + \theta( - \tau ) e^{i (\varepsilon_a - \varepsilon_j - \Oms ) \tau } \right] \\
=& { i \over 2\pi} \int d\omega e^{-i \omega 0^{+}} \left[ \frac{ W_{ij,ab}( \omega ) }{ ( \Oms - \omega ) - ( \varepsilon_b - \varepsilon_i) + i\eta }
+ \frac{ W_{ij,ab}( \omega ) }{ ( \Oms + \omega) - ( \varepsilon_a - \varepsilon_j ) + i\eta } \right]
\end{align*}
where we have introduced again a small damping $( \Oms \rightarrow \Oms + i\eta)$ of the excitation for convergency. This leads to
equation (11.8) by Strinati: \\
\begin{tcolorbox}[ width=\textwidth, colback={white}]
\begin{align*}
( \varepsilon_a - \varepsilon_i - \Omega_s ) & A_{ia}^{s}
+ \sum_{jb} A_{jb}^{s} \times v_{ai,bj} \stackrel{\text{TDA}}{=} { i \over 2 \pi} \sum_{jb} A_{jb}^{s} \int d\omega e^{-i \omega 0^+ } W_{ij,ab} \times \\
\hskip 1cm &\times \left[ \frac{1}{( \Oms - \omega) - (\varepsilon_b - \varepsilon_i) +i \eta } + \frac{1}{ (\Oms + \omega) - (\varepsilon_a - \varepsilon_j) + i\eta } \right]
\end{align*}
\end{tcolorbox}
\vskip .3cm \noindent or to a form similar to the static limit : \\
\begin{tcolorbox}[ width=\textwidth, colback={white}]
\begin{align}
( \varepsilon_a - \varepsilon_i - \Omega_s ) A_{ia}^{s}
+ \sum_{jb} \Big( v_{ai,bj} - \widetilde{W}_{ij,ab}(\Oms) \Big) A_{jb}^{s} \stackrel{\text{TDA}}{=} 0
\end{align}
\end{tcolorbox}
\vskip .3cm \noindent but with an effective dynamical screened Coulomb potential (see Pina eq. 24): \\
\begin{tcolorbox}[ width=\textwidth, colback={white}]
\begin{align}
\widetilde{W}_{ij,ab}(\Oms) &= { i \over 2 \pi} \int d\omega \; e^{-i \omega 0^+ } W_{ij,ab} \times \\
\hskip 1cm &\times \left[ \frac{1}{ (\Oms - \omega) - (\varepsilon_b - \varepsilon_i) +i \eta } + \frac{1}{ (\Oms + \omega) - (\varepsilon_a - \varepsilon_j) + i\eta } \right] \nonumber
\end{align}
\end{tcolorbox}
\vskip .3cm \noindent Taking now the spectral representation for the (time-ordered) $W_{ij,ab}( \omega )$, namely:
\begin{align*}
W_{ij,ab}(\omega) &= (ij|ab) + 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
& \times \Big( \frac{1}{ \omega-\Omega_m^{RPA} + i\eta } - \frac{1}{ \omega + \Omega_m^{RPA} - i\eta } \Big)
\end{align*}
the frequency integral becomes, closing in the lower half plane thanks to the $e^{-i\omega 0^+}$ factor.
\textcolor{red}{(This is the only place where the (+) in $3^+$ makes a difference as for the $GW$ operator ... check carefully origin ...)}:
\begin{align*}
\int d\omega e^{-i \omega 0^{+}} & \Big( \frac{ W_{ij,ab}( \omega ) }{ ( \Oms - \omega ) - ( \varepsilon_b - \varepsilon_i) + i\eta }
+ \frac{ W_{ij,ab}( \omega ) }{ ( \Oms + \omega) - ( \varepsilon_a - \varepsilon_j ) + i\eta } \Big)= \\
= & (-2i \pi) W_{ij,ab}( - \Oms + (\varepsilon_a - \varepsilon_j) ) + (-2i \pi) \times 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
& \times \Big( \frac{ 1 }{ ( \Oms - \Omega_m^{RPA} ) - ( \varepsilon_b - \varepsilon_i) + i\eta }
+ \frac{ 1 }{ ( \Oms + \Omega_m^{RPA} ) - ( \varepsilon_a - \varepsilon_j ) + i\eta } \Big)
\end{align*}
that is :
\begin{align*}
\widetilde{W}_{ij,ab}( \Oms ) = & (ij|ab) + 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
\times & \Big( \frac{ 1}{ -\Oms + (\varepsilon_a - \varepsilon_j) -\Omega_m^{RPA} + i\eta } - \frac{ 1}{ -\Oms + (\varepsilon_a - \varepsilon_j) + \Omega_m^{RPA} - i\eta } \\
+ & \frac{ 1 }{ ( \Oms - \Omega_m^{RPA} ) - ( \varepsilon_b - \varepsilon_i) + i\eta }
+ \frac{ 1 }{ ( \Oms + \Omega_m^{RPA} ) - ( \varepsilon_a - \varepsilon_j ) + i\eta } \Big)
\end{align*} \\
or
\begin{align*}
\widetilde{W}_{ij,ab}( \Oms ) = & (ij|ab) + 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
\times & \Big( \frac{ 1 }{ ( \Oms - \Omega_m^{RPA} ) - ( \varepsilon_b - \varepsilon_i) + i\eta } + \frac{ 1}{ ( \Oms - \Omega_m^{RPA} ) - (\varepsilon_a - \varepsilon_j) + i\eta } \\
+ & \frac{ 1 }{ ( \Oms + \Omega_m^{RPA} ) - ( \varepsilon_a - \varepsilon_j ) + i\eta } - \frac{ 1}{ ( \Oms +\Omega_m^{RPA} ) - (\varepsilon_a - \varepsilon_j) - i\eta } \Big)
\end{align*} \\
Keeping only \textbf{ the two resonant contributions : } (ACTUALLY the other ones in the last line cancel when $\eta \rightarrow 0$ !! check ...) one obtains Pina's equation (27) \\
\begin{tcolorbox}[ width=\textwidth, colback={white}]
\begin{align}
\widetilde{W}_{ij,ab}&( \Oms ) = (ij|ab) + 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
\times & \left[ \frac{ 1 }{ ( \Oms - \Omega_m^{RPA} ) - ( \varepsilon_b - \varepsilon_i) + i\eta } + \frac{ 1}{ ( \Oms - \Omega_m^{RPA} ) - (\varepsilon_a - \varepsilon_j) + i\eta }
\right] \nonumber
\end{align}
\end{tcolorbox}
\noindent This dynamical $\widetilde{W}( \Oms )$ kernel can be rewritten:
\begin{align}
\widetilde{W}_{ij,ab}( \Oms ) &= (ij|ab) + 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
& \times \left[ \frac{ 1 }{ \Omega_{ib}^{s} - \Omega_m^{RPA} + i\eta } + \frac{ 1}{ \Omega_{ja}^{s} - \Omega_m^{RPA} + i\eta }
\right] \nonumber
\end{align}
with e.g. $ \Omega_{ib}^{s} = \Oms - ( \varepsilon_b - \varepsilon_i) $. This formulation is similar to the spectral representation of $W_{ij,ab}(\omega)$ but with
$( \omega \rightarrow \Omega_{ib}^{s} )$ and $( \Omega \rightarrow \Omega_{ja}^{s} )$.
\vskip 1cm
\noindent {\textbf {Non-resonant contributions.}} We now add the $B_{jb}^{s}$ small contributions to the $\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) | N,s \rangle $ weight:
\begin{align*}
\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) | N,s \rangle^{\text{corr}} &= - \Big( e^{ -i \Omega_s t^{34} } \Big) \sum_{bj} \phi_j(x_3) \phi_b^*(x_4) B_{jb}^{s} \times \nonumber \\
\times & \Big( \theta( \tau_{34} ) e^{- i ( \varepsilon_j - \hOms ) \tau_{34} }
+ \theta( -\tau_{34} ) e^{ - i ( \varepsilon_b + \hOms) \tau_{34} } \Big)
\end{align*}
where ``corr" stands for correction. As a result, the integral associated with the screened Coulomb potential must be corrected by:
\begin{align*}
\int d34 \; \phi_a^*(x_3) & \phi_i(x_4)
\left[ \theta( \tau_{34} ) e^{i (\varepsilon_i + {\Omega_s \over 2} ) \tau_{34} } + \theta( - \tau_{34} ) e^{i (\varepsilon_a - {\Omega_s \over 2} ) \tau_{34} } \right] \times \\
& \times
e^{ i \Omega_s^{+} t^{34} } W(3^+,4) \langle N | T [{\hat \psi}(3) {\hat \psi}^{\dagger}(4) ] | N,s \rangle \times \theta( t^{34}_m - t_2 ) \\
= - \int d34 \; \phi_a^*(x_3) & \phi_i(x_4) \theta( t^{34}_m - t_2 ) W(3^+,4) \sum_{bj} \phi_j(x_3) \phi_b^*(x_4) B_{jb}^{s} e^{ - \delta t^{34} } \times \\
& \times \left[ \theta( \tau_{34} ) e^{i ( \varepsilon_i - \varepsilon_j + \Oms ) \tau_{34} } + \theta( - \tau_{34} ) e^{i (\varepsilon_a - \varepsilon_b - \Oms ) \tau_{34} } \right] \\
= - ({ \delta}^{-1} ) \sum_{jb} B_{jb}^{s} & \int d\tau \; W_{ib,aj}( {\tau}^{+} ) \times \\
& \times \left[ \theta( \tau ) e^{i ( \varepsilon_i - \varepsilon_j + \Oms ) \tau } + \theta( - \tau ) e^{i (\varepsilon_a - \varepsilon_b - \Oms ) \tau } \right]
\end{align*}
The time integration yields now:
\begin{align*}
\widetilde{W}_{ib,aj}(\Oms) =& \int d\tau \; W_{ib,aj}( \tau^{+} )
\times \; \left[ \theta( \tau ) e^{i (\varepsilon_i - \varepsilon_j+ \Oms ) \tau } + \theta( - \tau ) e^{i (\varepsilon_a - \varepsilon_b - \Oms ) \tau } \right] \\
=& \int { d\omega \over 2\pi } W_{ib,aj}( \omega ) \int d\tau \; e^{-i \omega ( \tau + 0^+) }
\left[ \theta( \tau ) e^{i (\varepsilon_i - \varepsilon_j+ \Oms ) \tau } + \theta( - \tau ) e^{i (\varepsilon_a - \varepsilon_b - \Oms ) \tau } \right] \\
=& { i \over 2\pi} \int d\omega e^{-i \omega 0^{+}} \left[ \frac{ W_{ib,aj}( \omega ) }{ ( \Oms - \omega ) - ( \varepsilon_j - \varepsilon_i) + i\eta }
+ \frac{ W_{ib,aj}( \omega ) }{ ( \Oms + \omega) - ( \varepsilon_a - \varepsilon_b ) + i\eta } \right]
\end{align*} \\
Plugging now the spectral representation of $W_{ib,aj}( \omega )$:
\begin{align*}
\widetilde{W}_{ib,aj}( \Oms ) = & (ib|aj) + 2 \sum_m^{OV} [ib|m] [aj|m] \times \\
\times & \Big( \frac{ 1 }{ ( \Oms - \Omega_m^{RPA} ) - ( \varepsilon_j - \varepsilon_i) + i\eta } + \frac{ 1}{ ( \Oms - \Omega_m^{RPA} ) - (\varepsilon_a - \varepsilon_b) + i\eta } \\
+ & \frac{ 1 }{ ( \Oms + \Omega_m^{RPA} ) - ( \varepsilon_a - \varepsilon_b ) + i\eta } - \frac{ 1}{ ( \Oms +\Omega_m^{RPA} ) - (\varepsilon_a - \varepsilon_b) - i\eta } \Big)
\end{align*}
where again the last line contribution cancel. \\
\noindent The correction to the bare Coulomb contribution stems from the correction to
$\langle N | T {\hat \psi}(5) {\hat \psi}^{\dagger}(5) | N,s \rangle$, namely :
$$
\langle N | T {\hat \psi}(5) {\hat \psi}^{\dagger}(5) | N,s \rangle^{\text{corr}} = - \Big( e^{ -i \Omega_s t_5 } \Big) \sum_{bj} \phi_j(x_5) \phi_b^*(x_5) B_{jb}^{s}
$$
leading to the following correction :
\begin{align*}
- \int d35 \; \phi_a^*(x_3) & \phi_i(x_3) e^{ i \Omega_s^{+} t_3 } v(3,5) \langle N | T [{\hat \psi}(5) {\hat \psi}^{\dagger}(5) ] | N,s \rangle \times \theta( t_5 - t_2 ) \\
& = \sum_{bj} B_{bj}^{s} v_{ia,bj} ( \delta^{-1} )
\end{align*}
As discussed above, there is no correction to the diagonal part of the BSE Hamiltonian so that one obtains beyond TDA so that: \\
\begin{tcolorbox}[ width=\textwidth, colback={white}]
\begin{align}
( \varepsilon_a - \varepsilon_i - \Omega_s ) A_{ia}^{s}
+& \sum_{jb} \Big( v_{ai,bj} - \widetilde{W}_{ij,ab}(\Oms) \Big) A_{jb}^{s} \nonumber \\
+ & \sum_{bj}\Big( v_{ai,jb} - \widetilde{W}_{ib,aj}(\Oms) \Big) B_{jb}^{s} = 0
\end{align}
\end{tcolorbox}
\noindent with:
\begin{align*}
\widetilde{W}_{ib,aj}( \Oms ) = & (ib|aj) + 2 \sum_m^{OV} [ib|m] [aj|m] \times \\
\times & \left[ \frac{ 1 }{ ( \Oms - \Omega_m^{RPA} ) - ( \varepsilon_j - \varepsilon_i) + i\eta } + \frac{ 1}{ ( \Oms - \Omega_m^{RPA} ) - (\varepsilon_a - \varepsilon_b) + i\eta } \right]
\end{align*}
%%%%%%%%%%
\end{document}

132
main.tex
View File

@ -1,132 +0,0 @@
\documentclass[journal=jctcce,manuscript=article]{achemso}
%%\usepackage{natbib}
%%\usepackage{notoccite}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{graphicx}% Include figure files:q
\usepackage{dcolumn}% Align table columns on decimal point
\usepackage{bm}% bold math
\usepackage{longtable}
\usepackage[usenames, dvipsnames]{color}
\usepackage{subfigure}
\usepackage{multirow}
\usepackage{dsfont}
\usepackage{ulem}
\usepackage{xcolor}
\DeclareMathOperator*{\argmax}{argmax}
\DeclareMathOperator*{\argmin}{argmin}
\newcommand\vari{{\varepsilon}_i}
\newcommand\vara{{\varepsilon}_a}
\newcommand\varj{{\varepsilon}_j}
\newcommand\varb{{\varepsilon}_b}
\newcommand\varn{{\varepsilon}_n}
\newcommand\varm{{\varepsilon}_m}
\newcommand\Oms{{\Omega}_s}
\newcommand\hOms{\frac{{\Omega}_s}{2}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\newcommand{\correction}[1]{\textcolor{blue}{#1}}
%\definecolor{myblue}{rgb}{0.0, 0.0, 1.0}
\author{Pierre-Fran\c{c}ois Loos}
\email{loos@irsamc.ups-tlse.fr}
\affiliation[LCPQ, Toulouse]{Laboratoire de Chimie et Physique Quantiques, Universit\'e de Toulouse, CNRS, UPS, France}
\author{Xavier Blase}
\email{xavier.blase@neel.cnrs.fr}
\affiliation[NEEL, Grenoble]{Universit\'e Grenoble Alpes, CNRS, Institut NEEL, F-38042 Grenoble, France}
\title{ Dynamical BSE }
\date{\today}
\begin{document}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{abstract}
\end{abstract}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Theory}
The Fourier components with respect to time $t_1$ of $iL_0(1, 4; 1', 3) = G(1, 3)G(4, 1')$ reads, dropping the (space/spin)-variables:
\begin{align*}
[iL_0]( \omega_1 ) = \frac{ 1 }{ 2\pi } \int d \omega \; G(\omega - \frac{\omega_1}{2} ) G( {\omega} + \frac{\omega_1}{2} ) e^{ i \omega \tau_{34} } e^{ i \omega_1 t^{34} }
\end{align*}
with $\tau_{34} = t_3 - t_4$ and
$t^{34} = (t_3 + t_4)/2$. Plugging now the 1-body Green's function Lehman representation, e.g.
$$
G(x_1,x_3 ; \omega) = \sum_n \frac{ \phi_n(x_1) \phi_n^*(x_3) } { \omega - \varepsilon_n + i \eta \text{sgn}(\varepsilon_n - \mu) }
$$
and projecting on $\phi_a^*(x_1) \phi_i(x_{1'})$, one obtains the $\omega_1= \Oms$ component
\begin{align*}
\int dx_1 dx_{1'} \; & \phi_a^*(x_1) \phi_i(x_{1'}) L_0(x_1,3;x_{1'},4; \Oms) = e^{i \Oms t^{34} } \times \\
& \frac{ \phi_a^*(x_3) \phi_i(x_4) } { \Oms - ( \vara - \vari ) + i \eta }
\Big( \theta( \tau ) e^{i ( \vari + \hOms) \tau }
+ \theta( - \tau ) e^{i (\vara - \hOms \tau } \Big)
\end{align*}
with $\tau = \tau_{34}$.
We further obtain the spectral representation of
$\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) | N,s \rangle$
expanding the field operators over a complete orbital basis creation/destruction operators:
\begin{align*}
\langle N | T {\hat \psi}(3) {\hat \psi}^{\dagger}(4) & | N,s \rangle = - \Big( e^{ -i \Omega_s t^{34} } \Big) \sum_{mn} \phi_m(x_3) \phi_n^*(x_4) \langle N | {\hat a}_n^{\dagger} {\hat a}_m | N,s \rangle \times \nonumber \\
\times & \Big( \theta( \tau ) e^{- i ( \varepsilon_m - \hOms ) \tau }
+ \theta( -\tau ) e^{ - i ( \varepsilon_n + \hOms) \tau } \Big)
\end{align*}
with $\tau = \tau_{34}$ and where the $ \lbrace \varepsilon_{n/m} \rbrace$ are proper addition/removal energies such that e.g.
$$
e^{ i H \tau } {\hat a}_m^{\dagger} | N \rangle = e^{ i (E_0^N + \varepsilon_m ) \tau } {\hat a} _m^{\dagger} | N \rangle
$$
Selecting (n,m)=(j,b) yields the largest components
$A_{jb}^{s} = \langle N | {\hat a}_j^{\dagger} {\hat a}_b | N,s \rangle $, while (n,m)=(b,j) yields much weaker
$B_{jb}^{s} = \langle N | {\hat a}_b^{\dagger} {\hat a}_j | N,s \rangle $ contributions. We used chemist notations with (i,j) indexing occupied orbitals and (a,b) virtual ones. Neglecting the $B_{jb}^{s}$ leads to the Tamm Dancoff approximation (TDA). Obtaining similarly the spectral representation of $ \langle N | T {\hat \psi}(1) {\hat \psi}^{\dagger}(1') | N,s \rangle$ ($t_{1'} = t_1^{+}$) projected onto $\phi_a^*(x_1) \phi_i(x_{1'})$,
one obtains after a few tedious manipulations (see Supplemental Information) the dynamical Bethe-Salpeter equation (DBSE) :
\begin{align}
( \varepsilon_a - \varepsilon_i - \Omega_s ) A_{ia}^{s}
&+ \sum_{jb} \Big( v_{ai,bj} - \widetilde{W}_{ij,ab}(\Oms) \Big) A_{jb}^{s} \\
&+ \sum_{bj} \Big( v_{ai,jb} - \widetilde{W}_{ib,aj}(\Oms) \Big) B_{jb}^{s}
= 0
\end{align}
with an effective dynamically screened Coulomb potential (see Pina eq. 24):
\begin{align}
\widetilde{W}_{ij,ab}(\Oms) &= { i \over 2 \pi} \int d\omega \; e^{-i \omega 0^+ } W_{ij,ab}(\omega) \times \\
\hskip 1cm &\times \left[ \frac{1}{ (\Oms - \omega) - ( \varb - \vari ) +i \eta } + \frac{1}{ (\Oms + \omega) - ( \vara - \varj ) + i\eta } \right] \nonumber
\end{align}
In the present study, we use the exact spectral representation of $W(\omega)$ at the RPA level:
\begin{align*}
W_{ij,ab}(\omega) &= (ij|ab) + 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
& \times \Big( \frac{1}{ \omega-\Omega_m^{RPA} + i\eta } - \frac{1}{ \omega + \Omega_m^{RPA} - i\eta } \Big)
\end{align*}
so that
\begin{align}
\widetilde{W}_{ij,ab}( \Oms ) &= (ij|ab) + 2 \sum_m^{OV} [ij|m] [ab|m] \times \\
& \times \left[ \frac{ 1 }{ \Omega_{ib}^{s} - \Omega_m^{RPA} + i\eta } + \frac{ 1}{ \Omega_{ja}^{s} - \Omega_m^{RPA} + i\eta }
\right] \nonumber
\end{align}
with e.g. $ \Omega_{ib}^{s} = \Oms - ( \varepsilon_b - \varepsilon_i) $. \textcolor{red}{Due to excitonic effects, the lowest BSE ${\Omega}_1$ excitation energy stands lower than the lowest $\Omega_m^{RPA}$ excitation energy, so that
e.g. $( \Omega_{ib}^{s} - \Omega_m^{RPA} )$ is strictly negative and cannot diverge. Further, $\Omega_{ib}^{s}$ and $\Omega_{ja}^{s}$ are necessarily negative for in-gap low lying BSE excitations, such that
$$
\left[ \frac{ 1 }{ \Omega_{ib}^{s} - \Omega_m^{RPA} + i\eta } + \frac{ 1}{ \Omega_{ja}^{s} - \Omega_m^{RPA} + i\eta }
\right]
<
\Big( \frac{1}{ \omega-\Omega_m^{RPA} + i\eta } - \frac{1}{ \omega + \Omega_m^{RPA} - i\eta } \Big) < 0
$$
in the limit $(\omega \rightarrow 0)$ of the standard adiabatic BSE . WELL, do we know the sign of
$[ij|m] [ab|m]$ ?? }
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\end{document}