note-name = , use-sort-key = false

# Exploring the Statically Screened $G3W2$ Correction to the $GW$ Self-Energy:
Charged Excitations and Total Energies of Finite Systems

Arno Förster<EMAIL_ADDRESS>Theoretical Chemistry, Vrije Universiteit,
De Boelelaan 1083, NL-1081 HV, Amsterdam, The Netherlands  Lucas Visscher
Theoretical Chemistry, Vrije Universiteit, De Boelelaan 1083, NL-1081 HV,
Amsterdam, The Netherlands

###### Abstract

Electron correlation in finite and extended systems is often described in an
effective single-particle framework within the $GW$ approximation. Here, we
use the statically screened second-order exchange (SOX) contribution to the
self-energy ($G3W2$) to calculate a perturbative correction to the $GW$ self-
energy. We use this correction to calculate total correlation energies of
atoms, relative energies, as well as charged excitations of a wide range of
molecular systems. We show that the second-order correction improves
correlation energies with respect to the RPA and also improves relative
energies for many, but not all considered systems. The dynamically screened
SOX term has previously been shown to consistently lower the highest occupied
molecular orbital (HOMO) quasiparticle (QP) energies and to increase the
lowest unoccupied molecular orbitals (LUMO) QP energies. We show here that the
statically screened $G3W2$ correction consistently increases the LUMO QP
energies, while no consistent trend can be observed for the HOMO levels. Also,
confirming previous results, the magnitude of the correction is much smaller
with the statically screened interaction, than with the dynamically screened
one. qsGW by itself is shown to be an excellent method for the calculation of
charged excitation of finite systems, and it can not consistently be improved
upon by the $G3W2$ correction. For range-separated hybrid starting points, the
description of fundamental gaps and HOMO QP energies is slightly worsened.
However, tremendous improvements upon the $GW$ LUMO energies leading to almost
perfect agreement with high-level coupled cluster reference methods are
observed. The evaluation of the statically screened $G3W2$ correction only
comes with small additional computational cost compared to $G_{0}W_{0}$ for
systems with up to 100 atoms and should therefore be suitable for practical
applications.

Perturbation Theory, $GW$, Quasiparticle

###### pacs:

##  I Introduction

Electron densities, total energies and spectral functions of an interacting
many-electron system can be extracted from its single-particle Green’s
function.Abrikosov  _et al._ (1975); Mattuck (1992); Martin  _et al._ (2016)
Instead of expanding the Green’s function directly in powers of the electron-
electron interaction, it is usually easier to expand the irreducible self-
energy and to calculate the interacting Green’s function via Dyson’s
equation.Dyson (1949) The difficulty is then to find an approximation which
captures the most important correlation effects and remains computationally
tractable.

In Hedin’s $GW$ approximation (GWA),Hedin (1965) the self-energy is obtained
as the first term of an expansion in terms of a screened electron-electron
interaction,Phillips (1961) where screening is (usually) calculated in the
random-phase approximation (RPA).Hubbard (1957) Compared to the Hartree-Fock
(HF) approximation, the first-order term in the expansion of the self-energy
in terms of the bare electron-electron interaction, it takes into account that
the electron-electron interaction at large distances is screened by the
presence of other electrons.Martin  _et al._ (2016) For properties which are
mainly dominated by long-range correlation effects, such as fundamental gaps
or ionization potentials in organic molecules, the GWA is therefore often very
accurate.Reining (2018); Golze  _et al._ (2019) This is especially true for
the most common computational approach, where the full solution of Dyson’s
equation is bypassed and the quasiparticle (QP) energies are calculated from a
non-interacting propagator.Hybertsen and Louie (1986) In the static long-range
limit, the errors introduced by the neglect of self-consistency and the
absence of higher order terms in the polarizability and self-energy cancel to
a large extentKotani  _et al._ (2007) by virtue of the Ward identity.Ward
(1950) When short-range correlation becomes important, the GWA is less
accurate.

From the $GW$ self-energy one can calculate correlation energies using the
Klein functional.Klein (1961) Evaluated with a Kohn-Sham (KS) Green’s
function, one obtains (particle-hole) RPA,Casida (1995); Dahlen  _et al._
(2006) often derived in the framework of the adiabatic connection (AC)
fluctuation-dissipation theoremRen  _et al._ (2012) or as a subset of terms in
the coupled cluster (CC)Coester (1958); Coester and Kümmel (1960); Čížek
(1966, 1969); Paldus  _et al._ (1972) singles and doubles (CCD)
expansion.Scuseria  _et al._ (2008, 2013) Since the effect of charge screening
dominates in this limit, the RPA describes long-range electron correlation
very accuratelyLangreth and Perdew (1977). At short electron-electron
distances, charge screening is less important and the RPA drastically
overestimates the magnitude of electron correlation.Singwi  _et al._ (1968) In
this limit it is crucial to take into account higher-order contributions
(vertex corrections) to the self-energy which can then dominate.Irmler  _et
al._ (2019) These contributions are fundamentally short-ranged and become less
and less important for large electron-electron distances.Jiang and Engel
(2007); Ren  _et al._ (2012)

Besides range-separation based approaches,Kurth  _et al._ (1999); Yan  _et
al._ (2000); Angyan  _et al._ (2005); Jiang and Engel (2007); Janesko  _et
al._ (2009a, b); Toulouse  _et al._ (2009); Zhu  _et al._ (2010); Toulouse
_et al._ (2010) the next-to-leading-order term, the second-order exchange
(SOX), is therefore the most obvious starting point to treat electron
correlation beyond the RPA. The most common correction is the second-order
screened exchange (SOSEX) correction.Freeman (1977) SOSEX can be derived from
the CC doubles (CCD) equations by neglecting ladder and ”crossed-ring”Scuseria
_et al._ (2013) diagrams.Scuseria  _et al._ (2008, 2013) SOSEX has been
applied to total energies of atoms, and molecules,Paier  _et al._ (2012); Ren
_et al._ (2013) solidsGrüneis  _et al._ (2009) and to the uniform electron
gasHummel  _et al._ (2019) and several extensions to it have been proposed as
well.Paier  _et al._ (2008); Toulouse  _et al._ (2011); Ángyán  _et al._
(2011); Hummel  _et al._ (2019) In another variant the SOX term is
renormalized by resummation of hole-hole ladder diagrams.Engel and Jiang
(2006); Jiang and Engel (2006) Typically, these approaches improve total
correlation energiesJiang and Engel (2007); Gould  _et al._ (2019) but
relative energies are not necessarily improved.Ren  _et al._ (2012)

The close connection between RPA correlation energies and the $GW$ method
suggests, that beyond-RPA methods might also be used to improve over the $GW$
method for charged excitations. For example, Ren et al. applied SOSEX to
moleculesRen  _et al._ (2015) and Maggio and KresseMaggio and Kresse (2017)
used a similar method to calculate HOMO QP energies for small molecules. It
is, however, more common to formulate beyond-$GW$ schemes by taking into
account vertex corrections in Hedin’s equations.Hedin (1965) Kutepov has
solved Hedin’s equations in a fully self-consistent fashion including several
vertex corrections.Kutepov (2016) He has applied the second-order self-energy
variant (coined $GW+G3W2$) to a wide range of metals, semiconductors, and
insulators,Kutepov and Kotliar (2017); Kutepov (2017, 2018, 2021a, 2021b,
2022) and observed major improvements over fully self-consistent $GW$ for
spectral properties using the $GW+G3W2$ approximation to the self-energy.
According to his work, the vertex corrected calculations should be performed
without any constraining approximations. This means, Hedin’s equations should
be solved fully self-consistently and the frequency dependence of all
quantities should be properly accounted for. As many other authors,Onida  _et
al._ (1995, 2002); Adragna  _et al._ (2003); Marini  _et al._ (2003); Bruneval
_et al._ (2005); Shishkin  _et al._ (2007); Romaniello  _et al._ (2009);
Sharma  _et al._ (2011); Chen and Pasquarello (2015); Hung  _et al._ (2016);
Cunningham  _et al._ (2018); Olsen  _et al._ (2019); Schmidt  _et al._ (2017);
Tal  _et al._ (2021); Cunningham  _et al._ (2021) he also emphasized the
importance of vertex corrections in the polarizabiliy.

On the other hand, the $GW+G3W2$ self-energy has also been applied by Grüneis
et al.Grüneis  _et al._ (2014) to calculate QP energies in solids in a more
approximate fashion. In contrast to Kutepov, Grüneis et al.Grüneis  _et al._
(2014) used a statically screened interaction and calculated perturbative
corrections to quasiparticle self-consistent $GW$ (qsGW)van Schilfgaarde  _et
al._ (2006); Kotani  _et al._ (2007) QP energies. They also combined this
method with a beyond RPA screening of the Coulomb interaction using a static
exchange-correlation kernel. As Kutepov, they reported major improvements over
the $GW$ method for the band structures of solids. Comparing their results to
Kutepov’s, the logical conclusion would be, that the errors introduced by the
static approximation, the use of a non-interacting Green’s function, and the
perturbative treatment of the $G3W2$ cancel to a large extent. A recent study
by Kutepov has confirmed this conjecture.Kutepov (2022) From a pragmatic point
of view this is of course convenient, since all of these approximations come
with drastically reduced computational cost compared to the rigorous, self-
consistent formalismKutepov (2016). In a recent study, the $GW+G3W2$ self-
energy has also been employed by Rinke, Ren and coworkersWang  _et al._ (2021)
to ionization potentials and electron affinities of molecules. Rinke, Ren and
coworkers used the dynamically (RPA) screened interaction and did not consider
any self-consistency in solving Dyson’s equation. Improvements over $GW$ were
found to be substantial, especially for electron affinities and the dependence
on the starting point was reduced.Wang  _et al._ (2021)

Encouraged by these results, we herein follow Grüneis et al.Grüneis  _et al._
(2014) and calculate perturbative $G3W2$ corrections to the $GW$ self-energy.
We follow their approach in using non-interacting Green’s functions throughout
and exclusively consider the statically screened $G3W2$ self-energy. However,
following Rinke, Ren and coworkersWang  _et al._ (2021), we do not consider
any beyond-RPA screening in our work. One should expect an expansion of the
polarizability in terms of $W$ to converge much faster for finite, than for
extended systems, since charge screening will be much weaker.Van Den Brink and
Sawatzky (2000) It is known for finite systems, that inclusion of the vertex
in the polarizability alone does not improveVerdozzi  _et al._ (1995);
Schindlmayr and Godby (1998) and might even deteriorate the accuracy of
charged excitations compared to $GW$.Lewis and Berkelbach (2019) While Grüneis
et al.Grüneis  _et al._ (2014) investigated periodic systems, we focus here on
the calculation of correlation energies as well as QP energies for a wide
range of finite systems: We show that massively improved absolute correlation
energies are obtained compared to the RPA and also that relative energies are
improved for the majority of the considered systems. For QP energies we show,
that the $G3W2$ correction provides improvements over $G_{0}W_{0}$ evaluated
with KS orbitals and energies from range-separated hybrid functionals,
especially for electron affinities. However, no consistent improvements over
qsGW are found which we show to be very accurate by itself.

The remainder of this work is organized as follows: In the next section, we
introduce the $G3W2$ self-energy and discuss how it can be used efficiently to
calculate corrections to $GW$ QP energies and RPA correlation energies. In
section IV, we first benchmark the performance of the $GW+G3W2$ self-energy
for absolute correlation energies as well as relative energies of several
organic molecules and then proceed to scrutinize the effect of the $G3W2$
self-energy correction to $G_{0}W_{0}$ on ionization potentials, electron
affinities and fundamental gaps of molecular systems. Finally, we summarize
and conclude this work in section ref. V.

##  II Theory

Dyson’s equation,Dyson (1949)

$G(1,2)=G^{0}(1,2)+G^{0}(1,\overline{3})\Sigma(\overline{3},\overline{4})G(\overline{4},2)\;,$
(1)

relates the interacting Green’s function $G$ of a many-electron system to the
Green’s function of a non-interacting reference system $G^{0}$, via the non-
Hermitian self-energy operator $\Sigma$. We have used the notation
$1=(\bm{r}_{1},\sigma_{1},\omega_{1})$ and integration over indices with upper
bars is implied. Following HubbardHubbard (1957), PhillipsPhillips (1961) and
Hedin,Hedin (1965) the self-energy can be expanded in terms of the screened
electron-electron interaction. The first two terms in this expansion are

$\Sigma^{GW+G3W2}(1,2)=\Sigma^{GW}(1,2)+\Sigma^{G3W2}(1,2)\;,$ (2)

with

$\Sigma^{GW}(1,2)=iG(1,2)W(1,2)\;,$ (3)

and

$\Sigma^{G3W2}(1,2)=-G(1,\overline{4})W(1,\overline{3})G(\overline{4},\overline{3})G(\overline{3},2)W(\overline{4},2)\;.$
(4)

$W$ is the screened Coulomb interaction which can be obtained from the
polarizability $P$ and the bare Coulomb interaction $V$ via a Dyson equation

$W(1,2)=V(1,2)+V(1,\overline{3})P(\overline{3},\overline{4})W(\overline{4},2)\;.$
(5)

As $\Sigma$, $P$ can be expanded in powers of the screened interaction where
to zeroth order in $W$, the RPA is obtained,

$P^{RPA}(1,2)=-iG(1,2)G(2,1)\;.$ (6)

Solving these Equations fully self-consistently is cumbersome due to the
complicated frequency dependence of the interacting Green’s function and the
second-order self-energy. Since we are interested in a computationally
efficient approach, we introduce a few approximations: First, all quantities
above are exclusively evaluated with non-interacting Green’s functions, which
implies that we need to approximate (1). Second, $\Sigma^{G3W2}$ is evaluated
perturbatively and the static limit for $W$ is used. We will discuss this term
in more detail in the next subsection. First, we outline our strategy on the
GWA level.

###  II.1 Partial self-consistency in the GWA

In a discrete basis of canonical molecular orbitals (MO) (we will label
general MOs by $p,q,r,\dots$, occupied orbitals by $i,j,k,\dots$, and virtual
orbitals by $a,b,c,\dots$), eq. (1) can be written in the form (see for
example ref. Hirata  _et al._ (2017))

$\sum_{r}\left\\{\epsilon_{pr}-\left[v_{xc}\right]_{pr}+\left[\Sigma_{xc}\right]_{pr}(\omega_{p})\right\\}U_{rq}(\omega_{p})=\omega_{p}U_{pq}(\omega_{p})\;.$
(7)

Here, $\epsilon$ is the Hamiltonian of a non-interacting reference system,
typically calculated within the HF approximation or KS-DFT, and is diagonal in
the MO basis. If the off-diagonal elements of $\Sigma_{xc}-v_{xc}$ are
neglected, (7) reduces to a set of independent, non-linear equations,

$\epsilon_{p}-\left[v_{xc}\right]_{pp}+\left[\Sigma_{xc}\right]_{pp}(\omega_{p})=\omega_{p}\;,$
(8)

for the QP corrections $\omega_{p}$ to the single-electron energy levels
$\epsilon_{p}$. If $\Sigma$ is calculated in the GWA, this approach is
referred to as $G_{0}W_{0}$. $G_{0}W_{0}$ can give very accurate results but
this heavily depends on the choice of the underlying exchange-correlation
functional.Bruneval  _et al._ (2006a, b); Bruneval and Marques (2013); Caruso
_et al._ (2016); Knight  _et al._ (2016) If the eigenvalues of the KS
Hamiltonian are too bad of an approximation to the QP energies (for typical
GGAs or the LDA, errors are typical of the order of a few eV), the
perturbative approach will fail.Bruneval  _et al._ (2006b); Caruso  _et al._
(2016) Therefore, it is mandatory to calculate the $G_{0}W_{0}$ correction
from a potential which is already close to the $GW$ self-energy.Atalla  _et
al._ (2013); Hellgren  _et al._ (2021)

In this work, we pursue two closely related and well-known strategies for
this: The first one is to perform a $G_{0}W_{0}$ calculation using KS
eigenvalues and orbitals from a range-separated hybrid calculation as input.
This is usually an excellent starting point since the screened Coulomb
potential in the RPA reduces to

$W(r\rightarrow\infty,\omega=0)=\frac{e^{-\lambda r}}{r}\;,$ (9)

with $\lambda$ proportional to the Fermi-momentum.Mattuck (1992) This is the
basic idea behind the construction of range-separated hybrid (RSH)
functionals. Solving the KS equations using such a potential emphasizes the
importance of screening of the Coulomb interaction for large electron-electron
distances and a perturbative correction with eq. (8) is usually rather
small.Kanchanakungwankul and Truhlar (2021) The range-separation parameter can
then also be adjusted in a way that the eigenvalues from the RSH calculation
gives very accurate QP energies by itself.Refaely-Abramson  _et al._ (2011,
2012); Kronik  _et al._ (2012)

The second approach we employ here is QP self-consistent $GW$ (qs$GW$):Kotani
_et al._ (2007); van Schilfgaarde  _et al._ (2006) For a given approximation
to the self-energy an optimized KS potential can be obtained by solving the
(linearized)Godby  _et al._ (1988) Sham-Schlüter equations (SSE)Sham and
Schl̈üter (1983); Casida (1995). With $\Sigma=GW$, the SSE can be solved as in
the exchange-only case if the self-energy is approximated as static.Sharp and
Horton (1953); Talman and Shadwick (1976); Engel, Eberhard and Dreizler (2013)
One then obtains a statically screened exchange (SEX) correction depending on
$Re\Sigma(\epsilon^{QP})$. The derivation via the SSE leaves some ambiguity in
the choice of the off-diagonal elements. As shown by Ismail-Beigi,Ismail-Beigi
(2017) one can also derive an optimized potential via minimization of the
gradient of the Klein functional,Klein (1961) leading to the following
Hermitian form which has already been suggested by van Schilfgaarde et al.van
Schilfgaarde  _et al._ (2006),

$\left[v^{OP}_{xc}\right]_{pq}=\frac{1}{2}\left\\{\text{Re}\left[\Sigma_{xc}\right]_{pq}(\epsilon_{p})+\text{Re}\left[\Sigma_{xc}\right]_{pq}(\epsilon_{q})\right\\}\;.$
(10)

Using this potential, the KS equations are then solved self-consistently. It
should be stressed that in fully self-consistent GW, an interacting Green’s
function is constructed in each iteration by solving (1) while in qs$GW$ a
non-interacting Green’s function is constructed in each iteration. qs$GW$
should thus rather be seen as an approach to optimize a non-interacting
Green’s function. It can also be combined with (8) ($G_{0}W_{0}@$qsGW),Klimeš
and Kresse (2014) However, this is typically not done in practice: In
situations where the single QP picture is valid, the self-energy will be
varying slowly around the QP positionPeter Kopietz (1997) and the effect of
the $G_{0}W_{0}$ correction will vanish.Lee  _et al._ (2020) qsGW is also
closely related to the RPA optimized effective potential (OEP) method when the
self-energy is approximated as static in the QP approximation.Klimeš and
Kresse (2014); Riemelmoser  _et al._ (2021) Note, however, that only the self-
energy is approximated as static, but not the screened interaction, as in
Hedin’s Coulomb hole (COH)SEX approximation.Hedin (1965)

###  II.2 Evaluation of the Self-Energy

We now proceed with the discussion of the $G3W2$ self-energy. To this end, we
first introduce the non-interacting Green’s function. Since we exclusively
work with non-interacting Green’s functions we will suppress the superscript
$(0)$ in the following. The time-ordered free propagator $G$ is diagonal in
the MO basis and can be expressed in terms of greater $G^{>}$ (particle) and
lesser $G^{<}$ (hole) propagators,

$G_{p}(t_{1},t_{2})=\Theta(t_{1}-t_{2})G_{p}^{>}(t_{1},t_{2})+\Theta(t_{2}-t_{1})G_{p}^{<}(t_{1},t_{2})\;,$
(11)

where

$\displaystyle G_{p}^{>}(t_{1},t_{2})=$ $\displaystyle
G_{p}^{>}(t_{1}-t_{2})=-i\overline{f}(\epsilon_{p})e^{-i\epsilon_{p}(t_{1}-t_{2})}$
(12) $\displaystyle G_{p}^{<}(t_{1},t_{2})=$ $\displaystyle
G_{p}^{<}(t_{1}-t_{2})=if(\epsilon_{p})e^{-i\epsilon_{p}(t_{1}-t_{2})}\;.$
(13)

Here, $f$ denotes the zero-temperature Fermi-function and $\overline{f}=1-f$,
and we have indicated that the Green’s function only depends on the time-
difference at equilibrium. The single-particle energies in (13) and (12) are
understood to be relative to the chemical potential which we place in the
middle of the HOMO-LUMO gap.

Figure 1: Diagrammatic representation of the lesser components of the $GW$
(top) and $G3W2$ (bottom) self-energy parts. Time-ordering from left to right
is implied. The signs at the vertices denote which branches of the Keldysh
contour are connected by $W$. Greater components are obtained by exchanging
plus and minus signs. Straight lines denote propagators and the double wiggly
line is the screened Coulomb interaction.

The structure of the self-energy terms has been scrutinized by Stefanucci,
Pavlyukh, van Leeuwen and coworkers.Stefanucci  _et al._ (2014); Pavlyukh  _et
al._ (2020) We follow their work and use the framework of the Keldysh
formalism to discuss the $G3W2$ self-energy.van Leeuwen  _et al._ (2015) In
the Keldysh formalism, we work on the contour
$\mathcal{C}=\mathcal{C}^{+}\cup\mathcal{C}^{-}$, with $\mathcal{C}^{+}$ being
the backward branch and $\mathcal{C}^{-}$ being the forward branch. The time-
ordered (anti-time-ordered) Green’s function $G^{--}$ ($G^{++}$) is built from
field-operators evolving on $\mathcal{C}^{-}$ ($\mathcal{C}^{+}$) while lesser
($G^{-+}=G^{<}$) and greater $G^{+-}=G^{>}$ Green’s function involve both
branches and describe propagation of holes and particles, respectively. In the
same way, the dynamically screened interaction can either connect both
different branches on the Keldysh contour ($W^{+-}=W^{>}$ and $W^{-+}=W^{<}$)
or not ($W^{--}$ and $W^{++}$).

The lesser components of both contributions to the self-energy in (2) are
given diagrammatically in figure 1. For the $GW$ self-energy each component
only consists of a single term (upper part of figure 1). For finite systems,
the $GW$ self-energy is most conveniently evaluated in a basis of localized
atomic orbitals and imaginary time,

$\Sigma^{GW,\lessgtr}_{\mu\nu}(i\tau)=-i\sum_{\mu^{\prime}\nu^{\prime}}G^{\lessgtr}_{\mu^{\prime}\nu^{\prime}}(i\tau)W^{\lessgtr}_{\mu\mu^{\prime}\nu\nu^{\prime}}(i\tau)\;,$
(14)

with $i\tau=i\tau_{1}-i\tau_{2}$. In our $GW$ implementation,Förster and
Visscher (2020) we transform the self-energy to the imaginary frequency axis
and analytically continue (either the diagonal elements in the MO basis in
$G_{0}W_{0}$ or the full self-energy matrix in qsGW) to the complex plane.
Greater and lesser component of the second-order contribution consist of four
terms each, since both intermediate vertices can connect both branches of the
Keldysh contour. Together, the four diagrams describe three distinct
scattering processes.Pavlyukh  _et al._ (2016, 2020) Among others, they are
responsible for spectral features which do not appear in fully self-consistent
$GW$, like the excitation of two plasmons and two particle hole pairs. In this
work, these terms should be of minor relevance only since we are interested in
improving QP energies and correlation energies. The last diagram in fig. 1,
however, only contributes to the $2h1p$ ($1h2p$) space and describes the
exchange of two final particles/holes.Pavlyukh  _et al._ (2020) In the same
way as the static SOX term it ensures the antisymmetry of the 4-point vertex
to first order in the electron-electron interaction.

Static screening comes with two distinct advantages. Obviously, it drastically
reduces the computational cost. Second, only the last diagram in the
expression in figure 1 remains. As an important side effect, the loss of the
positive definiteness of the spectral function Stefanucci  _et al._ (2014)
arising from straightforward inclusion of all four terms is avoided. The first
three diagrams vanish, since a static interaction can only connect points on
the same branch of the Keldysh contour.van Leeuwen  _et al._ (2015) To see
this clearly we first write down $\Sigma^{G3W2}$ explicitly,

$\displaystyle\Sigma^{G3W2}_{\mu\nu}(\tau_{12})=$
$\displaystyle-\sum_{\begin{subarray}{c}\kappa\kappa^{\prime}\lambda\lambda^{\prime}\\\
\mu^{\prime}\nu^{\prime}\end{subarray}}\int
d\tau_{3}d\tau_{4}G_{\mu^{\prime}\lambda}(\tau_{14})W_{\mu\mu^{\prime}\kappa\kappa^{\prime}}(\tau_{13})$
(15) $\displaystyle\times
G_{\lambda^{\prime}\kappa}(\tau_{43})G_{\kappa^{\prime}\nu^{\prime}}(\tau_{32})W_{\lambda\lambda^{\prime}\nu\nu^{\prime}}(\tau_{42})\;.$

Using a time-independent interaction, the integrals over the internal times
can be evaluated easily and transformation to frequency space using

$F(i\omega)=-i\int d\tau F(i\tau)e^{-i\omega\tau}$ (16)

gives

$\displaystyle\Sigma^{G3W2}_{\mu\nu}(i\omega)=$ $\displaystyle
i\sum_{\begin{subarray}{c}\kappa\kappa^{\prime}\lambda\lambda^{\prime}\\\
\mu^{\prime}\nu^{\prime}\end{subarray}}\int
d\tau_{12}e^{-i\omega\tau_{12}}G_{\mu^{\prime}\lambda}(\tau_{12})W_{\mu\mu^{\prime}\kappa\kappa^{\prime}}$
(17) $\displaystyle\times
G_{\lambda^{\prime}\kappa}(\tau_{21})G_{\kappa^{\prime}\nu^{\prime}}(\tau_{12})W_{\lambda\lambda^{\prime}\nu\nu^{\prime}}\;,$

which corresponds to the last diagram in figure 1 and can be expressed in
terms of greater and lesser components of $G$ only (plus and minus signs
alternate). Using (11), the integral has a simple analytic solution and
transforming to the MO basis gives

$\displaystyle\left[\Sigma^{G3W2}\right]_{pq}(\omega)=$
$\displaystyle\sum_{iab}\frac{W_{iapb}W_{ibqa}}{\epsilon_{a}+\epsilon_{b}-\epsilon_{i}-\omega}$
(18)
$\displaystyle-\sum_{ija}\frac{W_{iajp}W_{iqja}}{\epsilon_{a}-\epsilon_{i}-\epsilon_{j}+\omega}\;,$

which is the same expression already employed by Grüneis et al.Grüneis  _et
al._ (2014) The computational effort to evaluate these terms only grows as
$N^{3}$ with systems size if one is only interested in a small number of
diagonal elements. Calculation of the matrix elements of $W$ will generally
scale as $N^{4}$, but these will be available in some form when a $GW$
calculation is performed. The asymptotic scaling can be further reduced by
exploiting sparsity in the AO basis or in real space,Willow  _et al._ (2012,
2013); Mardirossian  _et al._ (2018) but we do not explore these possibilities
here.

For the self-energy beyond $GW$, we again rely on the assumption that $GW$
already gives rather accurate QP energies. We thus expand $\Sigma^{G3W2}$
around the $GW$ QP energies. At zeroth order, we obtain

$\epsilon^{GW+G3W2}_{p}=\epsilon^{GW}_{p}+\Sigma^{G3W2}(\epsilon^{GW}_{pp})\;.$
(19)

This is the most economical way to calculate the QP energy correction due to
$\Sigma^{G3W2}$. Note, that on the GWA level, the self-energy is often treated
to first order, which results in multiplication with a renormalization factor
$Z$. This is commonly done for solids, but it has been argued that the zeroth-
order treatment is to be preferred on the GWA level.Niquet and Gonze (2004) Of
course, evaluating (8) directly with the $\Sigma^{GW+G3W2}$ self-energy is
possible as well as has been done in ref. (73). Due to its much smaller
computational cost, all values presented in this work have been obtained using
equation (19).

###  II.3 Correlation Energies

For a given self-energy, total energies can be evaluated using the Klein
functional.Klein (1961) Following Dahlen et al.Dahlen  _et al._ (2006) we
obtain for the $\Sigma^{GW+G3W2}$ self-energy with a static $W$ in the second-
order term (Trace implies integration over spatial coordinates, spin and
frequency variables)

$\displaystyle E_{xc}=$ $\displaystyle
E_{x}[\phi^{KS}]+\frac{1}{2}\text{Tr}\left\\{\ln(1-P^{RPA}V)+P^{RPA}V\right\\}$
(20)
$\displaystyle+\frac{1}{4}\text{Tr}\left\\{G_{0}\Sigma^{G3W2}_{c}\right\\}\;.$

The first term is the HF exchange-correlation energy expression, evaluated
with KS orbitals, the second term is the usual expression for the RPA
correlation energyLangreth and Perdew (1977) and the third term is equivalent
to the SOX term in second-order Møller-Plesset perturbation theory (MP2) (see
for example eq. (A5) in ref. (15)), with the bare interaction lines replaces
by the statically screened interaction lines. When discussing correlation
energies in the following, we refer to the last term on the _r.h.s._ of (20)
as second-order statically screened Exchange (SOSSX). Correlation energies
evaluated with (20) is then called RPA+SOSSX, in analogy to RPA+SOSEX.

##  III Technical and Comptational Details

All expressions presented herein have been implemented in a locally modified
development version of the Amsterdam density functional (ADF) engine of the
Amsterdam modelling suite 2021.1 (AMS2021). The implementation of the GWA has
been outlined in previous work.Förster and Visscher (2020, 2021a, 2021b) In
all calculations, we expand 2-point correlation functions in correlation
consistent bases of Slater-type orbitals of triple- and quadruple-$\zeta$
quality (TZ3P and QZ6P, respectively).Förster and Visscher (2021a) Imaginary
time and imaginary frequency variables are discretized using non-uniform bases
$\mathcal{T}=\left\\{\tau_{\alpha}\right\\}_{\alpha=1,\dots N_{\tau}}$ and
$\mathcal{W}=\left\\{\omega_{\alpha}\right\\}_{\alpha=1,\dots N_{\omega}}$ of
sizes $N_{\tau}$ and $N_{\omega}$, respectively, tailored to each system. More
precisely, (16) is implemented as

$\displaystyle\overline{F}(i\omega_{\alpha})=$
$\displaystyle\Omega^{(c)}_{\alpha\beta}\overline{F}(i\tau_{\beta})$ (21)
$\displaystyle\underline{F}(i\omega_{\alpha})=$
$\displaystyle\Omega^{(s)}_{\alpha\beta}\underline{F}(i\tau_{\beta})\;,$ (22)

where $\overline{F}$ and $\underline{F}$ denote even and odd parts of $F$,
respectively. The transformation from imaginary frequency to imaginary time
only requires the (pseudo)inversion of $\Omega^{(c)}$ and $\Omega^{(s)}$,
respectively. Our procedure to calculate $\Omega^{(c)}$ and $\Omega^{(s)}$ as
well as $\mathcal{T}$ and $\mathcal{W}$ follows Kresse and coworkers.Kaltak
_et al._ (2014a, b); Liu  _et al._ (2016) The technical specifications of our
implementation have been described in the appendix of ref. (109). All 4-point
correlation functions are expressed in auxiliary basis sets of Slater type
functions which are usually 3 to 5 times larger than the primary bases.Förster
_et al._ (2020); Förster and Visscher (2020) The transformation between both
bases is implemented using the pair-atomic density fitting (PADF)
method.Baerends  _et al._ (1973) For an outline of the implementation of this
method, we refer to ref. (114).

We have implemented RPA and RPA+SOSSX correlation energies using (20). In
analogy to (18) and in the same way as in typical MP2 implementations in
quantum chemistry codes,Cooper  _et al._ (2017); Förster  _et al._ (2020) the
last term on the _r.h.s_ of (20) is evaluated directly in the MO basis. The
second term on the _r.h.s_ of (20) is evaluated in the space of auxiliary
basis functions, following our $GW$ implementation described in ref. (103).

###  III.1 Energies

The correlation energies can be calculated with different orbitals. In this
work, we use orbitals from HF, PBE0Ernzerhof and Scuseria (1999); Adamo and
Barone (1999), and the exact exchange (EXX) only OEPTalman and Shadwick
(1976); Engel, Eberhard and Dreizler (2013) in the implementation of Scuseria
and coworkersIzmaylov  _et al._ (2007) within the Krieger Lee Iafrate
(KLI)Krieger  _et al._ (1990) approximation.

####  III.1.1 Relative energies

For relative energies, we use imaginary frequency and time grids of 20 points.
Note, that 20 grid points are typically sufficient to converge total energies
to within $1e-6$ Hartree.Kaltak  _et al._ (2014a) Correlation energies are
calculated using,Helgaker  _et al._ (1997)

$E_{CBS}=E_{QZ}-E^{c}_{QZ}+\frac{E^{c}_{QZ}*4^{3}-E^{c}_{TZ}*3^{3}}{4^{3}-3^{3}}\;,$
(23)

where $E^{c}_{QZ}$ ($E^{c}_{TZ}$) denotes the correlation energy at the QZ6P
(TZ3P) level and $E_{QZ}$ is the total energy at the QZ6P level. The
extrapolation scheme has been shown to be suitable for correlation consistent
basis sets but can not be used for KS or HF contributions.Helgaker  _et al._
(1997); Jensen (2013) For the purpose of calculating relative energies, we
assume these contributions to be converged at the QZ6P level. We use the
_VeryGood_ numerical quality for integrals over real space and distance cut-
offs and the corresponding auxiliary basis set.Förster  _et al._ (2020)
Dependency thresholdsFörster and Visscher (2020) have been set to $5e^{-4}$.

####  III.1.2 Total Correlation Energies

We use 32 points for imaginary time and imaginary frequency each. We use here
the TZ3P and QZ6P basis sets with additional diffuse functions as described in
ref. (109) and additional tight 1s and 2p functions as will be described
below.See Supplemental Material at [URL will be inserted by publisher] for all
newly composed basis sets used in this work. Eq (23) is used for extrapolation
and numerical quality and auxiliary fit set quality are set to _excellent_. No
dependency thresholds have been set.

###  III.2 Quasiparticle Energies

Charged excitations in this work are calculated using (19). For qsGW, we set
the dependency threshold to $5e^{-3}$ and perform a maximum of 15 iterations
of the self-consistency cycle. Following the recommendations given in ref.
(110), we use _VeryGood_ numerical quality and the corresponding auxiliary
basis set. We perform $G_{0}W_{0}$ calculations using PBE, PBE0,
LRC-$\omega$PBEhRohrdanz  _et al._ (2009) and $\omega$B97-XChai and Head-
Gordon (2008) orbitals and eigenvalues. We use _Good_ numerical quality and
the corresponding auxiliary basis set for all $G_{0}W_{0}$ calculations and
set the dependency threshold to $5e^{-4}$. In ref. (103) these settings have
been shown to be appropriate. All QP energies are calculated using

$\epsilon_{n}^{CBS}=\epsilon_{n}^{QZ}-\frac{1}{N^{QZ}_{bas}}\frac{\epsilon_{n}^{QZ}-\epsilon_{n}^{TZ}}{\frac{1}{N^{QZ}_{bas}}-\frac{1}{N^{TZ}_{bas}}}\;,$
(24)

where $\epsilon_{n}^{QZ}$ ($\epsilon_{n}^{TZ}$) denotes the value of the QP
energy using QZ6P (TZ3P) and $N^{QZ}_{bas}$ and $N^{TZ}_{bas}$ denote the
respective numbers of basis functions (in spherical harmonics so that there
are 5 $d$ and 7 $f$ functions). This expression is commonly used for the
extrapolation of $GW$ QP energies to the complete basis set limit for
localized basis functionsVan Setten  _et al._ (2015). In ref. (109) we
demonstrated that using this extrapolation scheme with Slater type basis sets
good agreement with other codes is obtained for $GW$ QP energies.

##  IV Results

###  IV.1 Total Correlation Energies

Figure 2: Relative errors [in percent] of correlation energies calculates with
different methods compared to the exact values.Chakravorty  _et al._ (1993)
The RPA+SOX and RPA+RSOX values are taken from Engel and coworkersJiang and
Engel (2007) and the gRPA+ values are from Gould et al.Gould  _et al._ (2019)

To see how the inclusion of the $\Sigma^{G3W2}$ self-energy influences the
description of electron correlation effects, we first calculate the total
correlation energies of 16 atoms with in between 2 and 36 electrons. For all
systems with less or equal than 18 electrons, we compare the RPA and RPA+SOSSX
correlation energies to almost exact values by Froese-Fischer and
coworkersChakravorty  _et al._ (1993) and for the heavier elements we use the
CCSD(T) values by McCarthy and Thakkar as reference.McCarthy and Thakkar
(2011) For Argon, their CCSD(T) energy deviates from the value from ref. (127)
by only 0.01 %. We also compare them to different beyond-RPA approaches by
Jiang and EngelJiang and Engel (2007) (RPA+RSOX and RPA+SOX) and Gould et
al.Gould  _et al._ (2019) (gRPA+). To be consistent with ref. (27) and ref.
(46), we evaluate the correlation energies with EXX only OEP orbitals (EXX for
short), implemented within the KLI approximation.Krieger  _et al._ (1990);
Izmaylov  _et al._ (2007)

To obtain an idea about the numerical quality of our RPA correlation energies,
we compare them against the ones from Engel and coworkers which are free of
basis set errors.Jiang and Engel (2007) We find deviations between 3 and 15 %
for neutral atoms and much larger ones for cations. Clearly, our standard
basis sets are not compact enough for these systems (especially for the
cations) and do not capture the full correlation energy. Therefore, we augment
them with four tight 1s and two 2p functions each for TZ3P and QZ6P. Using
these basis sets, total energies for all atoms deviate to the ones from ref.
(27) by about 7 % on average. This is not perfect but accurate enough for a
qualitative comparison of the different beyond-RPA approaches.

The relative errors of correlation energies with respect to the reference
values are shown in figure 2.[note1]See Supplemental Material at [URL will be
inserted by publisher] for all raw data presented in this section. Simple
RPA@EXX overestimates the correlation energies by typically between 25 % and
100 %. In accordance with the expectation that the correlation energy is more
and more dominated by charge screening with increasing electron number, the
agreement with the exact values becomes better for larger atomic numbers.
RPA+SOSSX reduced the RPA correlation energy by 29 % on average, improving
agreement with the reference considerably. Especially for the systems with 2
or 4 electrons, the total correlation energies are still too high. For the
last three systems with 18 electrons or more, the agreement becomes much
better, but the correlation energy is still overestimated by a small amount.
As one can expect, RPA+SOX considerably underestimates the correlation
energies. RPA+RSOX shows a tendency to underestimate the correlation energies
and out of all assessed methods, the deviations to the reference energies are
clearly the smallest. For gRPA+, no clear trend in any direction can be
identified.

In summary, the results presented here show that the $GW+G3W2$ self-energy
improves the description of electron correlation over $GW$. This, however,
does not necessarily imply improvements for properties like relative energies
and charged excitations in realistic systems. For example, RPA is often very
accurate for relative energies since errors in the correlation energies tend
to cancel.Ren  _et al._ (2012) On the other hand, beyond-RPA methods, while
improving total correlation energies, often do not yield improvements for
relative energies. For this reason, we shortly assess the performance of
RPA+SOSSX in comparison to RPA for relative energies, before we move on to
charged excitations.

###  IV.2 Energy Differences

|  | RPA | RPA+SOSSX
---|---|---|---
dataset | $\overline{|\Delta E|}$ | PBE0 | HF | PBE0 | HF
ISO | 21.5 | 1.4 (14.3) | 2.1 (15.1) | 1.0 ( 9.1) | 1.9 (14.1)
FH51 | 31.0 | 1.6 (23.6) | 1.6 (21.2) | 1.6 (23.7) | 1.4 (24.0)
Table 1: Mean absolute deviations (MAD) and mean absolute percentage
deviations (MAPD) (in parantheses) for RPA and RPA+SOSSX using PBE0 and HF
orbitals for the ISO and FH51 datasets. The second column shows the average
value of the datapoints in each set. All values are in kcal/mol.

We now compare the accuracy of RPA+SOSSX against RPA for relative
energies.[note1] We consider HF and PBE0 starting points. Since we do not
evaluate the singles contribution to the total energy, the former is formally
more rigorous due to Brillouin’s theorem. on the other hand, RPA calculations
with PBE0 starting points (RPA@PBE0) are usually more accurate than for HF
(RPA@HF). Our benchmark is based on two datasets. The first one consists of 84
relative energies of different conformers of small and medium organic
molecules, including transition states of inversion and perycyclic reactions.
It is a compilation of the ISO34Grimme  _et al._ (2007); Goerigk  _et al._
(2017), INV24Goerigk and Sharma (2016); Goerigk  _et al._ (2017), and
BHPERIGoerigk and Grimme (2010); Goerigk  _et al._ (2017); Dinadayalane  _et
al._ (2002); Guner  _et al._ (2003); Ess and Houk (2005) datasets, with
reference values obtained using using the W1-F12Hylleraas (1929); Klopper and
Kutzelnigg (1987); Kutzelnigg and Klopper (1991) and W2-F12 protocolsMartin
and De Oliveira (1999); Parthiban and Martin (2001); Karton and Martin (2012),
except for the largest systems in INV24, where reference values have been
calculated on the DLPNO-CCSD(T)/_tight_ Riplinger  _et al._ (2013); Riplinger
and Neese (2013); Liakos  _et al._ (2015)/CBS level of theory. The second one
consists of 51 reaction energies of small organic and inorganic
moleculesGoerigk  _et al._ (2017); Friedrich and Hänchen (2013); Friedrich
(2015) and the reference values have been obtained on the CCSD(T)-F12/CBS
limit of theory.

For RPA@HF and RPA+SOSSX@HF, we report MADs and mean percentage deviations in
table 1. The results are promising. PBE0 reference orbitals generally leads to
better agreement with experiment, for both, RPA and RPA+SOSSX. Also
RPA+SOSSX@HF does not lead to improvements over RPA@HF but RPA+SOSSX@PBE0
leads to significant improvements over RPA@PBE0 for isomerization energies.
Relative errors are reduced from 14 % to 9 %. For FH51, no improvement is
found.

###  IV.3 Quasiparticle Energies

####  IV.3.1 Organic Acceptor Molecules

Figure 3: Ionization potentials (top), electron affinities (middle) and
fundamental gaps (bottom) from qsGW and qsGW + $G3W2$ for a dataset of 24
organic acceptor molecules. The dashed black diagonal lines are CCSD(T)
reference values. Experimental results are given for comparison as well. The
values in parantheses denote mean absolute deviations (MAD) with respect to
CCSD(T). All values are in eV.

We first assess the performance of the $G3W2$ QP energies in comparison to
qsGW for a set of 24 organic acceptor molecules.Richard  _et al._ (2016) The
reference values are of CC singles doubles and perturbative triples [CCSD(T)]
quality and have been extrapolated to the CBS limit. Comparing these
values[note1] to the experimental data shown in figure 3, we see that there
are sizable differences to the CC values, especially for the fundamental gaps.
Among the factors which might contribute to the discrepancies are errors in
the optimized geometries, missing zero point vibrational energy corrections,
and geometry relaxation after oxidation/reduction. For benzonitrile, the
authors of ref. (148) calculated the values of the latter two corrections to
be of the order of 0.18 and 0.14 eV, respectively. The errors in geometry or
thermodynamical contributions are more difficult to assess but can be sizable
as well: For example, the structure used in the calculations might not
correspond to a global minimum on the potential energy surface. Finally, we
note that the basis set extrapolation can also introduce some errors,
especially for the larger systems where no basis sets larger than QZ were
used.Knight  _et al._ (2016) We estimate the error of our own CBS limit
extrapolation to be of the order of 50 meV for IPs and EAs of medium organic
molecules.Förster and Visscher (2021a). Due to all these factors that affect
the direct comparison to experiment, we exclusively use the CCSD(T) reference
values for the following quantitative discussion and only show the
experimental values for comparison.

##### Performance of qsGW

For a variety of solids and metals it has been found that qsGW commonly
overestimates band gaps and IPs by about about 15-20 % when screening is
calculated within the RPA.van Schilfgaarde  _et al._ (2006); Shishkin  _et
al._ (2007); Kang and Hybertsen (2010); Svane  _et al._ (2010); Punya  _et
al._ (2011) In contrast, the qsGW IPs shown in figure 3a are in excellent
agreement with the CC reference values, with no systematic overestimation. The
fundamental gaps in figure 3c are even (a few exceptions aside) systematically
underestimated by qsGW. Overestimation of band gaps with qsGW issue is usually
explained by missing electron-hole interaction via vertex corrections in the
polarizability. Inclusion of an effective two-point kernel from time-dependent
(TD) DFT or the Bethe-Salpeter equation (BSE) has been demonstrated to
significantly improve the agreement of band gaps and IPs with
experiment,Shishkin  _et al._ (2007); Grüneis  _et al._ (2014); Chen and
Pasquarello (2015); Tal  _et al._ (2021); Cunningham  _et al._ (2018, 2021)
demonstrating the importance of beyond-RPA screening. For polar materials,
i.e. materials with strong longitudinal-optical (LO) and transverse-optical
(TO) phonon splitting,Botti and Marques (2013) electron-phonon coupling and
phonon contributions to the frequency-dependent screening can have a sizable
effect on the QP spectrum as well.Bechstedt  _et al._ (2005); Botti and
Marques (2013) For example, qsGW overestimates the experimental band gap of
$\text{V}_{2}\text{O}_{5}$ by about 100 %, which to a large extent is due to
LO-TO splitting.Bhandari  _et al._ (2015)

The systems we consider here are rather small and have a planar geometry. The
reduction of charge screening in low-dimensional materials has often been
emphasized, for example in comparative studies on bulk and layered
$\text{V}_{2}\text{O}_{5}$Bhandari and Lambrecht (2014); Bhandari  _et al._
(2015), $\text{MoS}_{2}$Molina-Sánchez and Wirtz (2011) or polythiophene.van
der Horst  _et al._ (2000) Anti-screening has been observed in a spin chainDe
Groot  _et al._ (1995) and also in finite conjugated systems.Van Den Brink and
Sawatzky (2000) Finally, LO-TO splitting will be absent entirely. These
qualitative differences most likely explain the much higher accuracy of qsGW
for the systems studied herein as compared to periodic systems.

The data clearly demonstrates qsGW to be an excellent first-principle method
for the description of charged excitations for these weakly correlated,
organic molecules. It is worthwhile to compare the performance of this method
for to previous benchmark results of different $GW$ methods. In ref. (83), the
accuracy of a large number of $GW$ methods has been assessed for the same
dataset. Out of all $GW$ based methods, the authors found
$G_{0}W_{0}@$LRC-$\omega$PBE to perform best (with optimized range separation
parameter), with a MAD of 0.13 eV for IPs, and 0.18 eV for EAs. qsGW seems to
be superior.

##### The effect of the statically screened $G3W2$ correction

Now we look at the effect of adding the $G3W2$ correction. As can be seen in
figure 3, on average it lowers the qsGW IPs (see also the MDs in table 2) by a
small amount and decreases the EAs by a relatively larger amount, implying
increasing fundamental gaps. For IPs, this slightly worsens the agreement with
the CC reference values, increasing the mean absolute deviation (MAD) from
0.09 eV to 0.12 eV. However, given the errors from the basis set limit
extrapolation, this difference is not significant. Also, the $G3W2$ correction
does not alter the MAD for the fundamental gaps. Since the statically screened
$G3W2$ correction tends to increase the QP energies, a method which
systematically overestimates IPs and EAs will be systematically improved by
the $G3W2$ correction. This is demonstrated in figure 4 for
$G_{0}W_{0}@\omega$B97-X and $G_{0}W_{0}@$LRC$\omega$PBEh. Here, the $G3W2$
correction slightly improves the MAD with respect to the CC reference values
for the IPs from 0.26 to 0.16 eV, and from 0.16 to 0.13 eV. For the EAs, the
improvements are tremendous, and the inclusion of the $G3W2$ term leads to
almost perfect agreement with the reference values. In table 2, these results
are summarized. Despite the great performance of the $G3W2$ correction for the
range-separated hybrids, the description of fundamental gaps is actually
deteriorated, which can be considered as a serious drawback of this method.

| $GW$ | $GW+G3W2$
---|---|---
| qsGW | $\omega B97X$ |  $\omega$PBEH | qsGW | $\omega B97X$ |  $\omega$PBEH
| Ionization Potentials
MAD | 0.09 | 0.26 | 0.16 | 0.12 | 0.16 | 0.13
MD | 0.08 | 0.26 | 0.16 | 0.00 | 0.14 | 0.05
MAX | 0.43 | 0.54 | 0.41 | 0.46 | 0.57 | 0.54
| Electron Affinities
MAD | 0.14 | 0.29 | 0.27 | 0.09 | 0.05 | 0.07
MD | 0.14 | 0.29 | 0.27 | -0.09 | -0.03 | -0.06
MAX | 0.25 | 0.38 | 0.36 | 0.23 | 0.15 | 0.20
| Fundamental Gaps
MAD | 0.13 | 0.11 | 0.16 | 0.13 | 0.19 | 0.16
MD | -0.06 | -0.03 | -0.12 | 0.09 | 0.17 | 0.11
MAX | 0.35 | 0.25 | 0.28 | 0.60 | 0.66 | 0.62
Table 2: Mean absolute deviations (MAD), mean signed deviattions (MD), and maximum errors (MAX) for IPs, EAs, and fundamental gaps for the 24 acceptor molecules for three different starting points for several $GW$ based methods plus the respective $G3W2$ corrections. All values are in eV. Figure 4: Ionization potentials (left) and electron affinities (right) from $G_{0}W_{0}@\omega$B97-X (+ $\Sigma^{G3W2}$) (top) and $G_{0}W_{0}@$LRC$\omega$-PBEh (+ $\Sigma^{G3W2}$) (bottom) and $G_{0}W_{0}@\omega$B97-X + $\Sigma^{G3W2}$, for a dataset of 24 organic acceptor molecules. The dashed black diagonal lines are CCSD(T) reference values and the values in parantheses denote mean absolute deviations (MAD) with respect to CCSD(T). All values are in eV. Starting Point | $G_{0}W_{0}$ |  $G_{0}W_{0}$+SOSEX | $G_{0}W_{0}+G3W^{(d)}2$ | $G_{0}W_{0}+G3W^{(s)}2$
---|---|---|---|---
| IP | EA | IP | EA | IP | EA | IP | EA
PBE | 0.56/0.66 | 0.44/0.60 | 0.33 | 0.08 | 0.28 | 0.06 | 0.67 | 0.09
PBE0 | 0.19/0.22 | 0.39/0.39 | 0.12 | 0.16 | 0.16 | 0.09 | 0.28 | 0.06
Table 3: Comparison of MADs with respect to the CC reference values for
different beyond-$GW$ method in comparison to $G_{0}W_{0}$ for different
starting points. All values are in eV. The subscript $(d)$ denote that the
screened interaction is dynamic, whereas $(s)$ denotes a static interaction.
For $G_{0}W_{0}$, the two numbers denote MAD with ADF/FHI-AIMS.

##### Comparison of dynamically and statically screened SOX

In table 3 we compare a few beyond-$GW$ approaches for the same dataset. These
are $GW$ \+ SOSEX from ref. (47) and $GW$ \+ dynamically screened $G3W2$ from
ref. (73) (denoted as $G_{0}W_{0}\Gamma^{(1)}_{0}$ by the authors of ref.
(73)). These methods have all been implemented perturbatively and only differ
in the way the electron-electron interaction in the SOX-term is screened (see
for example figure 1 in ref. (73)). The results from refs. (47) and ref. (73)
clearly demonstrate that especially methods with bad or only mediocre
performance, like $G_{0}W_{0}@$PBE or $G_{0}W_{0}@$PBE0, profit immensely from
vertex corrections with fully dynamical self-energy:
$G_{0}W_{0}\Gamma^{(1)}_{0}@$PBE0 is very accurate, with a MAD of 0.16 eV for
the IPs and 0.09 eV for EAs. Especially $G_{0}W_{0}\Gamma^{(1)}_{0}@$PBE
performs much better than $G_{0}W_{0}$PBE. For EAs, the MAD improve from 0.60
eV to 0.06 meV, however, with 0.28 eV for IPs,
$G_{0}W_{0}\Gamma^{(1)}_{0}@$PBE is not very accurate. For EAs, with a MAD of
0.06 eV, $G_{0}W_{0}@$PBE + SOSEX performs excellent, but it is considerably
less accurate for IPs (MAD = 0.33 eV). In summary,
$G_{0}W_{0}\Gamma^{(1)}_{0}@$PBE0 is the most accurate of these methods.
However, it can no beat qsGW and also not $G_{0}W_{0}@\omega$B-97-X + $G3W2$
and $G_{0}W_{0}@$LRC$\omega$PBEh + $G3W2$.

For the sake of a direct comparison of dynamically and statically screened
$G3W2$ corrections, we also calculated $G_{0}W_{0}@$PB + $G3W2$ or
$G_{0}W_{0}@$PBE0 + $G3W2$ with our implementation. First, we note the
reasonable agreement of the results obtained with our implementation and FHI-
AIMS on the $GW$ level, which allows for a qualitative comparison. The
differences mainly result from the different basis sets used in the
calculations.Ren  _et al._ (2015); Wang  _et al._ (2021) The MADs in table 3
clearly show, that the statically screened $G3W2$ correction does not give
good results for these starting points. For EAs, the performance of the
statically screened correction is comparable to the dynamical one. However,
the description of IPs is even worsened. Here, the methods with dynamical
screening are significantly better.

In ref. (73), it was found that the magnitude of the $G3W2$ correction was
much smaller, when the statically screened interaction instead of the
dynamically screened one was used. This is in line with our results. For the
HOMO QP energies, the correction often changes sign, i.e. the HOMO QP energy
is increased when the interaction is statically screened. Also in this case,
the correction for the LUMO level is typically much larger than the one for
the HOMO level, increasing the fundamental gaps. Due to the small magnitude of
the correction, the statically screened $G3W2$ correction works well if QP
energies are already well described on the $GW$ level. The dynamically
screened $G3W2$ term leads to a correction of larger magnitude and works best
for $GW$ methods which severely underestimate IPs and overestimate EAs. This
is the case for $G_{0}W_{0}@$PBE and (to a smaller extent) also for
$G_{0}W_{0}@$PBE0Knight  _et al._ (2016) and consequently, the addition of
SOSEX or dynamically screened $G3W2$ leads to large improvements.Wang  _et
al._ (2021) On the other hand, $G_{0}W_{0}@$HF underestimates the HOMO QP
levels and overestimates LUMO QP levels.Knight  _et al._ (2016); Wang  _et
al._ (2021) Therefore, the addition of dynamically screened SOX deteriorates
the results for this starting point.

####  IV.3.2 Ionization Potentials of Small Molecules

Starting Point | $GW$ | $GW+G3W2$
---|---|---
qsGW | 0.21 | 0.27
$\omega$B97-X | 0.20 | 0.27
LRC-$\omega$PBEh | 0.12 | 0.21
Table 4: MADs with respect to the EOM-CCSDT reference for the first ionization
potentials of a set of 40 small molecules. All values are in eV.

So far, the performance of the $G3W2$ self-energy correction has only been
assessed for a very specific type of molecules. We now also consider a second
database curated by Bartlett and coworkers.Ranasinghe  _et al._ (2019) They
calculated the IPs of 40 small organic and inorganic molecules, using equation
of motion (EOM) CCSD (EOM-CCSD)/cc-pVTZ & cc-pVQZ and EOM-CCSDT/cc-pVTZ. The
reference values we use here are obtained as follows: From the EOM-CCSD
results obtained with the cc-pVTZ and cc-pVQZ basis sets, we extrapolate the
IPs to the CBS limit with the formula by Helgaker et al.Helgaker  _et al._
(1997), eq. (23). Note, that we used eq.(24) to obtain our results.
Subsequently, we add the difference between the EOM-CCSD and EOM-CCSDT IPs to
the CBS limit extrapolated EOM-CCSD IPs. Thus, the reference values should be
close to EOM-CCSDT quality at the CBS limit. The MADs of the considered
methods with respect to the EOM-CCSDT reference values are shown in table 4.
The qsGW IPs are with a MAD of 0.21 eV still in reasonable agreement with the
reference values, but the agreement is worse than for the acceptor molecules.
For all tested $GW$ starting points, the $G3W2$ term worsens the IPs.

###  IV.4 Timings

| Azulene | Borrelidin
---|---|---
| QZ6P | TZ3P | QZ6P
$N_{atoms}$ | 18 | 78
$N_{freq}$ | 18 | 18
$N_{AO}$ | 782 | 1687 | 3184
$N_{aux}$ | 2296 | 9612 | 14257
$GW$ calc. [Core h] | 7.3 | 17 | 74
$G3W2^{*}$ correction [Core h] | 0.015 | 0.40 | 1.6

Table 5: Computational Timings (in core hours) for two selected molecules for
different basis sets. $N_{AO}$ and $N_{aux}$ denote the sizes of primary and
auxiliary basis, respectively.

∗Timing per state

Before concluding this work, we shortly comment on the computational timings
shown in table 5. Compared to SOSEX and $G3W2$ with dynamic $W$, our approach
comes with the advantage that it is computationally very cheap. For example,
the $G3W2$ correction for Azulene at the QZ6P level can be calculated in 2
seconds each for HOMO and LUMO on a 2.2 GHz intel Xeon (E5-2650 v4) node
(broadwell architecture) with 24 cores and 128 GB RAM each, while the
preceding $G_{0}W_{0}$ calculation takes 1100 seconds. We also calculated QP
energies of the Borrelidin molecule (structure taken from ref. (160)) with 78
atoms, 266 electrons and 3200 AOs. The $G_{0}W_{0}$ calculation alone took 74
core hours, while the $G3W2$ correction took 1.6 core hours per state.
Especially if only frontier orbitals are of interest as in many applications,
the $G3W2$ correction is thus computationally inexpensive. The comparison is
slightly flawed, since our low-scaling $GW$ implementation comes with a higher
prefactor than canonical implementations. The $GW$ implementation scales
quadratically, while the $G3W2$ corrections is naturally implemented with
$N^{4}$ scaling. Consequently, for much larger systems, The $G3W2$ correction
will become computationally more demanding. The corresponding expression for
the energy scales as $N^{5}$ (as opposed to our quadratic scaling RPA
implementation) but in canonical implementations MP2 energies can be routinely
calculated for systems with more than 100 atoms as well.Cooper  _et al._
(2017)

##  V Conclusion

In the GWA, the electron-electron self-energy in Dyons’s equation is expanded
in terms of a screened Coulomb interaction and truncated after first
order.Hedin (1965) In this work, we have analysed different aspects of the
statically screened second-order contribution to the self-energy ($G3W2$) and
applied it in a perturbative fashion to calculate correlation energies of
atoms, relative energies of chemical reactions, and IPs and EAs of a wide
range of molecules. The results we have presented herein shed some light on
many interesting aspects of the second-order correction and also raise many
questions which we hope to be able to address in future work.

For correlation energies, the second-order correction (SOSSX) reduces the
RPA@EXX correlation energies by around 25 % and brings them into much better
agreement with nearly exact reference values,Chakravorty  _et al._ (1993);
McCarthy and Thakkar (2011) but the correlation energies are still too large.
Generally, relative energies are also improved over RPA, but more detailed
investigations are needed to obtain a better understanding in what situations
RPA+SOSSX yields improvements over the RPA.

It would also be interesting to investigate the relation between the RSOX and
SOSSX terms further. SOSSX and RSOX are two different strategies to
renormalizes the SOX term. RPA+RSOX is based on the resummation of the
Epstein-Nesbet seriesSzabo, Attila and Ostlund (2012) of hole-hole ladder
diagramsEngel and Jiang (2006); Jiang and Engel (2006) while SOSSX results
from the resummation of ring diagrams, with the additional approximation of
static screening. Both classes of diagrams are contained in CCD which also
includes particle-particle ladder diagrams and a third class of diagrams,
coupling ladder and ring diagrams.Scuseria  _et al._ (2013) Thus, CCD can be
seen as a first-principle method to couple both renormalization strategies.
Unfortunately, CCD scales as $N^{6}$, as opposed to $N^{5}$. Since RPA+SOSSX
overestimates, and RSOX underestimates the correlation energies, it seems
worthwhile to look for alternative ways to combine SOSSX with RSOX while
retaining the $N^{5}$ scaling.

For charged excitations, our work reveals qualitative differences between the
dynamically screened $G3W2$ correction as recently tested in ref. (73) and the
statically screened one, benchmarked herein. The magnitude of the $G3W2$
correction becomes much smaller when the electron-electron interaction is
statically screened. The addition of the dynamically screened SOX term
consistently lowers EAs and increases IPs, while the statically screened one
often decreases the IPs.

Popular methods like $G_{0}W_{0}$ based on PBE and PBE0 starting points
predict too low ionization potentials and too high electron affinitiesWang
_et al._ (2021) and the addition of dynamically screened $G3W2$ correction
results in major improvements for these methods.Wang  _et al._ (2021) The
statically screened $G3W2$ correction, on the other hand, gives improvements
over $GW$ calculations which consistently underestimate QP energies. This is
the case for $G_{0}W_{0}@LRC\omega$PBEh and $G_{0}W_{0}@\omega$B97-X.
Especially for EAs, improvements are tremendous. For IPs and fundamental gaps,
the improvements are not consistent and seem to be system-specific. For
molecules with nearly 100 atoms, evaluating the perturbative $G3W2$ correction
for a few states does not come with significantly increased computational
costs compared to a $GW$ calculation. This is a big advantage over the
dynamically screened SOX correction.Wang  _et al._ (2021)

Also for qsGW, the statically screened $G3W2$ correction does not lead to
systematic improvements, even though both corrections improve fundamental
gaps. As an important byproduct of this work, we have shown that qsGW by
itself is an excellent method to calculate IPs and EAs for a set of medium
organic acceptor molecules. With MADs of 0.09 eV for IPs and 0.14 eV for EAs,
qsGW outperforms all other $GW$ methods, previously benchmarked for this
dataset.Knight  _et al._ (2016) Remarkably, fundamental gaps are
underestimated by qsGW compared to CCSD(T) reference values. This is different
from the common situation in extended systems and can be attributed to the
much weaker charge screening in finite systems. Also for a second set of 40
small molecules where no reference data for EAs was available, qsGW gives
excellent IPs. Charge screening can be expected to be stronger in larger
molecules and it would be interesting to see whether qsGW also performs well
for these or if beyond-RPA screening is then necessary.

Even though the SOX correction to the self-energy leads to improvements over
$GW$ in certain situations, our results confirm recent studiesLewis and
Berkelbach (2019); Bruneval  _et al._ (2021) which show that QP approximations
to $GW$ are difficult to improve upon diagrammatically, at least for molecular
systems. Systematic and reliable improvements over $GW$ are most likely only
possible starting from the fully self-consistent solution of the $GW$
equations.

###### Acknowledgements.

This research received funding from the Netherlands Organisation for
Scientific Research (NWO) in the framework of the Innovation Fund for
Chemistry and from the Ministry of Economic Affairs in the framework of the
_TKI/PPS-Toeslagregeling_.

## References

  * Abrikosov  _et al._ (1975) A. A. Abrikosov, P. L. Gorkov, and I. E. Dzyaloshinski,  _Methods of quantum field theory in statistical physics_ (Dover Publications INC. New York, 1975).
  * Mattuck (1992) R. D. Mattuck,  _A Guide to Feynman Diagrams in the Many-body Problem_ , 2nd ed. (Dover Publications INC. New York, 1992).
  * Martin  _et al._ (2016) R. M. Martin, L. Reining, and D. M. Ceperley,  _Interacting electrons_ (Cambridge University Press, 2016).
  * Dyson (1949) F. J. Dyson, Phys. Rev. 75, 1736 (1949).
  * Hedin (1965) L. Hedin, Phys. Rev. 139, A796 (1965).
  * Phillips (1961) J. C. Phillips, Phys. Rev. 123, 420 (1961).
  * Hubbard (1957) J. Hubbard, Proc. R. Soc. London. Ser. A. Math. Phys. Sci. 240, 539 (1957).
  * Reining (2018) L. Reining, Wiley Interdiscip. Rev. Comput. Mol. Sci. 8, e1344 (2018).
  * Golze  _et al._ (2019) D. Golze, M. Dvorak, and P. Rinke, Front. Chem. 7, 377 (2019), arXiv:1912.04893 .
  * Hybertsen and Louie (1986) M. S. Hybertsen and S. G. Louie, Phys. Rev. B 34, 5390 (1986).
  * Kotani  _et al._ (2007) T. Kotani, M. van Schilfgaarde, and S. V. Faleev, Phys. Rev. B 76, 165106 (2007), arXiv:0611002 [cond-mat] .
  * Ward (1950) J. C. Ward, Phys. Rev. 78, 182 (1950).
  * Klein (1961) A. Klein, Phys. Rev. 121, 950 (1961).
  * Casida (1995) M. E. Casida, Phys. Rev. A 51, 2005 (1995).
  * Dahlen  _et al._ (2006) N. E. Dahlen, R. van Leeuwen, and U. von Barth, Phys. Rev. A - At. Mol. Opt. Phys. 73, 012511 (2006).
  * Ren  _et al._ (2012) X. Ren, P. Rinke, C. Joas, and M. Scheffler, J. Mater. Sci. 47, 7447 (2012), arXiv:1203.5536 .
  * Coester (1958) F. Coester, Nucl. Phys. 7, 421 (1958).
  * Coester and Kümmel (1960) F. Coester and H. Kümmel, Nucl. Phys. 17, 477 (1960).
  * Čížek (1966) J. Čížek, J. Chem. Phys. 45, 4256 (1966).
  * Čížek (1969) J. Čížek, Adv. Chem. Physics, Vol. XIV XIV, 35 (1969).
  * Paldus  _et al._ (1972) J. Paldus, J. Čížek, and I. Shavitt, Phys. Rev. A 5, 50 (1972).
  * Scuseria  _et al._ (2008) G. E. Scuseria, T. M. Henderson, and D. C. Sorensen, J. Chem. Phys. 129, 231101 (2008).
  * Scuseria  _et al._ (2013) G. E. Scuseria, T. M. Henderson, and I. W. Bulik, J. Chem. Phys. 139 (2013), 10.1063/1.4820557, arXiv:1306.6360 .
  * Langreth and Perdew (1977) D. C. Langreth and J. P. Perdew, Phys. Rev. B 15, 2884 (1977).
  * Singwi  _et al._ (1968) K. S. Singwi, M. P. Tosi, R. H. Land, and A. Sjölander, Phys. Rev. 176, 589 (1968).
  * Irmler  _et al._ (2019) A. Irmler, A. Gallo, F. Hummel, and A. Grüneis, Phys. Rev. Lett. 123, 156401 (2019), arXiv:1903.05458 .
  * Jiang and Engel (2007) H. Jiang and E. Engel, J. Chem. Phys. 127, 184108 (2007).
  * Kurth  _et al._ (1999) S. Kurth, J. P. Perdew, and P. Blaha, Int. J. Quantum Chem. 75, 889 (1999).
  * Yan  _et al._ (2000) Z. Yan, J. P. Perdew, and S. Kurth, Phys. Rev. B - Condens. Matter Mater. Phys. 61, 16430 (2000).
  * Angyan  _et al._ (2005) J. G. Angyan, I. C. Gerber, A. Savin, and J. Toulouse, Phys. Rev. A 72, 012510 (2005).
  * Janesko  _et al._ (2009a) B. G. Janesko, T. M. Henderson, and G. E. Scuseria, J. Chem. Phys. 130, 081105 (2009a).
  * Janesko  _et al._ (2009b) B. G. Janesko, T. M. Henderson, and G. E. Scuseria, J. Chem. Phys. 131, 034110 (2009b).
  * Toulouse  _et al._ (2009) J. Toulouse, I. C. Gerber, G. Jansen, A. Savin, and J. G. Ángyán, Phys. Rev. Lett. 102, 096404 (2009), arXiv:0812.3302 .
  * Zhu  _et al._ (2010) W. Zhu, J. Toulouse, A. Savin, and J. G. Ángyán, J. Chem. Phys. 132, 244108 (2010).
  * Toulouse  _et al._ (2010) J. Toulouse, W. Zhu, J. G. Ángyán, and A. Savin, Phys. Rev. A - At. Mol. Opt. Phys. 82, 032502 (2010), arXiv:1006.2061 .
  * Freeman (1977) D. Freeman, Phys. Rev. B 15, 5512 (1977), arXiv:arXiv:1011.1669v3 .
  * Paier  _et al._ (2012) J. Paier, X. Ren, P. Rinke, G. E. Scuseria, A. Grüneis, G. Kresse, and M. Scheffler, New J. Phys. 14, 043002 (2012), arXiv:1111.0173 .
  * Ren  _et al._ (2013) X. Ren, P. Rinke, G. E. Scuseria, and M. Scheffler, Phys. Rev. B - Condens. Matter Mater. Phys. 88, 035120 (2013).
  * Grüneis  _et al._ (2009) A. Grüneis, M. Marsman, J. Harl, L. Schimka, and G. Kresse, J. Chem. Phys. 131, 154115 (2009).
  * Hummel  _et al._ (2019) F. Hummel, A. Grüneis, G. Kresse, and P. Ziesche, J. Chem. Theory Comput. 15, 3223 (2019), arXiv:1605.08643 .
  * Paier  _et al._ (2008) J. Paier, M. Marsman, and G. Kresse, Phys. Rev. B - Condens. Matter Mater. Phys. 78, 121201(R) (2008).
  * Toulouse  _et al._ (2011) J. Toulouse, W. Zhu, A. Savin, G. Jansen, and J. G. Ángyán, J. Chem. Phys. 135, 084119 (2011).
  * Ángyán  _et al._ (2011) J. G. Ángyán, R. F. Liu, J. Toulouse, and G. Jansen, J. Chem. Theory Comput. 7, 3116 (2011), arXiv:1404.1663 .
  * Engel and Jiang (2006) E. Engel and H. Jiang, Int. J. Quantum Chem. 106, 3242 (2006).
  * Jiang and Engel (2006) H. Jiang and E. Engel, J. Chem. Phys. 125, 1 (2006).
  * Gould  _et al._ (2019) T. Gould, A. Ruzsinszky, and J. P. Perdew, Phys. Rev. A 100, 022515 (2019).
  * Ren  _et al._ (2015) X. Ren, N. Marom, F. Caruso, M. Scheffler, and P. Rinke, Phys. Rev. B - Condens. Matter Mater. Phys. 92, 081104(R) (2015).
  * Maggio and Kresse (2017) E. Maggio and G. Kresse, J. Chem. Theory Comput. 13, 4765 (2017), arXiv:1706.01815 .
  * Kutepov (2016) A. L. Kutepov, Phys. Rev. B 94, 155101 (2016).
  * Kutepov and Kotliar (2017) A. L. Kutepov and G. Kotliar, Phys. Rev. B 96, 035108 (2017).
  * Kutepov (2017) A. L. Kutepov, Phys. Rev. B 95, 195120 (2017).
  * Kutepov (2018) A. L. Kutepov, arXiv:1809.06654v2 (2018), arXiv:1809.06654 .
  * Kutepov (2021a) A. L. Kutepov, Phys. Rev. B 104, 085109 (2021a), arXiv:2105.03770 .
  * Kutepov (2021b) A. L. Kutepov, Phys. Rev. Mater. 5, 083805 (2021b).
  * Kutepov (2022) A. L. Kutepov, Phys. Rev. B 105, 045124 (2022), arXiv:2109.01021 .
  * Onida  _et al._ (1995) G. Onida, L. Reining, R. W. Godby, R. Del Sole, and W. Andreoni, Phys. Rev. Lett. 75, 818 (1995).
  * Onida  _et al._ (2002) G. Onida, L. Reining, and A. Rubio, Rev. Mod. Phys. 74, 601 (2002).
  * Adragna  _et al._ (2003) G. Adragna, R. Del Sole, and A. Marini, Phys. Rev. B - Condens. Matter Mater. Phys. 68, 165108 (2003).
  * Marini  _et al._ (2003) A. Marini, R. Del Sole, and A. Rubio, Phys. Rev. Lett. 91, 256402 (2003), arXiv:0310495 [cond-mat] .
  * Bruneval  _et al._ (2005) F. Bruneval, F. Sottile, V. Olevano, R. Del Sole, and L. Reining, Phys. Rev. Lett. 94, 186402 (2005), arXiv:0504565 [cond-mat] .
  * Shishkin  _et al._ (2007) M. Shishkin, M. Marsman, and G. Kresse, Phys. Rev. Lett. 99, 246403 (2007).
  * Romaniello  _et al._ (2009) P. Romaniello, S. Guyot, and L. Reining, J. Chem. Phys. 131, 154111 (2009).
  * Sharma  _et al._ (2011) S. Sharma, J. K. Dewhurst, A. Sanna, and E. K. Gross, Phys. Rev. Lett. 107, 186401 (2011), arXiv:1107.0199 .
  * Chen and Pasquarello (2015) W. Chen and A. Pasquarello, Phys. Rev. B - Condens. Matter Mater. Phys. 92, 041115(R) (2015).
  * Hung  _et al._ (2016) L. Hung, F. H. da Jornada, J. Souto-Casares, J. R. Chelikowsky, S. G. Louie, and S. Öǧüt, Phys. Rev. B 94, 085125 (2016).
  * Cunningham  _et al._ (2018) B. Cunningham, M. Grüning, P. Azarhoosh, D. Pashov, and M. van Schilfgaarde, Phys. Rev. Mater. 2, 034603 (2018).
  * Olsen  _et al._ (2019) T. Olsen, C. E. Patrick, J. E. Bates, A. Ruzsinszky, and K. S. Thygesen, Nat. Comput. Mater. 5, 106 (2019).
  * Schmidt  _et al._ (2017) P. S. Schmidt, C. E. Patrick, and K. S. Thygesen, Phys. Rev. B 96, 205206 (2017), arXiv:1711.02922 .
  * Tal  _et al._ (2021) A. Tal, W. Chen, and A. Pasquarello, Phys. Rev. B 103, 161104 (2021).
  * Cunningham  _et al._ (2021) B. Cunningham, M. Gruening, D. Pashov, and M. van Schilfgaarde, arXiv:2106.05759v1 , 1 (2021), arXiv:2106.05759 .
  * Grüneis  _et al._ (2014) A. Grüneis, G. Kresse, Y. Hinuma, and F. Oba, Phys. Rev. Lett. 112, 096401 (2014).
  * van Schilfgaarde  _et al._ (2006) M. van Schilfgaarde, T. Kotani, and S. Faleev, Phys. Rev. Lett. 96, 226402 (2006), arXiv:0510408 [cond-mat] .
  * Wang  _et al._ (2021) Y. Wang, P. Rinke, and X. Ren, J. Chem. Theory Comput. 17, 5140 (2021), arXiv:2105.07386 .
  * Van Den Brink and Sawatzky (2000) J. Van Den Brink and G. A. Sawatzky, Europhys. Lett. 50, 447 (2000).
  * Verdozzi  _et al._ (1995) C. Verdozzi, R. W. Godby, and S. Holloway, Phys. Rev. Lett. 74, 2327 (1995).
  * Schindlmayr and Godby (1998) A. Schindlmayr and R. W. Godby, Phys. Rev. Lett. 80, 1702 (1998), arXiv:9710295 [cond-mat] .
  * Lewis and Berkelbach (2019) A. M. Lewis and T. C. Berkelbach, J. Chem. Theory Comput. 15, 2925 (2019), arXiv:1810.00456 .
  * Hirata  _et al._ (2017) S. Hirata, A. E. Doran, P. J. Knowles, and J. V. Ortiz, J. Chem. Phys. 147, 044108 (2017).
  * Bruneval  _et al._ (2006a) F. Bruneval, N. Vast, and L. Reining, Phys. Rev. B - Condens. Matter Mater. Phys. 74, 045102 (2006a).
  * Bruneval  _et al._ (2006b) F. Bruneval, N. Vast, L. Reining, M. Izquierdo, F. Sirotti, and N. Barrett, Phys. Rev. Lett. 97, 267601 (2006b).
  * Bruneval and Marques (2013) F. Bruneval and M. Marques, J. Chem. Theory Comput. 9, 324 (2013).
  * Caruso  _et al._ (2016) F. Caruso, M. Dauth, M. J. Van Setten, and P. Rinke, J. Chem. Theory Comput. 12, 5076 (2016), arXiv:1609.04665 .
  * Knight  _et al._ (2016) J. W. Knight, X. Wang, L. Gallandi, O. Dolgounitcheva, X. Ren, J. V. Ortiz, P. Rinke, T. Körzdörfer, and N. Marom, J. Chem. Theory Comput. 12, 615 (2016).
  * Atalla  _et al._ (2013) V. Atalla, M. Yoon, F. Caruso, P. Rinke, and M. Scheffler, Phys. Rev. B - Condens. Matter Mater. Phys. 88, 165122 (2013).
  * Hellgren  _et al._ (2021) M. Hellgren, L. Baguet, M. Calandra, F. Mauri, and L. Wirtz, Phys. Rev. B 103, 075101 (2021), arXiv:2005.14603 .
  * Kanchanakungwankul and Truhlar (2021) S. Kanchanakungwankul and D. G. Truhlar, J. Chem. Theory Comput. 17, 4823 (2021).
  * Refaely-Abramson  _et al._ (2011) S. Refaely-Abramson, R. Baer, and L. Kronik, Phys. Rev. B - Condens. Matter Mater. Phys. 84, 075144 (2011).
  * Refaely-Abramson  _et al._ (2012) S. Refaely-Abramson, S. Sharifzadeh, N. Govind, J. Autschbach, J. B. Neaton, R. Baer, and L. Kronik, Phys. Rev. Lett. 109, 226405 (2012), arXiv:1203.2357 .
  * Kronik  _et al._ (2012) L. Kronik, T. Stein, S. Refaely-Abramson, and R. Baer, J. Chem. Theory Comput. 8, 1515 (2012).
  * Godby  _et al._ (1988) R. W. Godby, M. Schlüter, and L. J. Sham, Phys. Rev. B 37, 10159 (1988).
  * Sham and Schl̈üter (1983) L. J. Sham and M. Schl̈üter, Phys. Rev. Lett. 51, 1888 (1983).
  * Sharp and Horton (1953) T. Sharp and G. Horton, Phys. Rev. 90, 317 (1953).
  * Talman and Shadwick (1976) J. D. Talman and W. F. Shadwick, Phys. Rev. A 14, 36 (1976).
  * Engel, Eberhard and Dreizler (2013) R. M. Engel, Eberhard and Dreizler,  _Density functional theory An Advanced Course_ (Springer, 2013).
  * Ismail-Beigi (2017) S. Ismail-Beigi, J. Phys. Condens. Matter 29, 385501 (2017), arXiv:arXiv:1406.0772v4 .
  * Klimeš and Kresse (2014) J. Klimeš and G. Kresse, J. Chem. Phys. 140, 054516 (2014).
  * Peter Kopietz (1997) Peter Kopietz,  _Bosonization of Interacting Fermions in Arbitrary Dimensions, Chapter 1_ (Springer, 1997).
  * Lee  _et al._ (2020) Y. Lee, T. Kotani, and L. Ke, Phys. Rev. B 101, 241409(R) (2020).
  * Riemelmoser  _et al._ (2021) S. Riemelmoser, M. Kaltak, and G. Kresse, J. Chem. Phys. 154, 154103 (2021), arXiv:2101.10009 .
  * Stefanucci  _et al._ (2014) G. Stefanucci, Y. Pavlyukh, A. M. Uimonen, and R. van Leeuwen, Phys. Rev. B 90, 115134 (2014), arXiv:1408.6163 .
  * Pavlyukh  _et al._ (2020) Y. Pavlyukh, G. Stefanucci, and R. van Leeuwen, Phys. Rev. B 102, 045121 (2020), arXiv:2004.05344 .
  * van Leeuwen  _et al._ (2015) R. van Leeuwen, N. E. Dahlen, G. Stefanucci, C. O. Almbladh, and U. Von Barth, in  _Time-Dependent Density Funct. Theory_, edited by M. A. Marques, C. A. Ullrich, F. Nogueira, A. Rubio, K. Burke, and E. K. Gross (Springer Heidelberg, 2015) pp. 185–217.
  * Förster and Visscher (2020) A. Förster and L. Visscher, J. Chem. Theory Comput. 16, 7381 (2020).
  * Pavlyukh  _et al._ (2016) Y. Pavlyukh, A. M. Uimonen, G. Stefanucci, and R. van Leeuwen, Phys. Rev. Lett. 117, 206402 (2016).
  * Willow  _et al._ (2012) S. Y. Willow, K. S. Kim, and S. Hirata, J. Chem. Phys. 137, 204122 (2012).
  * Willow  _et al._ (2013) S. Y. Willow, K. S. Kim, and S. Hirata, J. Chem. Phys. 138, 164111 (2013).
  * Mardirossian  _et al._ (2018) N. Mardirossian, J. D. McClain, and G. K. L. Chan, J. Chem. Phys. 148, 044106 (2018).
  * Niquet and Gonze (2004) Y. M. Niquet and X. Gonze, Phys. Rev. B - Condens. Matter Mater. Phys. 70, 245115 (2004).
  * Förster and Visscher (2021a) A. Förster and L. Visscher, J. Chem. Theory Comput. 17, 5080 (2021a).
  * Förster and Visscher (2021b) A. Förster and L. Visscher, Front. Chem. 9, 736591 (2021b).
  * Kaltak  _et al._ (2014a) M. Kaltak, J. Klimeš, and G. Kresse, J. Chem. Theory Comput. 10, 2498 (2014a).
  * Kaltak  _et al._ (2014b) M. Kaltak, J. Klimeš, and G. Kresse, Phys. Rev. B 90, 054115 (2014b).
  * Liu  _et al._ (2016) P. Liu, M. Kaltak, J. Klimeš, and G. Kresse, Phys. Rev. B 94, 165109 (2016), arXiv:1607.02859 .
  * Förster  _et al._ (2020) A. Förster, M. Franchini, E. van Lenthe, and L. Visscher, J. Chem. Theory Comput. 16, 875 (2020).
  * Baerends  _et al._ (1973) E. J. Baerends, D. E. Ellis, and P. Ros, Chem. Phys. 2, 41 (1973).
  * Cooper  _et al._ (2017) B. Cooper, S. Girdlestone, P. Burovskiy, G. Gaydadjiev, V. Averbukh, P. J. Knowles, and W. Luk, J. Chem. Theory Comput. 13, 5265 (2017).
  * Ernzerhof and Scuseria (1999) M. Ernzerhof and G. E. Scuseria, J. Chem. Phys. 110, 5029 (1999).
  * Adamo and Barone (1999) C. Adamo and V. Barone, J. Chem. Phys. 110, 6158 (1999).
  * Izmaylov  _et al._ (2007) A. F. Izmaylov, V. N. Staroverov, G. E. Scuseria, E. R. Davidson, G. Stoltz, and E. Canc̀s, J. Chem. Phys. 126, 084107 (2007).
  * Krieger  _et al._ (1990) J. B. Krieger, Y. Li, and G. J. Iafrate, Phys. Lett. A 146, 256 (1990).
  * Helgaker  _et al._ (1997) T. Helgaker, W. Klopper, H. Koch, and J. Noga, J. Chem. Phys. 106, 9639 (1997).
  * Jensen (2013) F. Jensen, Wiley Interdiscip. Rev. Comput. Mol. Sci. 3, 273 (2013).
  * (123) See Supplemental Material at [URL will be inserted by publisher] for all newly composed basis sets used in this work.
  * Rohrdanz  _et al._ (2009) M. A. Rohrdanz, K. M. Martins, and J. M. Herbert, J. Chem. Phys. 130, 054112 (2009).
  * Chai and Head-Gordon (2008) J. D. Chai and M. Head-Gordon, J. Chem. Phys. 128, 084106 (2008).
  * Van Setten  _et al._ (2015) M. J. Van Setten, F. Caruso, S. Sharifzadeh, X. Ren, M. Scheffler, F. Liu, J. Lischner, L. Lin, J. R. Deslippe, S. G. Louie, C. Yang, F. Weigend, J. B. Neaton, F. Evers, and P. Rinke, J. Chem. Theory Comput. 11, 5665 (2015).
  * Chakravorty  _et al._ (1993) S. J. Chakravorty, S. R. Gwaltney, E. R. Davidson, F. A. Parpia, and C. F. Fischer, Phys. Rev. A 47, 3649 (1993).
  * McCarthy and Thakkar (2011) S. P. McCarthy and A. J. Thakkar, J. Chem. Phys. 134, 044102 (2011).
  * (129) See Supplemental Material at [URL will be inserted by publisher] for all raw data presented in this section.
  * Grimme  _et al._ (2007) S. Grimme, M. Steinmetz, and M. Korth, J. Org. Chem. 72, 2118 (2007).
  * Goerigk  _et al._ (2017) L. Goerigk, A. Hansen, C. Bauer, S. Ehrlich, A. Najibi, and S. Grimme, Phys. Chem. Chem. Phys. 19, 32184 (2017).
  * Goerigk and Sharma (2016) L. Goerigk and R. Sharma, Can. J. Chem. 94, 1133 (2016).
  * Goerigk and Grimme (2010) L. Goerigk and S. Grimme, J. Chem. Theory Comput. 6, 107 (2010).
  * Dinadayalane  _et al._ (2002) T. C. Dinadayalane, R. Vijaya, A. Smitha, and G. N. Sastry, J. Phys. Chem. A 106, 1627 (2002).
  * Guner  _et al._ (2003) V. Guner, K. S. Khuong, A. G. Leach, P. S. Lee, M. D. Bartberger, and K. N. Houk, J. Phys. Chem. A 107, 11445 (2003).
  * Ess and Houk (2005) D. H. Ess and K. N. Houk, J. Phys. Chem. A 109, 9542 (2005).
  * Hylleraas (1929) E. A. Hylleraas, Zeitschrift für Phys. 54, 347 (1929).
  * Klopper and Kutzelnigg (1987) W. Klopper and W. Kutzelnigg, Chem. Phys. Lett. 134, 17 (1987).
  * Kutzelnigg and Klopper (1991) W. Kutzelnigg and W. Klopper, J. Chem. Phys. 94, 1985 (1991).
  * Martin and De Oliveira (1999) J. M. Martin and G. De Oliveira, J. Chem. Phys. 111, 1843 (1999), arXiv:9904038 [physics] .
  * Parthiban and Martin (2001) S. Parthiban and J. M. Martin, J. Chem. Phys. 114, 6014 (2001), arXiv:0101072 [physics] .
  * Karton and Martin (2012) A. Karton and J. M. Martin, J. Chem. Phys. 136, 124114 (2012).
  * Riplinger  _et al._ (2013) C. Riplinger, B. Sandhoefer, A. Hansen, and F. Neese, J. Chem. Phys. 139, 134101 (2013).
  * Riplinger and Neese (2013) C. Riplinger and F. Neese, J. Chem. Phys. 138, 034106 (2013).
  * Liakos  _et al._ (2015) D. G. Liakos, M. Sparta, M. K. Kesharwani, J. M. Martin, and F. Neese, J. Chem. Theory Comput. 11, 1525 (2015).
  * Friedrich and Hänchen (2013) J. Friedrich and J. Hänchen, J. Chem. Theory Comput. 9, 5381 (2013).
  * Friedrich (2015) J. Friedrich, J. Chem. Theory Comput. 11, 3596 (2015).
  * Richard  _et al._ (2016) R. M. Richard, M. S. Marshall, O. Dolgounitcheva, J. V. Ortiz, J. L. Brédas, N. Marom, and C. D. Sherrill, J. Chem. Theory Comput. 12, 595 (2016).
  * Kang and Hybertsen (2010) W. Kang and M. S. Hybertsen, Phys. Rev. B - Condens. Matter Mater. Phys. 82, 085203 (2010), arXiv:1006.4085 .
  * Svane  _et al._ (2010) A. Svane, N. E. Christensen, I. Gorczyca, M. van Schilfgaarde, A. N. Chantis, and T. Kotani, Phys. Rev. B - Condens. Matter Mater. Phys. 82, 115102 (2010).
  * Punya  _et al._ (2011) A. Punya, W. Lambrecht, and M. van Schilfgaarde, Phys. Rev. B - Condens. Matter Mater. Phys. 84, 165204 (2011).
  * Botti and Marques (2013) S. Botti and M. Marques, Phys. Rev. Lett. 110, 226404 (2013), arXiv:1304.7911 .
  * Bechstedt  _et al._ (2005) F. Bechstedt, K. Seino, P. H. Hahn, and W. G. Schmidt, Phys. Rev. B - Condens. Matter Mater. Phys. 72, 245114 (2005).
  * Bhandari  _et al._ (2015) C. Bhandari, W. Lambrecht, and M. van Schilfgaarde, Phys. Rev. B - Condens. Matter Mater. Phys. 91, 125116 (2015).
  * Bhandari and Lambrecht (2014) C. Bhandari and W. Lambrecht, Phys. Rev. B - Condens. Matter Mater. Phys. 89, 045109 (2014).
  * Molina-Sánchez and Wirtz (2011) A. Molina-Sánchez and L. Wirtz, Phys. Rev. B - Condens. Matter Mater. Phys. 84, 155413 (2011).
  * van der Horst  _et al._ (2000) J. van der Horst, P. A. Bobbert, P. de Jong, M. Michels, G. Brocks, and P. Kelly, Phys. Rev. B - Condens. Matter Mater. Phys. 61, 15817 (2000).
  * De Groot  _et al._ (1995) H. J. De Groot, P. A. Bobbert, and W. van Haeringen, Phys. Rev. B 52, 11000 (1995).
  * Ranasinghe  _et al._ (2019) D. S. Ranasinghe, J. T. Margraf, A. Perera, and R. J. Bartlett, J. Chem. Phys. 150, 074108 (2019).
  * Ma and Werner (2020) Q. Ma and H. J. Werner, J. Chem. Theory Comput. 16, 3135 (2020).
  * Szabo, Attila and Ostlund (2012) N. S. Szabo, Attila and Ostlund,  _Modern quantum chemistry: introduction to advanced electronic structure theory_ (Courier Corporation, 2012).
  * Bruneval  _et al._ (2021) F. Bruneval, N. Dattani, and M. J. van Setten, Front. Chem. 9, 749779 (2021).