CHAPTER 1 Introduction

The physical world is composed of interacting particles. In condensed-matter physics, electrons determine many material properties, such as electronic and thermal conductivity, magnetization, etc. Therefore, constructing an ab-initio theory that explains and also predicts electronic behavior in materials is very important. In general it is impossible to solve the Schrödinger equation due to the exponential growth of the dimensions of the many-body wavefunction with increasing number of electrons [START_REF] Montina | Exponential complexity and ontological theories of quantum mechanics[END_REF]. Approximations are therefore needed. Some of the most successful are Density Functional Theory (DFT) and Many-Body Perturbation Theory (MBPT). These theories avoid the use of the wavefunction and use simpler quantities, such as the electron density, and the one-body Green's function . By using a simplified framework one has to introduce effective potentials containing the complexity of the original problem.

In DFT [START_REF] Hohenberg | Inhomogeneous electron gas[END_REF][START_REF] Kohn | Self-consistent equations including exchange and correlation effects[END_REF], one replaces the interacting system with a fictitious system of non-interacting particles, the Kohn-Sham (KS) system, in which, under the influence of an effective potential (the Kohn-Sham potential), the electron density of the interacting system is reproduced. One then only has to solve a set of single-particle equations. These equations have to be solved self-consistently since the KS potential Introduction is a functional of the density. However, approximations to the Kohn-Sham potential are not intuitive, and, moreover, it is not always known how to extract the observable of interest from the density.

MBPT based on Green's function [START_REF] Bechstedt | Many-body approach to electronic excitations[END_REF][START_REF] Nakanishi | A general survey of the theory of the bethe-salpeter equation[END_REF][START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF], is the method of interest in this thesis. Standard approaches are based on the 1-GF or the 2-GF. The 1-GF is very useful to calculate photoemission spectra since the photoemission spectrum can be obtained in a straightforward way from the imaginary part of the 1-GF. Moreover, it gives access to many other properties, such as potential energy surfaces, electron densities, momentum distributions etc. Within the so called GW approximation [START_REF] Hedin | New method for calculating the one-particle green's function with application to the electron-gas problem[END_REF], MBPT has become the method of choice for the calculations of quasiparticle band structures [1,2,[START_REF] Ia Nechaev | Evidence for a direct band gap in the topological insulator bi 2 se 3 from theory and experiment[END_REF][START_REF] Vidal | Strong interplay between structure and electronic properties in CuIn(S, Se) 2 : A first-principles study[END_REF][START_REF] Vidal | Falsepositive and false-negative assignments of topological insulators in density functional theory and hybrids[END_REF][START_REF] Waroquiers | Band widths and gaps from the tran-blaha functional: Comparison with many-body perturbation theory[END_REF] and direct and inverse photo-emission spectra [START_REF] Athanasios | Quasiparticle self-consistent gw method applied to localized 4f electron systems[END_REF][START_REF] Sergey | All-electron selfconsistent gw approximation: Application to si, mno, and nio[END_REF][START_REF] Gatti | Understanding correlations in vanadium dioxide from first principles[END_REF][START_REF] Lim | Angle-resolved photoemission and quasiparticle calculation of zno: The need for d band shift in oxide semiconductors[END_REF][START_REF] Papalazarou | Valence-band electronic structure of v 2 o 3 : Identification of v and o bands[END_REF] of many materials improving substantially over the results provided by static mean-field electronic structure methods. However GW suffers from some fundamental shortcomings such as, self-screening [START_REF] Nelson | Self-interaction in green'sfunction theory of the hydrogen atom[END_REF][START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF], lack of full self-consistency (in extended systems) [START_REF] Caruso | Unified description of ground and excited states of finite systems: The self-consistent gw approach[END_REF][START_REF] Stan | Levels of selfconsistency in the GW approximation[END_REF][START_REF] Van Schilfgaarde | Quasiparticle self-consistent gw theory[END_REF], lack of size-consistency [START_REF] Erik Dahlen | Variational energy functionals of the green function and of the density tested on molecules[END_REF] and, in particular, it does not describe strong correlation (without imposing a magnetic ordering) [START_REF] Romaniello | Double excitations in finite systems[END_REF]. For example paramagnetic NiO, which is an insulator in nature, is described as a metal by GW . This is because GW , besides exchange, gives a classical description of the system and of its response to an additional electron or hole [START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF]. Moreover, it does not work well for spectral features which are related to other excitations of the system, such as e.g. plasmons [START_REF] Guzzo | Valence electron photoemission spectrum of semiconductors: Ab initio description of multiple satellites[END_REF].

The 2-GF is instead linked to neutral excitations. These excitation energies are obtained by solving the Bethe-Salpeter equation (BSE) [START_REF] Martin | Interacting Electrons: Theory and Computational Approaches[END_REF][START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF]. The BSE is a Dysonlike equation and relies on the 2-GF, specifically a part of it known as the two-body correlation function. The description of a single electron-hole pair is well-described by the most common and widely used approximation to the kernel of the BSE which is based on the GW approximation. This approximation to the BSE kernel is dynamical, meaning it is frequency dependent. However, in practical calculations, it is often considered static, meaning it becomes frequency-independent. Consequently, this treatment excludes the ability to describe double excitations, i.e., the simultaneous formation of two electron-hole pairs [START_REF] Romaniello | Double excitations in finite systems[END_REF][START_REF] Sangalli | Double excitations in correlated systems: A many-body approach[END_REF], and in general higher-order excitations.

In this thesis we will concentrate on the theoretical description of photoemission and optical spectroscopy. These two techniques are valuable tools to get information about the electronic structure, which is the cornerstone of our understanding of material properties.

Photoemission spectroscopy

Direct and inverse photoemission spectroscopy are the main techniques to study the electronic structure of solids. In direct Photoemission Spectroscopy (PES) [START_REF] Hüfner | Photoelectron Spectroscopy: Principles and Applications[END_REF][START_REF] Nordling | Precision method for obtaining absolute values of atomic binding energies[END_REF], high-energy photons, typically X-rays, are directed onto the surface of the material. These photons have enough energy to ionize the electrons of the atoms in the material. When an electron is ejected it creates a hole, namely an unoccupied energy level in the electronic structure. By measuring the kinetic energy of emitted electrons, information about the electronic structure is obtained, such as the valence band structure. Inverse Photoemission spectroscopy (IPES) [START_REF] Binnig | Tunneling spectroscopy and inverse photoemission: Image and field states[END_REF][START_REF] Cardona | Photoemission in Solids[END_REF] exploits the opposite process.

Here, an electron is absorbed into the material causing the emission of a photon. By measuring the photon energy, information about the conduction band structure is obtained.

For materials which are not strongly correlated it is possible to distinguish two kind of peaks in the spectrum: quasiparticles and satellites. Quasiparticle peaks are due to dressed one-particle excitations, i.e., it is possible to describe them as non-interacting particles with a finite lifetime and an effective mass that contains the effect of the interaction with all the particles of the system, such as nuclei and other electrons. Satellite peaks are due to the residual interaction and results in additional excitations. They cannot be treated in any independent-particle picture. These additional excitations of satellites peaks are usually due to the coupling of electrons with pairs of electrons or (quasi)bosons, such as plasmons (see for example in sodium [START_REF] Steiner | Simple metals[END_REF]), phonons, and magnons. At strong correlation this picture breaks down since quasiparticles and satellites mix. This regime is of particular interest since non-trivial or 'exotic' electronic properties of materials emerge and it represents a challenge for theoretical methods.

Figure 1.1 shows a schematic representation of the two techniques.

Absorption

Absorption occurs when a system undergoes a transition from a lower-energy state to a higher-energy state by absorbing energy, normally carried by photons. These photons do not possess enough energy to ionize the electrons, but they are sufficiently energetic to excite an electron from valence to a conduction band, typically in the form of visible, U V or X-rays photons. This process corresponds to the creation of an electron-hole pair. Figure 1.2 schematically represents the absorption process. Other excitations, or satellites, may also appear in the spectra. For example, it is well-known that the strong spatial localization may lead to the simultaneous formation of two eh pairs [START_REF] Sangalli | Double excitations in correlated systems: A many-body approach[END_REF]. For a complete analysis of photoemission spectroscopy and absorption, many different effects must be taken into account. For example, the incoming photon (electron) experiences surface effects before entering the bulk of the material; in the bulk, all electrons interact with each other and the nuclei, generating electron-phonon contributions. Moreover, thermal effects, defects, and impurities will modify the final spectra. In this thesis, we will focus our attention on the electron-electron interaction. Therefore, we will consider systems at zero temperature without defects or impurities. Furthermore, we will restrict the analysis to the Born-Oppenheimer approximation, where electrons interact with the nuclei only through a mean-field potential. Electronphonon interactions are not considered.

Introduction

Outline

The thesis is organized as follows:

• in chapter 2 we give a detailed introduction on the one-body and two-body Green's function. We also discuss some standard approximations that are used to solve the corresponding Dyson equations.

• In chapter 3 we analyze in detail the properties of the three-body Green's function . We show which physical information it contains, focusing in particular on the electron-electron-hole and hole-hole-electron channels. We derive an equation to retrieve the 1-GF from these channels of the 3-GF. We provide a similar analysis for the four-body Green's function.

• In chapter 4 we introduce the main result of this work: the multi-channel Dyson equation (MCDE). It is an equation that links Green's functions with different numbers of particles. First we show the general equation and then we focus on the MCDE that couples the one-body and three-body channels and the MCDE that couples the two-body and four-body channels. We also provide a general approach to obtain approximations to the MCDE self-energy that can be applied in practical calculations.

• In chapter 5 we apply the multi-channel Dyson equation that couples the onebody and three-body channels to the Hubbard dimer at 1/4 and 1/2 filling. We show that we obtain the exact result in both cases.

• In chapter 6 we draw our conclusions.

Background

CHAPTER 2

Background

The electronic Hamiltonian of an N -particle system is

Ĥ = N i=1 ĥ(r i ) + 1 2 N i,j=1 i̸ =j v(r i , r j ), (2.1)
where ĥ contains all the one-body effects, i.e., the kinetic energy and the electronnuclei interaction, and v is the two-body electron-electron interaction. In the secondquantization formalism, the Hamiltonian is

Ĥ = Ĥ0 + Ĥ1 = dx ψ † (x)h(r) ψ(x) + 1 2 dxdx ′ ψ † (x) ψ † (x ′ )v(r, r ′ ) ψ(x ′ ) ψ(x).
(2.2)

Here x means the set of space r and spin σ, and ψ and ψ † are annihilation and creation field operators, respectively. We will use atomic units ℏ = e = m e = 4πϵ 0 = 1 throughout this thesis. The one-body Hamiltonian is given by

h(r) = - 1 2 ∇ 2 + V (r), (2.3) 8
Background where ∇ 2 is the Laplacian operator, and V (r) is the electron-nuclei potential. The two-body term is defined in terms of the Coulomb interaction, i.e.,

v(r, r ′ ) = v(|r -r ′ |) = 1 |r -r ′ | . ( 2.4)
The time-independent Schrodinger equation contains all the information about a given system at equilibrium. It is given by Ĥ|Ψ⟩ = E|Ψ⟩, (2.5) where |Ψ⟩ is the many-body wave function, and E is the energy of the system. An exact solution to this eigenvalue problem is impossible to achieve due to the exponential growth of the calculation with the number of electrons. To solve the problem it is convenient to start from the non-interacting solution. Removing the electron-electron interaction in the Hamiltonian, the eigenvectors |Ψ⟩ reduce to single Slater determinants and the Schrodinger equation can be solved. From the solution of the noninteracting problem, it is possible to add perturbatively the electron-electron Coulomb interaction [START_REF] Chr | Note on an approximation treatment for manyelectron systems[END_REF]. The most used approximation to go beyond the non-interacting picture is the Hartree-Fock method (HF) [START_REF] Szabo | Modern quantum chemistry: introduction to advanced electronic structure theory[END_REF]. It aims to find the Slater determinant of single-electron wave functions that minimizes the total energy of a many-electron system. This is achieved by considering each electron's motion in an average field generated by all other electrons. Due to the complexity of the many-body waves function |Ψ⟩, other strategies to solve the many-body problem are used, for example coupled cluster (CC) [START_REF] Shavitt | Many-Body Methods in Chemistry and Physics: MBPT and Coupled-Cluster Theory[END_REF], configuration interaction (CI), Monte-Carlo (MC) [START_REF] Foulkes | Quantum monte carlo simulations of solids[END_REF], or functional theories. Functional theories are based on quantity simpler than the many-body waves function |Ψ⟩ but still with enough information to completely describe the property of interest, i.e., the density, the 1-density matrix, the one-body Green's function, etc.. In this thesis, we focus our attention on the Green's function theory.

The N-body Green's function

In many-body theory, the Green's function technique is one of the most widely used methods to study particle propagation. In this work, we focus on the non-relativistic, zero temperature, equilibrium case. The time-ordered n-body Green's function (n-GF) is defined as

G n (1, 2, ..., n, 1 ′ , 2 ′ , ..., n ′ ) = (-i) n ⟨Ψ N 0 |T [ ψ(1) ψ(2)... ψ(n) ψ † (n ′ )... ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N 0 ⟩, (2.6)
Background 9 with 1 = (x 1 , t 1 ), x 1 = (r 1 , σ 1 ), |Ψ N 0 ⟩ the ground-state many-body wavefunction of an N -particle system, and T the time ordering operator defined as T [ ψ(1)... ψ(n)] = p (-1) p θ(t p1 > t p2 )...θ(t pn-1 > t pn ) ψ(p 1 )... ψ(p n ), (2.7) where the sum runs over all the possible permutations, and p is the number of permutations with respect to the initial order t 1 ...t n . The field operators evolve in time in the Heisenberg picture ψ(1) = Û † (t 1 , 0) ψ(x 1 ) Û (t 1 , 0) = e i Ĥt1 ψ(x 1 )e -i Ĥt1 , (

where U (t, t ′ ) is the time evolution operator and we have assumed Ĥ time-independent.

While the Heisenberg picture is commonly used to describe time evolution, it is not the only approach. The interaction picture also plays an important role in the Green's function theory. The Hamiltonian can be divided in a one-body Ĥ0 and in a two body Ĥ1 term according to Eq. (2.2). In the interaction picture, the time evolution operator is written as Û (t, 0) = e -i Ĥ0t ÛI (t, 0) (

and the field operator ψI (1) = e i Ĥ0t1 ψ(x 1 )e -i Ĥ0t1 (2.10) evolves in time with Ĥ0 . The field operator in the Heisenberg picture is retrieved with ψ(1) = Û † I (t 1 , 0) ψI (1) ÛI (t 1 , 0). (2.11) Making use of the Gell-Mann and Low theorem [START_REF] Fetter | Quantum Theory of Many-Particle Systems[END_REF], the n-GF can be rewritten as

G n (1, ..., n, 1 ′ , ..., n ′ ) = (-i) n ⟨Ψ N 0 |T [ Ŝ ψI (1) ψI (2)... ψI (n) ψ † I (n ′ )... ψ † I (2 ′ ) ψ † I (1 ′ )]|Ψ N 0 ⟩ ⟨Ψ N 0 |T [ Ŝ]|Ψ N 0 ⟩ , (2.12
) where Ŝ is the so-called scattering operator

Ŝ = e -i ∞ -∞ dt Ĥ1(t) , (2.13)
and Ĥ1 (t) is the time-dependent electron-electron Hamiltonian given by Ĥ1 (t) = 1 2 dxdx ′ ψ † I (x, t) ψ † I (x ′ , t)v(r, r ′ ) ψI (x ′ , t) ψI (x, t).

(2.14)

The definition of the n-GF in Eq. (2.12) is of great interest since it moves all the complexity inside the scattering operator. By expanding the exponential in Eq. (2.13), Background it is possible to analyze the n-GF at each order in the electron-electron Coulomb interaction. Utilizing the Wick theorem [START_REF] Wick | The evaluation of the collision matrix[END_REF] after this expansion led to the development of the Feynman diagrams technique. In sections 2.2.6 and 2.3.4, we will explore it.

Before starting the description of the one-, two-, and three-body Green's functions, we report the Martin-Schwinger (MS) hierarchy [START_REF] Martin | Theory of many-particle systems[END_REF]. It is an important relation that links a general n-GF with the n±1-GF as follows

i d dt k -ĥ(k) G n (1, ..., n, 1 ′ , ..., n ′ ) = n j=1
(-1) j+k δ(k, j ′ )

G n-1 (1, ..., k, ..., n, 1 ′ , ..., j′ , ..., n ′ ) -i d 1v(k, 1)G n+1 (1, ..., n, 1+ , 1 ′ , ...n ′ , 1++ ).

(2.15)

Here, the indices with a tilde are removed from the Green's function, 1 + = (x 1 , t + 1 ) = (x 1 , t 1 + η) with η → 0 + , and we have utilized the notation

v(1, 2) = v(r 1 , r 2 )δ(t 1 -t 2 ).
(2.16)

Through the Martin-Schwinger hierarchy, it is possible to link the 1-GF to higherorders GFs. However, dealing with high-order Green's functions can become exceedingly complex and computationally demanding. Therefore, we will adopt wellestablished strategies that allow us to systematically eliminate higher-order Green's functions from the MS equation (2.15) thanks to the definition of a functional called self-energy [START_REF] Csanak | Green's function technique in atomic and molecular physics[END_REF][START_REF] Strinati | Application of the green's functions method to the study of the optical properties of semiconductors[END_REF].

The one-body Green's function

The one-body Green's function is the most widely used propagator. Its importance is due to the simple link with observables of interest, such as the density, the total energy, and photoemission spectra.

Following Eq. (2.6), it is defined as

G 1 (1, 1 ′ ) = -i⟨Ψ N 0 |T [ ψ(1) ψ † (1 ′ )]|Ψ N 0 ⟩. (2.17)
The time-ordered product can be written explicitly as

T [ ψ(1) ψ † (1 ′ )] = ψ(1) ψ † (1 ′ ) for t 1 > t 1 ′ -ψ † (1 ′ ) ψ(1) for t 1 ′ > t 1 , (2.18)
from which

G 1 (1, 1 ′ ) = -iΘ(t 1 -t 1 ′ )⟨Ψ N 0 | ψ(1) ψ † (1 ′ )|Ψ N 0 ⟩ + iΘ(t 1 ′ -t 1 )⟨Ψ N 0 | ψ † (1 ′ ) ψ(1)|Ψ N 0 ⟩, (2.19)
where Θ(t) is the Heaviside step function. The first term on the right-hand side can be interpreted as the probability amplitude for an electron to be added to an N -electron system at time t 1 ′ at position r 1 ′ with spin σ 1 ′ , and to be found at position r 1 with spin σ 1 at a later time t 1 > t 1 ′ . The other term describes the probability amplitude for the propagation of a hole. For this reason, it is convenient to define

G e 1 (1, 1 ′ ) = -iΘ(t 1 -t 1 ′ )⟨Ψ N 0 | ψ(1) ψ † (1 ′ )|Ψ N 0 ⟩, (2.20)
G h 1 (1, 1 ′ ) = iΘ(t 1 ′ -t 1 )⟨Ψ N 0 | ψ † (1 ′ ) ψ(1)|Ψ N 0 ⟩. (2.21)
The hole propagator G h is linked to direct photoemission spectroscopy, where an electron is removed from the system, which is equivalent to adding a hole, while, the electron propagator G e describes inverse photoemission spectroscopy, where an electron is added to the system. It is interesting to start the analysis from the simplest case: the non-interacting 1-GF.

Non-interacting Green's function

In the non-interacting regime, the Hamiltonian reduces to the one-body part, i.e., Ĥ0 . With a canonical transformation of the field operators

ψ(x, t) = p ϕ p (x)ĉ p (t) ψ † (x, t) = p ϕ * p (x)ĉ † p (t), (2.22)
where ϕ p (x) is s basis set of one-electron spinorbitals, the Hamiltonian becomes

Ĥ0 = p ϵ p ĉ † p ĉp , (2.23)
and the creation and annihilation operators evolve in time according to

ĉp (t) = e -iϵpt ĉp ĉ † p (t) = e iϵpt ĉ † p (2.24)
with the condition ĉp (0) = ĉp . The energies of the occupied states are below the Fermi level. Consequently, the expectation value of a creation and an annihilation operator on the ground state of a non-interacting system is

⟨Ψ N 0 |ĉ † q ĉp |Ψ N 0 ⟩ = δ pq f p ⟨Ψ N 0 |ĉ p ĉ † q |Ψ N 0 ⟩ = δ pq (1 -f p ) (2.25)
where

f p = 0 if ϵ p > E F 1 if ϵ p < E F , (2.26)
Background where E F is the Fermi energy. With these relations, the non-interacting Green's function is easily found to be

G 0 1 (1, 1 ′ ) = -i p e -iϵp(t1-t 1 ′ ) ϕ p (x 1 )ϕ * p (x 1 ′ )[θ(t 1 -t 1 ′ )(1 -f p ) -θ(t 1 ′ -t 1 )f p ].
(2.27)

The G 0 1 can hence be written as

G 0 1 (1, 1 ′ ) = im G 0 im (t 1 -t 1 ′ )ϕ i (x 1 )ϕ * m (x 1 ′ ) (2.28)
where G 0 im is given by (2.29) This expression allows us to find a very useful relation that we will use later in this thesis. Let us start the analysis with the electron contribution

G 0 im (t 1 -t 1 ′ ) = G 0(e) im (t 1 -t 1 ′ ) + G 0(h) im (t 1 -t 1 ′ ) = -i p e -iϵp(t1-t 1 ′ ) [θ(t 1 -t 1 ′ )(1 -f p )δ ip δ pm -θ(t 1 ′ -t 1 )f p δ mp δ pi ].
G 0(e) im (t 1 -t 1 ′ ) = -i p θ(t 1 -t 1 ′ )(1 -f p )e -iϵp(t1-t 1 ′ ) δ ip δ pm = -i pql θ(t 1 -t ′′ + t ′′ -t 1 ′ )(1 -f p )(1 -f q )e -iϵp(t1-t ′′ ) e -iϵq(t ′′ -t 1 ′ ) × δ ip δ pl δ lq δ qm = i l G 0(e) il (t 1 -t ′′ )G 0(e) lm (t ′′ -t 1 ′ ) (2.30)
where we have added two completeness relations I = a |a⟩⟨a|, and used the relation

(1 -f p ) = (1 -f p ) 2 . Finally, we applied the restriction t 1 > t ′′ > t 1 ′ to simplify the relation θ(t 1 -t ′′ + t ′′ -t 1 ′ ) = θ(t 1 -t ′′ )θ(t ′′ -t 1 ′
). Following similar steps but for the hole part, we obtain

G 0(h) im (t 1 -t 1 ′ ) = -i l G 0(h) il (t 1 -t ′′ )G 0(h) lm (t ′′ -t 1 ′ ) (2.31)
with the condition t 1 ′ > t ′′ > t 1 . Notice that merging two non-interacting electron or hole Green's functions gives a difference in sign. Considering the electron and hole part together

G 0 im (t 1 -t 1 ′ ) = i l G 0(e) il (t 1 -t ′′ )G 0(e) lm (t ′′ -t 1 ′ ) -i l G 0(h) il (t 1 -t ′′ )G 0(h) lm (t ′′ -t 1 ′ ).
(2.32)

Background 13
It is interesting for subsequent analysis to report the relation in position space

G 0 1 (1, 1 ′ ) = i dyG 0(e) 1 (x, y, t 1 -t ′′ )G 0(e) 1 (y, x ′ , t ′′ -t 1 ′ ) -i dyG 0(h) 1 (x, y, t 1 -t ′′ )G 0(h) 1 (y, x ′ , t ′′ -t 1 ′ ). (2.33)
These last relations are also valid for a mean field theory, such as Hartree-Fock. Before moving on to describe the interacting theory, it is instructive to see how the Martin-Schwinger hierarchy reduces when applied to the non-interacting one-body Green's function. Since v(1, 2) = 0 and G n-1 = I when n = 1, we simply have

i d dt 1 -ĥ(1) G 0 1 (1, 1 ′ ) = δ(1, 1 ′ ). (2.34)
From the definition of the inverse

d2[G 0 1 ] -1 (1, 2)G 0 1 (2, 1 ′ ) = d2G 0 1 (1, 2)[G 0 1 ] -1 (2, 1 ′ ) = δ(1, 1 ′ ) (2.35)
we obtain the well-known relation

[G 0 1 ] -1 (1, 2) = i d dt 1 -ĥ(1) δ(1, 2).
(2.36)

The Lehmann representation

In this section, we will rewrite the 1-GF to highlight its connection with the charged excitation energies of the system. This representation is commonly known as the Lehmann or spectral representation. By introducing the resolution of the identity [START_REF] Viktor | Application of quantum field theory methods to the many body problem[END_REF], where |Ψ M k ⟩ indicates the k-th eigenstate of the M -electron system, we obtain

M k |Ψ M k ⟩⟨Ψ M k | = I in Eq. (2.
G 1 (1, 1 ′ ) = -iΘ(t 1 -t 1 ′ ) m ⟨Ψ N 0 | ψ(1)|Ψ N +1 m ⟩⟨Ψ N +1 m | ψ † (1 ′ ) Ψ N 0 ⟩ + iΘ(t 1 ′ -t 1 ) n ⟨Ψ N 0 | ψ † (1 ′ )|Ψ N -1 n ⟩⟨Ψ N -1 n | ψ(1) Ψ N 0 ⟩, (2.37)
where we have used the fact that field operators only select states with N +1 particles if t 1 > t

Background

We thus obtain the following expression for the 1-GF

G 1 (x 1 , x 1 ′ ; τ ) = -iΘ(τ ) m ⟨Ψ N 0 | ψ(x 1 )|Ψ N +1 m ⟩⟨Ψ N +1 m | ψ † (x 1 ′ )|Ψ N 0 ⟩e i(E N 0 -E N +1 m )τ + iΘ(-τ ) n ⟨Ψ N 0 | ψ † (x 1 ′ )|Ψ N -1 n ⟩⟨Ψ N -1 n | ψ(x 1 ) Ψ N 0 ⟩e -i(E N 0 -E N -1 n )τ .
(2.39)

where τ = t 1 -t 1 ′ . We notice that the 1-GF depends only on the time difference τ . This is due to the system's time-independent Hamiltonian, resulting in a ground state that is not affected by the time of particle addition. What truly matters is the duration of the particle's propagation, described by the time difference τ . For a more compact notation, we define

f m (x 1 ) = ⟨Ψ N 0 | ψ(x 1 )|Ψ N +1 m ⟩, (2.40)
and

g n (x 1 ) = ⟨Ψ N -1 n | ψ(x 1 ) Ψ N 0 ⟩. (2.41)
They are generally referred to as the Feynman-Dyson amplitudes. Therefore, Eq. (2.39) takes the form

G 1 (x 1 , x 1 ′ ; τ ) = -iΘ(τ ) m f m (x 1 )f * m (x 1 ′ )e i(E N 0 -E N +1 m )τ + iΘ(-τ ) n g n (x 1 )g * n (x 1 ′ )e -i(E N 0 -E N -1 n )τ . (2.42)
Note that the amplitudes f m as well as g n are, in general, not orthogonal nor linearly independent. One can only prove that the total set is complete:

m f m (x 1 )f * m (x 1 ′ ) + n g n (x 1 )g * n (x 1 ′ ) = δ(x 1 -x 1 ′ ). (2.43)
Using the Fourier transform of the Heaviside step function

∞ -∞ dt Θ(±t)e -iαt e iωt = lim η→0 + ±i ω -α ± iη , (2.44)
we can obtain the Fourier transform of Eq. (2.42). We arrive at the Lehmann representation in frequency space for the 1-GF a

G 1 (x, x ′ ; ω) = lim η→0 + m f m (x)f * m (x ′ ) ω -(E N +1 m -E N 0 ) + iη + n g n (x)g * n (x ′ ) ω -(E N 0 -E N -1 n ) -iη .
(2.45)

a we use the following definition of Fourier transform of a general function f : f (ω) = dtf (t)e iωt ; and for the reverse Fourier transform: ) energies of the system. Figure 2.1 depicts the poles of the 1-GF in the complex frequency space.

f (t) = 1/(2π) dω f (ω)e -iωt .
We observe that in the case of a non-interacting system or more generally with wave functions described by single Slater determinants, the 1-GF can be expressed as follows

G 0 1 (x 1 , x 1 ′ ; ω) = i ϕ i (x 1 )ϕ * i (x 1 ′ ) ω -ϵ 0 i + sign(ϵ 0 i -µ)iη . (2.46)
With the change of basis given in Eq. (2.22), it becomes

G 0 ij (ω) = δ ij ω -ϵ 0 i + sign(ϵ 0 i -µ)iη . (2.47)
Usually the non interacting energies and wave-functions are replaced by the Kohn-Sham ones obtained from DFT calculations.

Physical properties

The 1-GF is of great interest due to its straightforward connection to several essential properties of the system. It allows the calculation of:

1. the ground state expectation value of any single-particle operator; Background 2. the ground state energy 3. the one-particle charged excitation energies of the system.

Let us discuss these points in detail.

In the second-quantization formalism, any one-particle operator Ô = i o(x i ) is given by Ô = dx ψ † (x)o(x) ψ(x).

(2.48)

The ground-state expectation value of Ô can thus be evaluated from the 1-GF according to

⟨Ψ N 0 | Ô |Ψ N 0 ⟩ = dx lim x ′ →x o(x)⟨Ψ N 0 | ψ † (x ′ ) ψ(x) |Ψ N 0 ⟩ = -i dx lim x ′ →x lim t ′ →t + o(x)G 1 (x, x ′ ; t -t ′ ). (2.49)
Thus the 1-GF allows us to calculate the ground-state expectation value of any oneparticle operator. One significant example of it is the density operator. It is defined as

ρ(r) = N i=1 δ(r -r i ), (2.50)
using Eq. (2.49) we obtain

ρ(r) = -i lim τ →0 - G 1 (r, r; τ ). (2.51)
The link to the ground-state energy is described by the Galitskii-Migdal formula [START_REF] Viktor | Application of quantum field theory methods to the many body problem[END_REF]

E = ⟨Ψ N 0 | Ĥ |Ψ N 0 ⟩ = - i 2 dx 1 lim t ′ 1 →t + 1 lim x ′ 1 →x1 i ∂ ∂t 1 + h(x 1 ) G 1 (1, 1 ′ ).
(2.52)

The one-particle excitation spectrum can be obtained from the spectral function which is defined as

A(x 1 , x 1 ′ ; ω) = 1 π sign(µ -ω)ImG 1 (x 1 , x 1 ′ ; ω). (2.53)

The self-energy and the Dyson equation

So far, we have discussed what can be calculated using the 1-GF. However, the equations presented so far do not provide a direct method to compute the 1-GF. In this section, we address the challenge of finding a way to calculate the 1-GF. The main idea is to start with what is possible to calculate, namely the non-interacting 1-GF, and construct a closed equation to add interaction to it. All the information related to interactions is encapsulated in an operator called self-energy. This quantity is intricate, being non-local and frequency-dependent, and it is not possible to calculate it exactly (except for some simple model systems). Hence, approximations are required.

In the following, we will explore some common approximations and their limitations, especially concerning strongly correlated materials.

To define the self-energy, we write the Martin-Schwinger hierarchy for the 1-GF

i d dt k -ĥ(k) G 1 (1, 1 ′ ) = δ(1, 1 ′ ) -i d2v(1, 2)G 2 (1, 2 + , 1 ′ , 2 ++ ). (2.54)
Applying the identity

δ(1, 1 ′ ) = d2G -1 1 (1, 2)G 1 (2, 1 ′ ) (2.55)
on the right-hand side of Eq. (2.54), and using Eqs. (2.35) and (2.36) we obtain the so-called Dyson equation

G 1 (1, 1 ′ ) = G 0 1 (1, 1 ′ ) + d2d2 ′ G 0 1 (1, 2)Σ(2, 2 ′ )G 1 (2 ′ , 1 ′ ), (2.56)
where the self-energy Σ is defined as

Σ(1, 1 ′ ) = -i d2d2 ′ v(1, 2)G 2 (1, 2 + , 2 ′ , 2 ++ )G -1 1 (2 ′ , 1 ′ ). ( 2

.57)

A simple approximation to the self-energy is the Hartree-Fock approximation. Being a mean-field approximation, it is sufficient to consider an independent particle picture to derive it. This implies that the 2-GF in Eq. (2.57) is simply defined as the antisymmetric product of two 1-GFs, without considering interaction between them, namely

G 2 (1, 2, 1 ′ , 2 ′ ) = G 1 (1, 1 ′ )G 1 (2, 2 ′ ) -G 1 (1, 2 ′ )G 1 (2, 1 ′ ). (2.58)
We thus obtain the HF self-energy,

Σ HF (1, 1 ′ ) = δ(1, 1 ′ ) d2v(1, 2)ρ(2) + iv(1, 1 ′ )G 1 (1, 1 ′ ) = δ(1, 1 ′ )v H (1) + iv(1, 1 ′ )G 1 (1, 1 ′ ). (2.59)
The first term on the right-hand side of Eq. (2.59) describes the local interaction with the electron density, i.e., the Hartree term, while the second term is non-local in the interaction, i.e., the Fock or exchange term.

Background

Making use of the frequency dependence of the 1-GF, as given in Eq. (2.45), it becomes possible to express the Dyson equation as a function of the frequency

G 1 (x 1 , x 1 ′ ; ω) = G 0 1 (x 1 , x 1 ′ ; ω) + d2d2 ′ G 0 1 (x 1 , x 2 ; ω)Σ(x 2 , x 2 ′ ; ω)G 1 (x 2 ′ , x 1 ′ ; ω).
(2.60) Applying the inverse of G from the right and the inverse of G 0 from the left, the Dyson equation becomes

G -1 1 (x 1 , x 1 ′ ; ω) = [G 0 1 ] -1 (x 1 , x 1 ′ ; ω) -Σ(x 1 , x 1 ′ ; ω). (2.61)
The above equation is interesting because it shows more clearly that a static selfenergy just shift the poles of G 0 1 , whereas a dynamical self-energy also creates extra poles.

It is instructive to analyse the 1-GF in the diagonal approximation, i.e., G ij = 0 if i ̸ = j. In this approximation, the 1-GF projected in a basis takes the form

G ii (ω) = ω -ϵ 0 i -Σ ii (ω) -1 . (2.62)
where ϵ 0 i is the eigenvalue of the non-interacting system. Writing the spectral function (2.53) in the diagonal approximation yields

A ii (ω) = 1 π |ImΣ ii (ω)| [ω -ϵ 0 i -ReΣ ii (ω)] 2 + ImΣ ii (ω) 2 , (2.63)
where ReΣ and ImΣ indicates the real and imaginary parts, respectively, of the selfenergy. This representation allows us to better understand the role of the self-energy. The imaginary part of Σ is associated with the broadening of the Lorentzian peak. Physically, this quantity can be interpreted as the inverse of the particle's lifetime. On the other hand, the real part of Σ shifts the non-interacting peaks according to the interactions present in the system. If the peak remains sharp, namely ImΣ → 0, the zeros of the real part are at

ϵ i = ϵ 0 i + ReΣ ii (ω = ϵ 0 i ), (2.64)
in this picture, there is a shift in the peak's position, but it maintains a one-toone correspondence with the non-interacting particle peak. In this case, we are in the so-called quasi-particle limit [START_REF] Nozieres | Theory Of Quantum Liquids[END_REF][START_REF] Pines | Elementary Excitations in Solids[END_REF]. Another approach to achieve the quasiparticle representation is to consider the static self-energy. This limit is well known in literature and gives good results for materials where quasi-particles are stable. We refer to this systems as weakly correlated materials. If the broadening is significant, i.e., when the ImΣ is large, the one-particle picture is no longer correct for interpreting the system. The presence of this additional significant term in the denominator leads to the creation of extra peaks in the spectral function. These additional peaks are known as satellites. As the spectral function is normalized, these new structures take weight away from the quasi-particle peak. When the satellite structures become relevant in the total spectra, we enter the strongly correlated regime. Studying this scenario becomes more challenging.

Hedin's equations and GW approximation

In 1965, Lars Hedin formally derived a closed set of five equations to calculate the 1-GF [START_REF] Hedin | New method for calculating the one-particle green's function with application to the electron-gas problem[END_REF]. They are defined as

G 1 (1, 1 ′ ) = G 0 1 (1, 1 ′ ) + d2d2 ′ G 0 1 (1, 2)Σ(2, 2 ′ )G 1 (2 ′ , 1 ′ ), (2.65) Σ(1, 2) = v H (1)δ(1, 2) + i d3d4G 1 (1, 4)W (3, 1) Γ(4, 2; 3), (2.66) W (1, 2) = v(1, 2) + d3d4v(1, 3) χ(3, 4)W (4, 2), (2.67)
χ(1, 2) = -i d3d4G 1 (1, 3)G 1 (4, 1) Γ(3, 4; 2), (2.68)
Γ(1, 2; 3) =δ(1, 2)δ(2, 3) + d4d5d6d7 δΣ xc (1, 2) δG 1 (4, 5) G 1 (4, 6)G 1 (7, 5) Γ(6, 7; 3), (2.69)

Here, v H is the Hartree potential, W is the dynamically screened Coulomb potential, χ is the irreducible polarizability, and Γ is the irreducible vertex function. Equation (2.66) defines the exchange-correlation part of the self-energy, denoted as Σ xc , which excludes the Hartree term. It is worth noting that the first four equations are integral equations, while the last equation is an integro-differential equation. Solving the last equation involves a functional derivative of the exchange-correlation part of the selfenergy with respect to the 1-GF, which is related to the 2-GF, as we will discuss in more detail in the next chapter. This makes the last equation more complex and challenging. As a first approximation, it is common to consider the vertex function at zero order in the potential, thus simplifying the equations. By setting Γ(1, 2; 3) = δ(1, 2)δ(2, 3) the exchange-correlation part of the self-energy takes the form On the left there is the selfconsistent procedure that is difficult to do in practical calculation. On the right there is the Hedin pentagon in the GW approximation. In this later case, the equations to solve reduce to four, and accordingly the pentagon reduces to a quadrangle.

Σ GW xc (1, 2) = iG(1, 2)W (2, 1), ( 2
where the superscript RP A stands for random-phase approximation. In the RPA, the irreducible polarizability is simplified to the product of two Green's functions, capturing the electron-electron interaction effects in a mean-field-like manner.

This approximation is commonly referred to as the GW approximation. Notice that Σ GW xc is very similar to the exchange self-energy in Eq. (2.59), the only difference is that the potential in GW is dynamically screened.

Hedin's equations have to be solved self-consistently, as shown in Figure 2.2.

In the GW scheme, one starts from the non-interacting 1-GF to calculate the irreducible polarizability (Eq. (2.71)), and from it the screened Coulomb potential (Eq. (2.67)). From the screened Coulomb potential, one then computes the selfenergy (Eq. (2.70)). Finally, with the updated self-energy, the 1-GF is calculated (Eq. (2.65)). These steps are repeated iteratively until convergence is achieved.

However, performing all these calculations for real systems is computationally demanding. As a result, a common approach is to start from the best available 1-GF, often obtained from Hartree-Fock or Kohn-Sham calculations, and perform the steps of the Hedin pentagon only once. The obtained self-energy is then used to solve the Dyson equation. This approach is called one-shot GW or G 0 W 0 .

The frequency dependence of the self-energy introduces further complexity to the calculations. To address the computational challenges, a common strategy is to employ a static approximation for the self-energy. This strategy allows for a good description of quasiparticle peaks. However, the static approximation completely misses the satellite structure.

Feynman diagram analysis

In this section, we will provide a general introduction to the Feynman diagram technique. It is a way to visualize long equations, and we will use it in section 4. First, we will explore the creation of Feynman diagrams from the definition of the 1-GF. Next, we will examine the general properties of these diagrams. Finally, we will show how the Hedin equations are represented diagrammatically.

We begin by defining the 1-GF using equation (2.12)

G 1 (1, 1 ′ ) = -i ⟨Ψ N 0 |T [ Ŝ ψI (1) ψ † I (1 ′ )]|Ψ N 0 ⟩ ⟨Ψ N 0 |T [ Ŝ]|Ψ N 0 ⟩ . (2.72)
Next, we expand the scattering operator:

Ŝ = I - i 2 d1d2 ψ † I (1) ψ † I (2)v(1, 2) ψI (2) ψI (1) - 1 4 d1d2d3d4 ψ † I (1) ψ † I (2)v(1, 2) ψI (2) ψI (1) ψ † I (3) ψ † I (4)v(3, 4) ψI (4) ψI (3) + ... (2.73)
where 1 = (x 1 , t 1 ) and v(1, 2) given in Eq. (2.16). At order n, the prefactor is (-i) n 2 n , with 2n creation and 2n annihilation operators, and n bare Coulomb interactions.

Let us start the analysis by considering the first order of the numerator (N ) in Eq (2.72). It reads

N (1, 1 ′ ) = - 1 2 d2d3v(2, 3)⟨Ψ N 0 |T [ ψ † I (2) ψ † I (3) ψI (3) ψI (2) ψI (1) ψ † I (1 ′ )]|Ψ N 0 ⟩. (2.74)
To obtain the Feynman diagrams, we will use Wick's theorem, which allows us to couple the six field operators in all possible ways. For instance

⟨Ψ N 0 |T [ ψ † I (2) ψ † I (3) ψI (3) ψI (2) ψI (1) ψ † I (1 ′ )]|Ψ N 0 ⟩ = -⟨Ψ N 0 |T [ ψI (3) ψ † I (3)]|Ψ N 0 ⟩⟨Ψ N 0 |T [ ψI (1 ′ ) ψ † I (2)]|Ψ N 0 ⟩⟨Ψ N 0 |T [ ψI (2) ψ † I (1 ′ )]|Ψ N 0 ⟩ +⟨Ψ N 0 |T [ ψI (1) ψ † I (2)]|Ψ N 0 ⟩⟨Ψ N 0 |T [ ψI (2) ψ † I (3)]|Ψ N 0 ⟩⟨Ψ N 0 |T [ ψI (3) ψ † I (1 ′ )]|Ψ N 0 ⟩ -⟨Ψ N 0 |T [ ψI (2) ψ † I (3)]|Ψ N 0 ⟩⟨Ψ N 0 |T [ ψI (3) ψ † I (2)]|Ψ N 0 ⟩⟨Ψ N 0 |T [ ψI (1) ψ † I (1 ′ )]|Ψ N 0 ⟩ + ... =iG 0 1 (3, 3)G 0 1 (1, 2)G 0 1 (2, 1 ′ ) -iG 0 1 (1, 2)G 0 1 (2, 3)G 0 1 (3, 1 ′ ) +iG 0 1 (2, 3)G 0 1 (3, 2)G 0 1 (1, 1 ′ ) + ... (2.75)
By inserting these combinations into Eq. (2.74), we obtain the first-order Feynman diagrams, depicted in Fig. 2.3. To interpret the diagrams, we associate a non-interacting 1-GF with each arrow line, an interaction v with each wiggly line, and we remember that each internal point is integrated over space and time and summed over spin. The first term contains G 0 (3, 3) = G 0 (3, 3 + ) which represents the interaction of the added particle with the density of a non-interacting system (see Eq. (2.51)). This term corresponds to the Hartree contribution. The second term represents the first-order non-local term, which can be identified as the Fock term. To achieve the self-consistent Hartree-Fock solution, one must use the corresponding self-energy, i.e., Eq. (2.59). From this diagrammatic analysis, the HF self-energy is obtained by removing the two external G 0 lines, since G = G 0 + G 0 Σ HF G 0 + .... With this self-energy, the Dyson equation (2.56) has to be solved iteratively until convergence.

The last term in Fig. 2.3 is a vacuum diagram, meaning that there are fermion lines disconnected from the others. In this example, the line that connects 1 and 1 ′ is disconnected from the other lines. It can be shown (see, for example, [START_REF] Fetter | Quantum Theory of Many-Particle Systems[END_REF]) that such vacuum diagrams do not contribute to the 1-GF since they cancel out when divided by the denominator of Eq.(2.72). Henceforth, whenever we mention Feynman diagrams, we implicitly refer only to connected diagrams.

The other possible Wick contractions in Eq. (2.74) give the same connected result as in Fig. (2.3), since the diagrams have the same topology and the internal indices are integrated. In total, there are two combinations that give the HF terms. The sum of these terms cancels the 1/2 in front of Eq. (2.74).

To obtain higher-order Feynman diagrams, one has to follow the same procedure at each order of the expansion of the scattering operator Eq. (2.67). To simplify the process, Feynman introduced rules to create diagrams at order n with the correct sign and prefactor [START_REF] Fetter | Quantum Theory of Many-Particle Systems[END_REF]. They are the following ones:

1. Draw all topologically distinct connected diagrams with n interaction lines v and 2n + 1 non-interacting 1-GFs;

2. Affix a prefactor (-1) F i n to each diagram, where F is the number of closed fermion loops;

3. Treat a Green's function with equal time variable as G 0 (1, 1 + ) (as we have done above for the Hartree contribution).

Since in practice the 1-GF is obtained through the Dyson equation (2.56), where an approximation to the self-energy is needed, the diagrammatic analysis often focuses directly on the self-energy. Two main definitions of the self-energy are particularly important 1. proper self-energy: it is the sum of all self-energy diagrams that cannot be disconnected into two separate self-energy contributions by removing one G 0 line.

2. skeleton self-energy: it is a self-energy diagram that contains no other self-energy insertions except for itself. Having a proper self-energy is crucial when solving the Dyson equation. Using an improper self-energy leads to a double-counting problem, where certain diagrams enter the 1-GF more times than they should. This issue becomes evident when examining the Dyson equation order by order, namely

G = G 0 + G 0 ΣG 0 + G 0 ΣG 0 ΣG 0 + ....
The concept of a skeleton self-energy is significant in simplifying the understanding of diagrams, as it allows us to focus on the essential contributions while neglecting some higher-order effects. To illustrate this point, let us consider the GW self-energy Σ = iGW . It relies on the 1-GF both directly and through the irreducible polarizability present in W , as given in equations (2.67) and (2.68). When considering this 1-GF as non-interacting, only skeleton diagrams are present. Non-skeleton diagrams Background are obtained by improving the 1-GF within the self-energy. This occurs by solving the Hedin pentagon self-consistently. In contrast, more skeleton diagrams cannot be obtained in this manner; they can only be obtained by improving the self-energy. For instance, in the case of the Hedin equations in the GW approximation, new skeleton diagrams are obtained by adding vertex corrections to the self-energy.

The diagrammatic representation of Hedin's equations is shown in Fig. 2.5. When considering the GW approximation, the diagrams for the self-energy and irredicible polarizability become those depicted in Fig. 2.6. In this approximation, the diagrams for the screened Coulomb potential are sometimes referred to as bubble diagrams, because of their shape, as shown in Fig. 2.7.

The two-body Green's function

The object of interest in this section is the 2-GF. In the context of electronic structure, the 2-GF is widely used to study absorption spectra and it has the task to efficiently describe the excitonic effects by improving the RPA irreducible polarizability. It is also used to describe particle-particle interactions, which are relevant, for example, in the study of superconductivity, where two electrons interact with each other through a phonon, forming the so-called Cooper pairs [START_REF] Bardeen | Theory of superconductivity[END_REF].

General properties

The two-body Green's function is defined as

G 2 (1, 2, 1 ′ , 2 ′ ) = -⟨Ψ N 0 |T [ ψ(1) ψ(2) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N 0 ⟩. (2.76)
Adapting Eq. (2.12) to the two-body case and making use of Wick's theorem, it is possible to show (see, for example Ref. [START_REF] Fetter | Quantum Theory of Many-Particle Systems[END_REF]) that the 2-GF can be decomposed as

G 2 (1, 2, 1 ′ , 2 ′ ) = G 1 (1, 1 ′ )G 1 (2, 2 ′ ) -G 1 (1, 2 ′ )G 1 (2, 1 ′ ) + δG 2 (1, 2, 1 ′ , 2 ′ ), (2.77)
where the two combinations of GG describe the independent propagation of two dressed particles, while the δG 2 terms include all the ways in which these two particles interact with each other. In the context of the 2-GF, we refer to the terms GG as non-interacting terms. Here, "interaction" specifically denotes the interaction between the two particles.

To begin our analysis, we will focus our attention on the time dependency. The 2-GF given in Eq. (2.76) depends on four distinct times, which can be reduced to three time differences (if the Hamiltonian is time independent, which is the case of this theses). In the majority of cases, it is used to study the propagation of two particles (electrons or holes) or one electron and one hole. As a result, only the time required to describe this propagation is relevant. The other two time differences account for the creation and annihilation of the electron-hole pair or the particle pair. In most cases, these processes are considered to be instantaneous, and wherefore, these times are set to zero.

The eh 2-GF

To describe the electron-hole contribution, we set the times t 1 ′ = t + 1 and t 2 ′ = t + 2 in Eq. (2.76). The eh 2-GF hence becomes

G eh 2 (1, 2, 1 ′ , 2 ′ ) = -⟨Ψ N 0 |T [( ψ(x 1 ) ψ † (x 1 ′ )) t1 ( ψ(x 2 ) ψ † (x 2 ′ )) t2 ]|Ψ N 0 ⟩, (2.78)
where ( ψ(x 1 ) ψ † (x 1 ′ )) t1 means that the two field operators are set at the same time t 1 .

It is important to note that the other choice of ordering the field operators, namely

( ψ(x 1 ) ψ † (x 2 ′ )) t1 ( ψ(x 2 ) ψ † (x 1 ′ )) t2
, is also valid. This choice effectively describes the propagation of an electron-hole pair, as evident from its spectral representation. The calculation of the electron-hole channel is well-established and it can be found in literature (see, e.g., Refs. [START_REF] Csanak | Green's function technique in atomic and molecular physics[END_REF][START_REF] Strinati | Application of the green's functions method to the study of the optical properties of semiconductors[END_REF]). We will present here only the final result. In frequency space, it reads

G eh 2 (x 1 , x 2 , x 1 ′ , x 2 ′ ; ω) = -i lim η→0 + n χ n (x 1 , x 1 ′ ) χn (x 2 , x 2 ′ ) ω -(E N n -E N 0 ) + iη - χn (x 1 , x 1 ′ )χ n (x 2 , x 2 ′ ) ω + (E N n -E N 0 ) -iη (2.79) Background 27
where the electron-hole amplitudes have been defined as

χ n (x 1 , x 1 ′ ) = ⟨Ψ N 0 | ψ(x 1 ) ψ † (x 1 ′ )|Ψ N n ⟩, (2.80) χn (x 1 , x 1 ′ ) = ⟨Ψ N n | ψ(x 1 ) ψ † (x 1 ′ )|Ψ N 0 ⟩. (2.81)
The poles of Eq. (2.79) clearly describe neutral excitation (or deexcitation) energies.

Considering the non-interacting term, i.e.,

G 0 2 (1, 2, 1 ′ , 2 ′ ) = G 1 (1, 1 ′ )G 1 (2, 2 ′ ) - G 1 (1, 2 ′ )G 1 (2, 1 ′
) with this time choice, only one term depends on time since

G 1 (1, 1 ′ )G 1 (2, 2 ′ ) = G 1 (x 1 , x 1 ′ ; 0 + )G 1 (x 2 , x 2 ′ ; 0 + ), (2.82) G 1 (1, 2 ′ )G 1 (2, 1 ′ ) = G 1 (x 1 , x 2 ′ ; t 1 -t 2 )G 1 (x 2 , x 1 ′ ; t 2 -t 1 ). (2.83)
Therefore, only the spectral representation of Eq. ( 2.83) has poles. It reads as follows (considering for simplicity the non-interacting 1-GF (2.46))

[G 0 1 (x 1 , x 2 ′ )G 0 1 (x 2 , x 1 ′ )](ω) = = -i v,c ψ c (x 1 )ψ * c (x 2 ′ )ψ v (x 2 )ψ * v (x 1 ′ ) ω -(ϵ 0 c -ϵ 0 v ) + iη - ψ v (x 1 )ψ * v (x 2 ′ )ψ c (x 2 )ψ * c (x 1 ′ ) ω + (ϵ 0 c -ϵ 0 v ) -iη , (2.84)
where v (c) stands for valence (conduction) and the notation [G 0 1 G 0 1 ](ω) is used for convolution in frequency space. This combination has poles at the neutral excitation energies. For this reason, to study the eh channel, the so-called two-particle correlation function is defined

L(1, 2, 1 ′ , 2 ′ ) = -G eh 2 (1, 2, 1 ′ , 2 ′ ) + G 1 (1, 1 ′ )G 1 (2, 2 ′ ), (2.85)
in which the static combination of G 1 G 1 in Eq.(2.82) is removed from G eh 2 . This corresponds to remove the term n = 0 in Eq. (2.79)

The pp 2-GF

To describe the propagation of two particles, the times in Eq. (2.76) are chosen as t 2 = t + 1 and t 1 ′ = t + 2 ′ , and the particle-particle 2-GF is given by

G pp 2 (1, 2, 1 ′ , 2 ′ ) = ⟨Ψ N 0 |T [( ψ(x 1 ) ψ(x 2 )) t1 ( ψ † (x 1 ′ ) ψ † (x 2 ′ )) t 1 ′ ]|Ψ N 0 ⟩. (2.86)
The spectral representation is

G pp 2 (x 1 , x 2 , x 1 ′ , x 2 ′ ; ω) = i lim η→0 + n ζ n (x 1 , x 2 )ζ * n (x 2 ′ , x 1 ′ ) ω -(E N +2 n -E N 0 ) + iη - ζn (x 1 ′ , x 2 ′ ) ζ * n (x 2 , x 1 ) ω + (E N -2 n -E N 0 ) -iη (2.87)

Background

where the pp amplitudes are defined as

ζ n (x 1 , x 2 ) = ⟨Ψ N 0 | ψ(x 1 ) ψ(x 2 )|Ψ N +2 n ⟩, (2.88) ζn (x 1 ′ , x 2 ′ ) = ⟨Ψ N 0 | ψ † (x 1 ′ ) ψ † (x 2 ′ )|Ψ N -2 n ⟩, (2.89)
and ζ * n , and ζ * n ) are the complex conjugate. The poles of Eq. (2.87) clearly describe the excitation energy of a system where two particles (electrons or holes) are removed. In this case, both non-interacting terms in Eq. (2.77) are time-dependent. Indeed we get

G 1 (1, 1 ′ )G 1 (2, 2 ′ ) = G 1 (x 1 , x 1 ′ ; t 1 -t 1 ′ )G 1 (x 2 , x 2 ′ ; t 1 -t 1 ′ ) (2.90) G 1 (1, 2 ′ )G 1 (2, 1 ′ ) = G 1 (x 1 , x 2 ′ ; t 1 -t 1 ′ )G 1 (x 2 , x 1 ′ ; t 1 -t 1 ′ ). (2.91)
Their spectral representation (where again for simplicity we consider the non-interacting

1-GF) is [G 0 1 (x 1 , x 1 ′ )G 0 1 (x 2 , x 2 ′ )](ω) = -i c,c ′ ψ c (x 1 )ψ * c (x 1 ′ )ψ c ′ (x 2 )ψ * c ′ (x 2 ′ ) ω -(ϵ 0 c + ϵ 0 c ′ ) + iη + i v,v ′ ψ v (x 1 )ψ * v (x 1 ′ )ψ v ′ (x 2 )ψ * v ′ (x 2 ′ ) ω -(ϵ 0 v + ϵ 0 v ′ ) -iη , (2.92)
[G 0 1 (x 1 , x 2 ′ )G 0 1 (x 2 , x 1 ′ )](ω) = -i c,c ′ ψ c (x 1 )ψ * c (x 2 ′ )ψ c ′ (x 2 )ψ * c ′ (x 1 ′ ) ω -(ϵ 0 c + ϵ 0 c ′ ) + iη + i v,v ′ ψ v (x 1 )ψ * v (x 2 ′ )ψ v ′ (x 2 )ψ * v ′ (x 1 ′ ) ω -(ϵ 0 v + ϵ 0 v ′ ) -iη . (2.93)
Here it is clear that both combinations have poles at the energies corresponding to the removal of two particles. Therefore, the non-interacting term for this channel is composed of both combinations.

To solve the 2-GF in practice, a recursive equation with a structure similar to the Dyson equation (2.56) must be defined. These equations are already known in the literature. For the eh channel, the Dyson-like equation is known as the Bethe-Salpeter equation (BSE) [START_REF] Csanak | Green's function technique in atomic and molecular physics[END_REF][START_REF] Salpeter | A relativistic equation for bound-state problems[END_REF], while for the pp channel, it is sometimes referred to as pp BSE [START_REF] Csanak | Green's function technique in atomic and molecular physics[END_REF]. In this work, we will only briefly introduce a common strategy to find them. Instead, we will focus our attention on an intuitive and alternative way to derive them based on the Feynman diagram technique. Therefore, before analyzing these two equations, we will study the diagrammatic structure of the 2-GF.

Feynman diagram analysis

In this section, we will analyze the Feynman diagrams of the interacting part of the 2-GF, namely δG 2 as given in Eq. (2.77). We will focus only on the interacting part, as the non-interacting part has diagrams identical to those of the 1-GF, which has been previously considered.

Let us begin by applying Eq. (2.12) to the 2-GF case

G 2 (1, 2, 1 ′ , 2 ′ ) = - ⟨Ψ N 0 |T [ Ŝ ψI (1) ψI (2) ψ † I (2 ′ ) ψ † I (1 ′ )]|Ψ N 0 ⟩ ⟨Ψ N 0 |T [ Ŝ]|Ψ N 0 ⟩ . (2.94)
The process of obtaining diagrams for the 2-GF follows the same principles as those used for the 1-GF, as discussed in section 2.2.6. The linked-cluster theorem [START_REF] Wen | Quantum Field Theory of Many-Body Systems:From the Origin of Sound to an Origin of Light and Electrons: From the Origin of Sound to an Origin of Light and Electrons[END_REF] ensures that the vacuum diagrams are canceled by the denominator. This holds for every n-GF. Focusing only on the interacting part δG 2 implies that in the Wick contractions every field operator of G 0 2 must be contracted with one from the expansion of the scattering operator Ŝ in such a way that there is at least one interaction that links the two non-interacting 1-GFs. This approach effectively captures the essential interactions between the particles within the 2-GF framework, allowing us to study the correlated dynamics of the system in a diagrammatic manner.

The rules to create the diagrams of the pp 2-GF in Eq. ( 2.86) or of the two-particle correlation function L in Eq. ( 2.85) at order n in the interaction are as follows 1. Draw all topologically distinct connected diagrams with n interaction lines v and 2n + 2 non-interacting 1-GF;

2. Add a prefactor (-1) F (-1) l i n to each diagram, where F is the number of closed fermion loops, and l is the number of interchanges of two incoming or outgoing fermion lines.

To better understand the meaning of l, let us begin by analyzing the particleparticle (pp) channel. In Fig. 2.8, we can see two topologically distinct diagrams that differ by a minus sign due to the interchange of two lines. Figure 2.9 illustrates some examples of second order diagrams in this channel. It is important to note that, for the pp channel, the arrows of the two particles are always in the same direction.

In the analysis of the electron-hole (eh) channel, the two arrows always point in opposite directions, which creates some differences in the possibilities of creating diagrams. Figure 2.10 illustrates the two first-order topologically distinct connected diagrams for the eh channel. Figure 2.11 displays second-order diagrams created by repeating the same first-order contribution. These diagrams will be analyzed in more detail in the next section. Figure 2.12 shows more second-order diagrams that cannot be obtained by repeating first-order contributions.

After this general introduction to 2-GF Feynman diagrams, let us analyze how to construct an equation to calculate the 2-GF. Figure 2.12: Two second-order diagrams that are not created by repeating first-order contributions.

Background

eh channel and BSE

Standard theory

The Bethe-Salpeter equation (BSE) is a Dyson-like equation that describes the eh channel of the 2-GF, or, to be precise, the two-particle correlation function L given by Eq. (2.85). Its main purpose is to incorporate the interaction between the two particles considered. The most important field of application of the BSE is optical absorption. The standard approach to introduce the BSE utilizes a functional derivative approach [START_REF] Salpeter | A relativistic equation for bound-state problems[END_REF][START_REF] Strinati | Application of the green's functions method to the study of the optical properties of semiconductors[END_REF]. One studies how the 1-GF is modified by an external time-dependent potential φ which is set to zero in the end. The two-particle correlation function is thus defined as follows

L(1, 2, 1 ′ , 2 ′ ) ≡ δG(1, 1 ′ ; [φ]) δφ(2 ′ , 2) φ=0 = -G 2 (1, 2, 1 ′ , 2 ′ ) + G 1 (1, 1 ′ )G 1 (2, 2 ′ ). (2.95)
L is hence a generalized response function. It is strictly related to the reducible polarizability or standard response function

χ(1, 2) ≡ δρ(1; [φ]) δφ(2) φ=0 = -iL(1, 2, 1 + , 2 + ), (2.96)
Background that is linked to the irreducible polarizability in Eq. (2.68) through a Dyson-like equation [START_REF] Strinati | Application of the green's functions method to the study of the optical properties of semiconductors[END_REF] b

χ(1, 2) = χ(1, 2) + d3d4 χ(1, 3)v(3, 4)χ(4, 2). (2.97)
From the definition of the two-particle correlation function L (Eq. (2.95)), the BSE can be obtained [START_REF] Csanak | Green's function technique in atomic and molecular physics[END_REF][START_REF] Strinati | Application of the green's functions method to the study of the optical properties of semiconductors[END_REF]. It reads

L(1, 2, 1 ′ , 2 ′ ) = L 0 (1, 2, 1 ′ , 2 ′ ) + L 0 (1, 3, 1 ′ , 3 ′ )Ξ(3 ′ , 4 ′ , 3, 4)L(4, 2, 4 ′ , 2 ′ ) (2.98)
where repeated indices are integrate,

L 0 (1, 2, 1 ′ , 2 ′ ) = G(1, 2 ′ )G(2, 1 ′ ), and
Ξ(1 ′ , 2 ′ , 1, 2) ≡ δΣ(1 ′ , 1) δG(2, 2 ′ ) = -iδ(1, 1 ′ )δ(2, 2 ′ )v(1, 2) + δΣ xc (1 ′ , 1) δG(2, 2 ′ ) , (2.99)
is the kernel of the BSE. Diagrammatically the BSE is represented as

. (2.

100)

To obtain an expression suitable for practical calculations, we need to make approximations to Σ xc . At first order, it corresponds to the exchange term Σ

F xc (1 ′ , 1) = iv(1 ′ , 1)G(1 ′ , 1). Consequently, the kernel Ξ takes the form Ξ HF (1 ′ , 2 ′ , 1, 2) = -iδ(1, 1 ′ )δ(2, 2 ′ )v(1, 2) + iδ(1, 2 ′ )δ(2, 1 ′ )v(1, 2).
(2.101)

The first term on the right-hand side is called eh exchange term and results from the functional derivative of the Hartree potential, while the second is the eh attraction. Figure 2.13 shows its diagrammatic structure.

To include screening in the kernel, we consider the exchange correlation self-energy in the GW approximation. This yields

Ξ GW (1 ′ , 2 ′ , 1, 2) = -iδ(1, 1 ′ )δ(2, 2 ′ )v(1, 2) + iδ(1, 2 ′ )δ(2, 1 ′ )W (1, 2), ( 2

Diagrammatic analysis

Let us now analyze how to obtain the BSE with the GW kernel in Eq. (2.102) from a diagrammatic analysis. First, we need to define the non-interacting terms. Considering the definition of the eh 2-GF in Eq. (2.78), and the analysis of the time dependence in Eqs. (2.82) and ( 2

.83), it is natural to choose L 0 (1, 2, 1 ′ , 2 ′ ) = G 1 (1, 2 ′ )G 1 (2, 1 ′ )
as starting point. This combination of GG terms includes all the poles of the noninteracting term. Notice that this is the combination of GG with a minus sign in the definition of G 0 2 . To add interactions, we aim to construct a Dyson-like equation, namely L = L 0 + L 0 ΞL. To construct the kernel, it is instructive to examine the first-order term in the interaction, namely L 0 ΞL 0 . At this order, the only topologically distinct connected diagrams are represented in Fig. 2.10. Following the Feynman rules for the prefactor, it is natural to define the kernel for the first-order term as the HF kernel (2.101).

We are now going to analyze why, to improve it, we have to screen only the attraction term and not the exchange term. To do that, let us consider as starting (2,[START_REF] Bechstedt | Many-body approach to electronic excitations[END_REF].

Background

L H (1, 2, 1 ′ , 2 ′ ) = G 1 (1, 2 ′ )G 1 (2, 1 ′ ) + (-i)G 1 (1, 3)G 1 (3, 1 ′ )v(3, 4)G 1 (4, 2 ′ )G 1 (2, 4) + (-i) 2 G 1 (1, 3)G 1 (3, 1 ′ )v(3, 4)G 1 (4, 4 ′ )G 1 (4 ′ , 4)v(4 ′ , 5)G 1 (5, 2 ′ )G 1 (2, 5) + ... = G 1 (1, 2 ′ )G 1 (2, 1 ′ ) + G 1 (1, 3)G 1 (3, 1 ′ )W RP A (3, 4)G 1 (4, 2 ′ )G 1
(2.103)

Its diagrammatic representation is shown in Fig. 2.14. This equation clearly shows that the BSE naturally screened the eh exchange interaction. Therefore, if the Hartree potential is screened in the kernel, e.g., Ξ

H (1 ′ , 2 ′ , 1, 2) = -iδ(1, 1 ′ )δ(2, 2 ′ )W RP A (1, 2),
we would double count screening contributions.

By doing the same analysis with only the attraction term, the result is different. Figure (2.15) shows that the BSE is not able to improve the eh. The only way to improve this term is to improve the bare Coulomb interaction directly in BSE kernel. The more natural way is to screen it. Notice that this screening will not lead to doublecounting problems. This analysis shows an intuitive way to obtain the standard BSE kernel in the GW approximation only by looking at the diagrams.

It is important to notice that, considering the full kernel Ξ = -iv + iW , starting from third-order in v, diagrams like the one in Fig. 2.16 will appear. This shows that for the eh exchange interaction screening beyond the RPA level is naturally included. Another important conclusion is that for improving the screening of the eh exchange interaction one has to improve the part of the kernel created by δΣ xc /δG.

How to solve the BSE in practice

The two-particle correlation function is linked to experiments through the response function (Eq. (2.96)). To have a simpler relation between the two quantities L = -iL is defined. The BSE for L is

L(1, 2, 1 ′ , 2 ′ ) = L0 (1, 2, 1 ′ , 2 ′ ) + L0 (1, 3, 1 ′ , 3 ′ ) Ξ(3 ′ , 4 ′ , 3, 4) L(4, 2, 4 ′ , 2 ′ ), (2.104)
with the GW kernel equal to

Ξ(1 ′ , 2 ′ , 1, 2) = δ(1, 1 ′ )δ(2, 2 ′ )v(1, 2) -δ(1, 2 ′ )δ(2, 1 ′ )W (1, 2), (2.105)
and the RPAx kernel equal to

Ξ(1 ′ , 2 ′ , 1, 2) = δ(1, 1 ′ )δ(2, 2 ′ )v(1, 2) -δ(1, 2 ′ )δ(2, 1 ′ )v(1, 2). (2.106)
The BSE for L in Eq. (2.104) is the one commonly solved.

With the following change of basis

Lij;mo (ω) = dx 1 dx 2 dx 1 ′ dx 2 ′ L(x 1 , x 2 , x 1 ′ , x 2 ′ ; ω)ϕ * i (x 1 )ϕ j (x 1 ′ )ϕ m (x 2 ′ )ϕ * o (x 2 ),
(2.107) the non-interacting term L0 is given by

L0 ij;mo = δ im δ jo (f i -f j ) ϵ i -ϵ j -ω + iηsign(f i -f j ) , (2.108)
where f i is the spinorbital occupation (f i = 1 if ϵ i ≤ µ and 0 otherwise, with µ the chemical potential) and we used the quasiparticle approximation for the 1-GF [START_REF] Martin | Interacting Electrons: Theory and Computational Approaches[END_REF][START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF].

The BSE becomes Lij;mo (ω) = L0 ij;mo (ω) + i ′ j ′ m ′ o ′ L0 ij;i ′ j ′ (ω) Ξi ′ j ′ ;m ′ o ′ Lm ′ o ′ ;mo (ω), (2.109) 36

Background

where the kernel Ξ is considered static. The BSE expressed in Eq. (2.109) can be formally solved [START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF] giving

Lij;mo = H 2p -Iω -1 ij;mo (f o -f m ), (2.110)
where we defined a two-particle Hamiltonian

H 2p ij;mo = (ϵ i -ϵ j )δ im δ jo + (f j -f i ) Ξij;mo . (2.111)
The presence of (f o -f m ) in Eq. (2.110) selects only a part of the Hamiltonian H 2p , which we refer to as excitonic Hamiltonian H 2p,exc [START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF],

H 2p,exc ij;mo = H 2p,reso vc;v ′ c ′ Ξcoupling vc;c ′ v ′ -[ Ξcoupling cv ;v ′ c ′ ] * -[H 2p,reso cv ;c ′ v ′ ] * . (2.112)
Here we have indicated the occupied states with v, v ′ , and the unoccupied ones with c, c ′ . The matrix H 2p,reso only involves positive energy transitions and it is given by,

H 2p,reso vc;v ′ c ′ = (ϵ c -ϵ v )δ v,v ′ δ c,c ′ + Ξvc;v ′ c ′ , (2.113)
with Ξvc;v ′ c ′ = v vcc ′ v ′ -W vcv ′ c ′ .
(2.114)

The matrix elements of the so-called electron-hole exchange v vcc ′ v ′ and of the screened electron-hole interaction W vcv ′ c ′ are given by,

v vcc ′ v ′ = dxdx ′ ϕ * v (x)ϕ * c (x ′ ) 1 |r -r ′ | ϕ c ′ (x ′ )ϕ v ′ (x), (2.115)
W vcv ′ c ′ = dxdx ′ ϕ * v (x)ϕ * c (x ′ )W (x, x ′ )ϕ v ′ (x ′ )ϕ c ′ (x).
(2.116)

The lower-right block -[H 2p,reso vcv ′ c ′ ] * is the anti-resonant part of the excitonic Hamiltonian which only involves negative energy solutions, while the off-diagonal terms couple the eigenvalue equations of the resonant and anti-resonant Hamiltonians. By using the spectral representation of the inverse excitonic Hamiltonian,

H 2p,exc -I ω -1 ij;mo = λ,λ ′ A ij λ S -1 λλ ′ A * mo λ ′ E λ - ω , (2.117)
where A ij λ and E λ are the eigenvectors and eigenvalues, respectively, of the excitonic Hamiltonian

H 2p,exc ij;mo A mo λ = E λ A ij λ , (2.118)
Background 37 and S λλ ′ is the overlap matrix

S λλ′ = ij A * ij λ A ij λ ′ , (2.119)
which differs from the identity due to the non-orthogonality of the eigenvectors of a non-Hermitian matrix, we arrive at

Lij;mo = λ,λ ′ A ij λ S -1 λλ ′ A * mo λ ′ E λ -ω (f o -f j ).
(2.120)

pp channel

The pp channel of the 2-GF describes the propagation of two particles (electrons or holes).We introduce it because it will be used in chapter 4.

The standard procedure to obtain the pp BSE is through the functional derivative of the anomalous Green's function [START_REF] Somà | Ab initio self-consistent gorkov-green's function calculations of semi-magic nuclei: Numerical implementation at second order with a two-nucleon interaction[END_REF]. We propose two alternative approaches to derive it: 1) from a diagrammatic analysis; 2) from the Martin-Schwinger hierarchy.

Diagrammatic analysis

As for the eh BSE, we start the analysis from the non-interacting term. Both combinations of GG depend on a time difference (see Eqs. (2.90) and (2.91)). Therefore, both terms have to be taken into account in G

0(pp) 2 G 0(pp) 2 (1, 2, 1 ′ , 2 ′ ) = G 1 (1, 1 ′ )G 1 (2, 2 ′ ) -G 1 (1, 2 ′ )G 1 (2, 1 ′ ).
(2.121)

Following the structure of a Dyson-like equation we define the pp BSE as

G pp 2 (1, 2, 1 ′ , 2 ′ ) = G 0(pp) 2 (1, 2, 1 ′ , 2 ′ ) + G 0(pp) 2 (1, 2, 3 ′ , 4 ′ )K(3 ′ , 4 ′ , 3, 4)G pp 2 (3, 4, 1 ′ , 2 ′ ), ( 2
.122) To construct the kernel K we iterate the pp BSE and compare the first order in the interaction G 0 2 KG 0 2 with the first order diagrams, which are given in Fig. 2.8. This comparison suggests the kernel

K(1 ′ , 2 ′ , 1, 2) = i 2 δ(1, 1 ′ )δ(2, 2 ′ )v(1, 2)
, where the prefactor 1/2 prevents the double counting of diagrams. Notice that, the two diagrams in Fig. 2.17 are topologically equivalent and therefore we can equivalently rewrite

G pp 2 (1, 2, 1 ′ , 2 ′ ) = G 0(pp) 2 (1, 2, 1 ′ , 2 ′ ) + G 1 (1, 3 ′ )G 1 (2, 4 ′ ) K(3 ′ , 4 ′ , 3, 4)G pp 2 (3, 4, 1 ′ , 2 ′ ), (2.123) with K(1 ′ , 2 ′ , 1, 2) = iδ(1, 1 ′ )δ(2, 2 ′ )v(1, 2). (2.124) Background Figure 2
.17: Two first order diagrams in the pp channel. They are topologically equivalent since in both cases the arrow move from 1 (2) to 1 ′ (2 ′ ).

Analyzing the next orders in the interaction, it becomes clear that this kernel behaves like the attractive part of the eh kernel. Therefore, the pp BSE behaves similarly as it is shown in Fig. 2.15, namely it is not able to improve the bare Coulomb interaction between the particle pairs. Adding screening to the bare Coulomb potential is the most natural way to improve the kernel (2.124). Notice that the screening does not create double-counting problems.

pp BSE from the Martin-Schwinger hierarchy

An alternative way to obtain a Dyson-like equation for G pp 2 is from the Martin-Schwinger (MS) hierarchy (Eq. (2.15)). The relations we will use are general, therefore in the following we drop the superscript "pp". Let us start by constructing two equations from the MS hierarchy. In the first one we apply the time derivative with respect to t 1 , and in the second one we apply the time derivative with respect to t 2

G 2 (1, 2, 1 ′ , 2 ′ ) = G 0 1 (1, 1 ′ )G 1 (2, 2 ′ ) -G 0 1 (1, 2 ′ )G 1 (2, 1 ′ ) + iG 0 1 (1, 3 ′ )v(3 ′ , 4 ′ , 3, 4)G 3 (2, 3, 4, 1 ′ , 2 ′ , 4 ′ ), G 2 (1, 2, 1 ′ , 2 ′ ) = G 0 1 (2, 2 ′ )G 1 (1, 1 ′ ) -G 0 1 (2, 1 ′ )G 1 (1, 2 ′ ) -iG 0 1 (2, 3 ′ )v(3 ′ , 4 ′ , 3, 4)G 3 (1, 3, 4, 1 ′ , 2 ′ , 4 ′ ), (2.125) where v(1, 2, 1 ′ , 2 ′ ) = δ(1, 1 ′ )δ(2, 2 ′ )v(1, 2)
. By applying Eq. (A.3) from the right of Eq. (2.125), we get

G 2 (1, 2, 1 ′ , 2 ′ ) = G 1 (1, 1 ′ )G 1 (2, 2 ′ ) -G 1 (1, 2 ′ )G 1 (2, 1 ′ ) + iG 1 (1, 3 ′ )v(3 ′ , 4 ′ , 3, 4) × [G 3 (2, 3, 4, 1 ′ , 2 ′ , 4 ′ ) -G 2 (3, 4, 5 ′ , 4 ′ )G -1 1 (5 ′ , 5)G 2 (5, 2, 1 ′ , 2 ′ )], G 2 (1, 2, 1 ′ , 2 ′ ) = G 1 (1, 1 ′ )G 1 (2, 2 ′ ) -G 1 (1, 2 ′ )G 1 (2, 1 ′ ) -iG 1 (2, 3 ′ )v(3 ′ , 4 ′ , 3, 4) × [G 3 (1, 3, 4, 1 ′ , 2 ′ , 4 ′ ) -G 2 (3, 4, 5 ′ , 4 ′ )G -1 1 (5 ′ , 5)G 2 (5, 1, 1 ′ , 2 ′ )].
(2.126)

We now introduce GG -1 in the left-hand side of the square brackets, and we apply Eq.(A.9) from the right in order to have a recursive equation. We obtain

G 2 (1, 2, 1 ′ , 2 ′ ) = G 0 2 (1, 2, 1 ′ , 2 ′ ) + iG 1 (1, 3 ′ )G 1 (2, 6 ′ )G -1 1 (6 ′ , 6)v(3 ′ , 4 ′ , 3, 4) ×[G 3 (6,3,4,7 ′ ,8 ′ ,4 ′ )-G 2 (3,4,5 ′ ,4 ′ )G -1 1 (5 ′ ,5)G 2 (5,6,7 ′ , 8 ′ )] 1 2 G -1 2 (7 ′ , 8 ′ , 7,8)G 2 (7,8,1 ′ , 2 ′ ) G 2 (1, 2, 1 ′ , 2 ′ ) = G 0 2 (1, 2, 1 ′ , 2 ′ ) -iG 1 (2, 3 ′ )G 1 (1, 6 ′ )G -1 1 (6 ′ , 6)v(3 ′ , 4 ′ , 3, 4) ×[G 3 (6,3,4,7 ′ ,8 ′ ,4 ′ )-G 2 (3,4,5 ′ ,4 ′ )G -1 1 (5 ′ ,5)G 2 (5,6,7 ′ ,8 ′ )] 1 2 G -1 2 (7 ′ , 8 ′ , 7,8)G 2 (7,8,1 ′ ,2 ′ ). (2.127)
By summing the two equations, we arrive at

G 2 (1, 2, 1 ′ , 2 ′ ) = G 0 2 (1, 2, 1 ′ , 2 ′ ) + i 2 [G 0 2 (1, 2, 3 ′ , 6 ′ )]G -1 1 (6 ′ , 6)v(3 ′ , 4 ′ , 3, 4) ×[G 3 (6,3,4,7 ′ ,8 ′ ,4 ′ )-G 2 (3,4,5 ′ ,4 ′ )G -1 1 (5 ′ ,5)G 2 (5,6,7 ′ ,8 ′ )] 1 2 G -1 2 (7 ′ ,8 ′ ,7,8)G 2 (7,8,1 ′ ,2 ′ ) = G 0 2 (1, 2, 1 ′ , 2 ′ ) + G 0 2 (1, 2, 3 ′ , 6 ′ )K pp (3 ′ , 6 ′ , 7, 8)G 2 (7, 8, 1 ′ , 2 ′ ), (2.128)
where we defined the kernel as

K pp (3 ′ , 6 ′ , 7, 8) = i 4 G -1 1 (6 ′ , 6)v(3 ′ , 4 ′ , 3, 4)[G 3 (6, 3, 4, 7 ′ , 8 ′ , 4 ′ ) -G 2 (3, 4, 5 ′ , 4 ′ )G -1 1 (5 ′ , 5)G 2 (5, 6, 7 ′ , 8 ′ )]G -1 2 (7 ′ , 8 ′ , 7, 8) = i 4 G -1 1 (6 ′ , 6)v(3 ′ , 4 ′ , 3, 4)R(6, 3, 4, 7 ′ , 8 ′ 4 ′ )G -1 2 (7 ′ , 8 ′ , 7, 8) (2.129) with R = G 3 -G 2 G -1
1 G 2 the so-called pph response function [START_REF] Ethofer | The six-point Green function and the application of Faddeev's theory to the many body problem[END_REF][START_REF] Schuck | Equation of motion method for strongly correlated fermi systems and extended rpa approaches[END_REF]. We associate this equation with the pp channel due to the structure of its non-interacting term. To compare this kernel to the one defined by the diagrammatic analysis (Eq. (2.124)) we use the HF approximation for R and G -1 2 . The HF kernel then reads

K HF (1 ′ , 2 ′ , 1, 2) = i 4 G -1 1 (2 ′ , 5)v(1 ′ , 3 ′ , 3, 4)G 1 (5, 3 ′ ) [G 1 (4, 5 ′ )G 1 (3, 4 ′ ) -G 1 (4, 4 ′ )G 1 (3, 5 ′ )][G 0 2 ] -1 (4 ′ , 5 ′ , 1, 2) = i 4 v(1 ′ , 2 ′ , 3, 4)G 0 2 (3, 4, 4 ′ , 5 ′ )[G 0 2 ] -1 (4 ′ , 5 ′ , 1, 2) = i 4 v(1 ′ , 2 ′ , 3, 4)[δ(1, 3)δ(2, 4) -δ(1, 4)δ(2, 3)] = i 4 [v(1 ′ , 2 ′ , 1, 2) -v(1 ′ , 2 ′ , 2, 1)] = i 4 [δ(1, 1 ′ )δ(2, 2 ′ )v(1, 2) -δ(2, 1 ′ )δ(1, 2 ′ )v(1, 2)], (2.130)
40

Background

where we used Eq. (A.5) for [G 0 2 ] -1 . By exploiting symmetry relations, it is easy to show that Eq. (2.128) with the HF kernel in Eq. (2.130) is equal to Eq. (2.123) with first order kernel in Eq. (2.124).

We have shown a method to analytically obtain the pp kernel (Eq. (2.129)). The analysis of this kernel is very complicated and outside the scope of this thesis. However, we showed that in an independent particle picture it coincides with the kernel obtained from the diagrammatic analysis Eq. (2.124).

As for the eh channel, also in the pp case it is instructive to project the equation on the space that diagonalize G 0 1 , Eq. (2.47). By using the change of basis of Eq. (2.22), the matrix elements of the pp 2-GF are

G pp ij,mo (ω) = dx 1 dx 1 ′ dx 2 dx 2 ′ ϕ * i (x 1 )ϕ * j (x 2 )G pp (x 1 , x 2 , x 1 ′ , x 2 ′ ; ω)ϕ m (x 1 ′ )ϕ o (x 2 ′ ). ( 2
G 0(pp) ij,mo (ω) = [G 0 im G 0 jo ](ω) -[G 0 io G 0 jm ](ω) = -i (δ im δ jo -δ io δ jm )(1 -f i -f j ) ω -(ϵ 0 i + ϵ 0 j ) + iηsign(ϵ 0 i -µ)
.

(2.132)

Making use of the anti-commutation relations of the creation and annihilation operator the following symmetry relations hold

G pp ij,mo = -G pp ji,mo = -G pp ij,om = G pp ji,om .
(2.133)

The above equation shows redundant information, to avoid them the restrictions i > j, and m > o are used. With these restrictions, the condition δ io δ jm in the non-interacting term (Eq. (2.132)) is never fulfilled. Therefore, only the combination G 0 im G 0 jo is non-zero, and the pp BSE in the HF approximation reduces to

G pp ij,mo (ω) = [G 0 im G 0 jo ](ω) + [G 0 ii ′ G 0 jj ′ ](ω)K HF i ′ j ′ ,m ′ o ′ G pp m ′ o ′ ,mo (ω) (2.134) with K HF ij,mo = iv ijom -iv ijmo ≡ iv ijom (2.135) and v ijom = dx 1 dx 2 ϕ * i (x 1 )ϕ * j (x 2 )v(r 1 , r 2 )ϕ o (x 2 )ϕ m (x 1 ). (2.136)
We finally express the HF kernel (2.135) in position space

K HF (1, 2, 1 ′ , 2 ′ ) = i[δ(1, 1 ′ )δ(2, 2 ′ )v(1, 2) -δ(2, 1 ′ )δ(1, 2 ′ )v(1, 2)],
(2.137) the above kernel is shown diagrammatically in Fig. 2.18, and it is the particle-particle kernel we will use in section 4. Notice that, like for the eh BSE, it is possible to recast the equation collecting a -i, i.e., by using the following definitions: Gpp = -iG pp and K = -iK.

CHAPTER 3

Beyond G 1 and G 2 : the 3-body and 4-body GF

To study direct and inverse photoemission spectroscopy the one-body Green's function is the most natural quantity since its spectral function is directly linked to electron removal and addition. However, the approximations needed to calculate it, yield spectra that are not always accurate. In particular, the description of the satellite structure is very difficult. In this chapter, we show that it is possible to calculate the spectra also starting from the three-body Green's function. We note that this is a general strategy: the more information the fundamental quantity contains the less information is required in the effective potential, i.e. the self-energy in our case, to describe the relevant many-body effects. Indeed, we will show that already at the level of the non-interacting 3-GF there is information related to satellites. Therefore, as we will show, a static self-energy (3-body self-energy) is sufficient to obtain both quasiparticles and satellites in the photoemission spectra.

The three-body Green's function

According to Eq. (2.6) the 3-GF is defined by

G 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = i⟨Ψ N 0 |T [ ψ(1) ψ(2) ψ(3) ψ † (3 ′ ) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N 0 ⟩. (3.1) 44
The 3-GF and the 4-GF

To understand its composition we use Wick's theorem in Eq. (2.12). We can rewrite the 3-GF as

G 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = G 0 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) + G 1 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) + δG 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) (3.2)
where the superscript 0, 1 refers to the number of interactions between particles.

The first term on the right-hand side of Eq. (3.2) is composed of all the possible combinations of G 1 G 1 G 1 without interaction between them. We refer to it as the non-interacting term. It is equal to

G 0 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = G 1 (1, 1 ′ )G 1 (2, 2 ′ )G 1 (3, 3 ′ ) + G 1 (1, 2 ′ )G 1 (2, 3 ′ )G 1 (3, 1 ′ ) + G 1 (1, 3 ′ )G 1 (2, 1 ′ )G 1 (3, 2 ′ ) -G 1 (1, 1 ′ )G 1 (2, 3 ′ )G 1 (3, 2 ′ ) -G 1 (1, 2 ′ )G 1 (2, 1 ′ )G 1 (3, 3 ′ ) -G 1 (1, 3 ′ )G 1 (2, 2 ′ )G 1 (3, 1 ′ ). (3.3)
The second term, G 1 3 , is composed of all the possible combinations of G 1 δG 2 , namely a free particle and δG 2 , which corresponds to an interacting two-particle system, that propagate without interaction between them. It is equal to

G 1 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = G 1 (1, 1 ′ )δG 2 (2, 3, 2 ′ , 3 ′ ) -G 1 (1, 2 ′ )δG 2 (2, 3, 1 ′ , 3 ′ ) + G 1 (1, 3 ′ )δG 2 (2, 3, 1 ′ , 2 ′ ) -G 1 (2, 1 ′ )δG 2 (1, 3, 2 ′ , 3 ′ ) + G 1 (2, 2 ′ )δG 2 (1, 3, 1 ′ , 3 ′ ) -G 1 (2, 3 ′ )δG 2 (1, 3, 1 ′ , 2 ′ ) + G 1 (3, 1 ′ )δG 2 (1, 2, 2 ′ , 3 ′ ) -G 1 (3, 2 ′ )δG 2 (1, 2, 1 ′ , 3 ′ ) + G 1 (3, 3 ′ )δG 2 (1, 2, 1 ′ , 2 ′ ) (3.4)
The last term, δG 3 , represents the propagation of three interacting particles. Since the Coulomb interaction in the Hamiltonian (2.1) is a two-body term, δG 3 starts to be different from zero at second order in the interaction. An equation for this term is shown in App. B.

The 3-GF depends on six times or five time differences when the Hamiltonian is time independent, and the total number of permutations of the field operators in Eq. (3.1) due to the T operator is 6! = 720. Depending on the order of the field operators (and therefore of the times) the 3-GF yields different information. In general, it describes the propagation of three particles (electrons or holes) and the 3-GF can therefore be split into four components:

G eee 3 , G hhh 3 , G eeh 3
and G hhe 3 . In order to make this separation explicit one can rewrite the six time-ordered field operators in Eq. (3.1) as a sum of products of two terms each containing three time-ordered field operators

G 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = = i n ⟨Ψ N 0 |T [ ψ(1) ψ(2) ψ(3)]|Ψ N +3 n ⟩⟨Ψ N +3 n |T [ ψ † (3 ′ ) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N 0 ⟩ -i n ⟨Ψ N 0 |T [ ψ † (3 ′ ) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N -3 n ⟩⟨Ψ N -3 n |T [ ψ(1) ψ(2) ψ(3)]|Ψ N 0 ⟩ -i n ⟨Ψ N 0 |T [ ψ(1) ψ(2) ψ † (3 ′ )]|Ψ N +1 n ⟩⟨Ψ N +1 n |T [ ψ(3) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N 0 ⟩ + i n ⟨Ψ N 0 |T [ ψ(3) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N -1 n ⟩⟨Ψ N -1 n |T [ ψ(1) ψ(2) ψ † (3 ′ )]|Ψ N 0 ⟩ + other orderings. ( 3

The eeh and hhe 3-GF

To study photoemission spectra we are not interested in the eee and hhh parts. Therefore, we focus our attention on G hhe 3 and G h 3 terms considered in Eq. (3.5) to show that they contain information about removal and addition energies. We follow a similar procedure that Csanak et al. [START_REF] Csanak | Green's function technique in atomic and molecular physics[END_REF] used to find the eh/he part of the two-particle Green's function.

Case 1: t 1 , t 2 , t 3 ′ > t 3 , t 1 ′ , t 2 ′ We set t 1 , t 2 , t 3 ′ > t 3 , t 1 ′ , t 2 ′
without fixing the order of t 1 , t 2 and t 3 ′ , and of t 3 , t 1 ′ and t 2 ′ . Then, for this time ordering, we have

G e 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = = -i n ⟨Ψ N 0 |T [ ψ(1) ψ(2) ψ † (3 ′ )]|Ψ N +1 n ⟩⟨Ψ N +1 n |T [ ψ(3) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N 0 ⟩ = i n χ n (1, 2, 3 ′ ) χn (1 ′ , 2 ′ , 3). (3.6)
The 3-GF and the 4-GF

The electron-electron-hole amplitudes have been defined as

χ n (1, 2, 3 ′ ) = ⟨Ψ N 0 |T [ ψ(1) ψ(2) ψ † (3 ′ )]|Ψ N +1 n ⟩ χn (1 ′ , 2 ′ , 3) = ⟨Ψ N +1 n |T [ ψ † (1 ′ ) ψ † (2 ′ ) ψ(3)]|Ψ N 0 ⟩. (3.7)
Making explicit the times in the Heisenberg representation of the field operators, these amplitudes can be rewritten as

χ n (1, 2, 3 ′ ) = exp[i/3(t 1 + t 2 + t 3 ′ )(E N 0 -E N +1 n )]X n (x 1 , x 2 , x 3 ′ ; τ 12 , τ 23 ′ ) χn (1 ′ , 2 ′ , 3) = exp[-i/3(t 1 ′ + t 2 ′ + t 3 )(E N 0 -E N +1 n )] Xn (x 1 ′ , x 2 ′ , x 3 ; τ 1 ′ 2 ′ , τ 2 ′ 3 ) (3.8)
where X n and Xn are defined by

X n (x 1 , x 2 , x 3 ′ ; τ 12 , τ 23 ′ ) = = i̸ =j̸ =k=1,2,3 ′ (-1) P θ(τ ij )θ(τ jk ) exp[ i 3 (E N 0 (2τ ij +τ jk )+E N +1 n (2τ jk +τ ij ))] × ⟨Ψ N 0 |Υ(x i )e -iHτij Υ(x j )e -iHτ jk Υ(x k )|Ψ N +1 n ⟩ (3.9) Xn (x 1 ′ , x 2 ′ , x 3 ; τ 1 ′ 2 ′ , τ 2 ′ 3 ) = = i̸ =j̸ =k=1 ′ ,2 ′ ,3 (-1) P θ(τ ij )θ(τ jk ) exp[ i 3 (E N +1 n (2τ ij + τ jk ) + E N 0 (2τ jk + τ ij ))] × ⟨Ψ N +1 n |Υ(x i )e -iHτij Υ(x j )e -iHτ jk Υ(x k )|Ψ N 0 ⟩ (3. 10
)
where P is the number of permutations with respect to the initial ordering i = 1,

j = 2, k = 3 ′ or i = 1 ′ , j = 2 ′ , k = 3. Finally, Υ(x i ) is given by Υ(x i ) = ψ(x i ) if i = 1, 2, 3 ψ † (x i ) if i = 1 ′ , 2 ′ , 3 ′ . ( 3
G e 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = i n exp[iτ (E N 0 -E N +1 n )] × X n (x 1 , x 2 , x 3 ′ ; τ 12 , τ 23 ′ ) Xn (x 1 ′ , x 2 ′ , x 3 ; τ 1 ′ 2 ′ , τ 2 ′ 3 ) (3.12)
where we defined

τ = 1 3 (t 1 + t 2 + t 3 ′ ) - 1 3 (t 3 + t 1 ′ + t 2 ′ ) τ ij = t i -t j . (3.13)
The important point is that Eq. (3.12) depends on τ only through an exponential factor in which it multiplies the electron addition energies. In this way, Eq. (3.12)

Theory 47 has a similar structure as the electron part of the 1-GF.

Case 2:

t 3 , t 1 ′ , t 2 ′ > t 1 , t 2 , t 3 ′
For this ordering of times we obtain

G h 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = = i n ⟨Ψ N 0 |T [ ψ(3) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N -1 n ⟩⟨Ψ N -1 n |T [ ψ(1) ψ(2) ψ † (3 ′ )]|Ψ N 0 ⟩ = -i n ζn (1 ′ , 2 ′ , 3)ζ n (1, 2, 3 ′ ). (3.14)
The hole-hole-electron amplitudes have been defined as

ζ n (1, 2, 3 ′ ) = ⟨Ψ N -1 n |T [ ψ(1) ψ(2) ψ † (3 ′ )]|Ψ N 0 ⟩ ζn (1 ′ , 2 ′ , 3 ′ ) = ⟨Ψ N 0 |T [ ψ † (1 ′ ) ψ † (2 ′ ) ψ(3)]|Ψ N -1 n ⟩. (3.15)
As before, we make explicit the times in the Heisenberg representation of the field operators, and we arrive at

ζ n (1, 2, 3 ′ ) = exp[-i/3(t 1 + t 2 + t 3 ′ )(E N 0 -E N -1 n )]Z n (x 1 , x 2 , x 3 ′ ; τ 12 , τ 23 ′ ) ζn (1 ′ , 2 ′ , 3) = exp[i/3(t 1 ′ + t 2 ′ + t 3 )(E N 0 -E N -1 n )] Zn (x 1 ′ , x 2 ′ , x 3 ; τ 1 ′ 2 ′ , τ 2 ′ 3 ) (3.16)
where Z and Z are defined by

Z n (x 1 , x 2 , x 3 ′ ; τ 12 , τ 23 ′ ) = = i̸ =j̸ =k=1,2,3 ′ (-1) P θ(τ ij )θ(τ jk ) exp[ i 3 (E N 0 (2τ jk + τ ij ) + E N -1 n (2τ ij + τ jk ))] × ⟨Ψ N -1 n |Υ(x i )e -iHτij Υ(x j )e -iHτ jk Υ(x k )|Ψ N 0 ⟩ (3.17) Zn (x 1 ′ , x 2 ′ , x 3 ; τ 1 ′ 2 ′ , τ 2 ′ 3 ) = = i̸ =j̸ =k=1 ′ ,2 ′ ,3 (-1) P θ(τ ij )θ(τ jk ) exp[ i 3 (E N -1 n (2τ jk + τ ij ) + E N 0 (2τ ij + τ jk ))] × ⟨Ψ N 0 |Υ(x i )e -iHτij Υ(x j )e -iHτ jk Υ(x k )|Ψ N -1 n ⟩. ( 3
G h 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = -i n exp[-iτ (E N 0 -E N -1 n )] × Zn (x 1 ′ , x 2 ′ , x 3 ; τ 1 ′ 2 ′ , τ 2 ′ 3 )Z n (x 1 , x 2 , x 3 ′ ; τ 12 , τ 23 ′ ). ( 3

.19)

The 3-GF and the 4-GF

In this case the time difference τ in the exponential function multiplies electron removal energies. By a similar analysis, one can show that the other time orderings do not have factorizable exponential functions in terms of τ . Therefore, we can write G 3 as follows

G 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = G e 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ )θ(τ + F (τ 12 , τ 3 ′ 1 , τ 1 ′ 2 ′ , τ 31 ′ )) +G h 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ )θ(-τ + F (τ 1 ′ 2 ′ , τ 31 ′ , τ 12 , τ 3 ′ 1 )) + other orderings, (3.20)
where the time orderings of the two cases described above are ensured by the Heaviside functions, and F is defined as

F (τ 12 , τ 3 ′ 1 , τ 1 ′ 2 ′ , τ 31 ′ ) = i̸ =j̸ =k=1,2,3 ′ 1 3 (τ ij -τ ki )θ(τ jk )θ(τ ki ) - i̸ =j̸ =k=1 ′ ,2 ′ ,3 1 3 (τ ij -τ ki )θ(τ jk )θ(τ ij ). (3.21)
In any of the terms in other orderings, it is impossible to factorize an exponential function of the form exp[±iτ (

E N 0 -E N ±1 n )]
, with τ defined in Eq. (3.13). To obtain the spectral representation of the 3-GF, we Fourier transform with respect to τ , which yields

G e+h 3 (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ , x 3 ′ ;τ 12 , τ 23 ′ , τ 1 ′ 2 ′ , τ 2 ′ 3 , ω) = = - n e -i[ω-(E N +1 n -E N 0 )]F (τ12,τ 3 ′ 1 ,τ 1 ′ 2 ′ ,τ 31 ′ ) × X n (x 1 , x 2 , x 3 ′ ; τ 12 , τ 23 ′ ) Xn (x 1 ′ , x 2 ′ , x 3 ; τ 1 ′ 2 ′ , τ 2 ′ 3 ) ω -(E N +1 n -E N 0 ) + iη - n e -i[ω-(E N 0 -E N -1 n )]F (τ 1 ′ 2 ′ ,τ 31 ′ ,τ12,τ 3 ′ 1 ) × Zn (x 1 ′ , x 2 ′ , x 3 ; τ 1 ′ 2 ′ , τ 2 ′ 3 )Z n (x 1 , x 2 , x 3 ′ ; τ 12 , τ 23 ′ ) ω -(E N 0 -E N -1 n ) -iη . (3.22)
All the terms in other orderings are nonsingular at frequencies equal to electron removal or addition energies. Therefore, their contribution to the final spectra is null. From Eq. (3.22) it is clear that the poles of the eeh and the hhe parts of the 3-GF are the same as the poles of the 1-GF, see Eq. (2.53). In analogy to the time difference of the 1-GF, it is possible to interpret τ as the time of the combined propagation of the added particle (electron or hole) and the electron-hole pair. The four remaining time differences correspond to the following physical processes: 1) the time between the added particle and the creation of the electron-hole pair; 2) the Theory time needed to create the electron-hole pair; 3) the time needed to recombine the electron-hole pair; 4) the time between the recombination and the removal of the particle. Which time difference corresponds to which process depends on the order of the times. For the description of (inverse) photoemission spectroscopy all four processes can be considered instantaneous. Therefore, we can take the limit τ ij → 0 for each of the four time differences. It is convenient to choose a time ordering that is coherent with the chronology of the (inverse) photoemission process. For example, in direct photoemission spectroscopy first an electron is emitted from the system leading to the creation of electron-hole pairs. After a time τ the electron-hole pairs recombine and finally an electron is added. This corresponds to the following order of the field operators ψ † ψ † ψ ψ † ψ ψ that act on |Ψ N 0 ⟩. This order of the field operators can be obtained with the following choice for the time differences,

τ 12 = 0 -, τ 23 ′ = 0 -, τ 1 ′ 2 ′ = 0 + , τ 2 ′ 3 = 0 + . (3.23)
We notice that, due to the presence of the Heaviside step functions in Eqs. (3.9), (3.17), (3.10), and (3.18), only one term in the sum remains after fixing the time differences. We note that other choices for the time differences are possible in order to obtain the same order of creation and annihilation operators mentioned above. The final result does not depend on this choice. With the time differences given in Eq. (3.23) we obtain the following expression for G e+h

3 G e+h 3 (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ , x 3 ′ ; ω) = = n X n (x 1 , x 2 , x 3 ′ )X * n (x 1 ′ , x 2 ′ , x 3 ) ω -(E N +1 n -E N 0 ) + iη + n Z * n (x 1 ′ , x 2 ′ , x 3 )Z n (x 1 , x 2 , x 3 ′ ) ω -(E N 0 -E N -1 n ) -iη , (3.24)
where the electron-electron-hole and hole-hole-electron amplitudes, X n and Z n , respectively, are defined as

X n (x 1 , x 2 , x 3 ′ ) = ⟨Ψ N 0 | ψ † (x 3 ′ ) ψ(x 2 ) ψ(x 1 )|Ψ N +1 n ⟩, (3.25)
Z n (x 1 , x 2 , x 3 ′ ) = ⟨Ψ N -1 n | ψ † (x 3 ′ ) ψ(x 2 ) ψ(x 1 )|Ψ N 0 ⟩. (3.26)
For completeness, we also give here the explicit expressions of the complex conjugates of these amplitudes,

X * n (x 1 ′ , x 2 ′ , x 3 ) = ⟨Ψ N +1 n | ψ † (x 1 ′ ) ψ † (x 2 ′ ) ψ(x 3 )|Ψ N 0 ⟩, (3.27)
Z * n (x 1 ′ , x 2 ′ , x 3 ) = ⟨Ψ N 0 | ψ † (x 1 ′ ) ψ † (x 2 ′ ) ψ(x 3 )|Ψ N -1 n ⟩. (3.28)
The representation of G 3 in Eq. (3.24) is similar to the Lehmann representation of G 1 , i.e., the poles are the same but the amplitude corresponding to each pole is different.

The 3-GF and the 4-GF

We can thus define the spectral function

A 3 (ω) corresponding to G e+h 3 (ω) according to A 3 (ω) = 1 π sign(µ -ω)ImG e+h 3 (ω), (3.29)
where, for notational convenience, the spin-position arguments are omitted. It can be verified that G e+h

3

(ω) can be retrieved from A 3 (ω) according to

G e+h 3 (ω) = µ -∞ dω ′ A 3 (ω ′ ) ω -ω ′ -iη + +∞ µ dω ′ A 3 (ω ′ ) ω -ω ′ + iη . (3.30)
By comparing the above expression to Eq. (3.24) we see that A 3 (ω) can be written as

A 3 (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ , x 3 ′ ; ω) = n X n (x 1 , x 2 , x 3 ′ )X * n (x 1 ′ , x 2 ′ , x 3 )δ(ω -(E N +1 n -E N 0 )) + n Z n (x 1 , x 2 , x 3 ′ )Z * n (x 1 ′ , x 2 ′ , x 3 )δ(ω -(E N 0 -E N -1 n )). (3.31)
It is easy to show that A 3 (ω) is hermitian and positive definite. However, it is not normalized. Our goal is to calculate the spectral function of G 1 . We therefore require an equation that yields G 1 from G e+h

3

. Such a relation can be obtained by contracting the position-spin variables of the field operators that correspond to electron-hole pairs followed by integration over the contracted variables, i.e.

dx 2 dx 3 dx 2 ′ dx 3 ′ δ(x 2 -x 3 ′ )δ(x 3 -x 2 ′ )G e+h 3 (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ , x 3 ′ ; ω) = = dx 2 dx 3 G e+h 3 (x 1 , x 2 , x 3 , x 1 ′ , x 3 , x 2 ; ω) = dx 2 dx 3 n ⟨Ψ N 0 | ψ † (x 2 ) ψ(x 2 ) ψ(x 1 )|Ψ N +1 n ⟩⟨Ψ N +1 n | ψ † (x 1 ′ ) ψ † (x 3 ) ψ(x 3 )|Ψ N 0 ⟩ ω -(E N +1 n -E N 0 ) + iη + ⟨Ψ N 0 | ψ † (x 1 ′ ) ψ † (x 3 ) ψ(x 3 )|Ψ N -1 n ⟩⟨Ψ N -1 n | ψ † (x 2 ) ψ(x 2 ) ψ(x 1 )|Ψ N 0 ⟩ ω -(E N 0 -E N -1 n ) -iη = n N 2 ⟨Ψ N 0 | ψ(x 1 )|Ψ N +1 n ⟩⟨Ψ N +1 n | ψ † (x 1 ′ )|Ψ N 0 ⟩ ω -(E N +1 n -E N 0 ) + iη + (N -1) 2 ⟨Ψ N 0 | ψ † (x 1 ′ )|Ψ N -1 n ⟩⟨Ψ N -1 n | ψ(x 1 )|Ψ N 0 ⟩ ω -(E N 0 -E N -1 n ) -iη = N 2 G e 1 (x 1 , x 1 ′ ; ω) + (N -1) 2 G h 1 (x 1 , x 1 ′ ; ω) (3.32)
where we used that

dx ψ † (x) ψ(x)|Ψ N n ⟩ = N |Ψ N n ⟩, (3.33)
Theory 51 which leads to

G e 1 (x 1 , x 1 ′ , ω) = 1 N 2 dx 2 dx 3 G e 3 (x 1 , x 2 , x 3 , x 1 ′ , x 3 , x 2 ; ω) (3.34) G h 1 (x 1 , x 1 ′ , ω) = 1 (N -1) 2 dx 2 dx 3 G h 3 (x 1 , x 2 , x 3 , x 1 ′ , x 3 , x 2 ; ω). (3.35)
These relations show that the 1-GF is included in G e+h

3

. We notice that, the normalization prefactor depends on the choice of time differences in Eq. (3.23) As for the 1-GF, the definition of the 3-GF in Eq. (3.1) is not useful for practical calculations since its evaluation requires the knowledge of the N -body ground state wave function. Similarly, the expression of G e+h . Therefore, in the same spirit as for the 1-GF, we introduce a self-energy Σ 3 that is defined by the following Dyson equation

G e+h 3 (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ , x 3 ′ ; ω) = = G 0(e+h) 3 (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ , x 3 ′ ; ω) + G 0(e+h) 3 (x 1 , x 2 , x 6 , x 4 ′ , x 5 ′ , x 3 ′ ; ω) × Σ 3 (x 4 ′ , x 5 ′ , x 6 ′ , x 4 , x 5 , x 6 ; ω)G e+h 3 (x 4 , x 5 , x 3 , x 1 ′ , x 2 ′ , x 6 ′ ; ω). (3.36)
To solve the above Dyson equation, the starting point is the non-interacting 3-GF. Starting from equation (3.3), taking into account the choice for the time differences in Eq. (3.23), and performing a Fourier transform with respect to τ given in Eq. (3.13), we obtain the following expression for G

0(e+h) 3 G 0(e+h) 3 (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ , x 3 ′ ; ω) = = dω ′ dω ′′ (2π) 2 G 0 1 (x 1 , x 1 ′ ; ω + ω ′ -ω ′′ )G 0 1 (x 2 , x 2 ′ ; ω ′′ )G 0 1 (x 3 , x 3 ′ ; ω ′ ) +G 0 1 (x 1 , x 2 ′ ; ω)G 0 1 (x 2 , x 3 ′ )G 0 1 (x 3 , x 1 ′ ) + G 0 1 (x 1 , x 3 ′ )G 0 1 (x 2 , x 1 ′ ; ω)G 0 1 (x 3 , x 2 ′ ) -G 0 1 (x 1 , x 1 ′ ; ω)G 0 1 (x 2 , x 3 ′ )G 0 1 (x 3 , x 2 ′ ) -G 0 1 (x 1 , x 3 ′ )G 0 1 (x 2 , x 2 ′ ; ω)G 0 1 (x 3 , x 1 ′ ) - dω ′ dω ′′ (2π) 2 G 0 1 (x 1 , x 2 ′ ; ω + ω ′ -ω ′′ )G 0 1 (x 2 , x 1 ′ ; ω ′′ )G 0 1 (x 3 , x 3 ′ ; ω ′ ), (3.37)
The 3-GF and the 4-GF where G 0 1 (ω) is defined in equation (2.46), and

G 0 1 (x 1 , x 1 ′ ) = G 0 1 (x 1 , x 1 ′ , τ → 0 -) = iγ(x 1 , x 1 ′ ) = i v ϕ v (x 1 )ϕ * v (x 1 ′ ). (3.38)
One can recognize two types of contributions on the right-hand side of Eq. (3.37). The first type contains a product of three non-interacting 1-GFs of which only one depends on the frequency. From Eqs. (2.46) and (3.38) it becomes clear that these contributions correspond to quasi-particles since their poles correspond to a single eigenenergy. The second type contains two convolutions. Let us work out one of these contributions. We obtain

dω ′ dω ′′ (2π) 2 G 0 1 (x 1 , x 1 ′ ; ω + ω ′ -ω ′′ )G 0 1 (x 2 , x 2 ′ ; ω ′′ )G 0 1 (x 3 , x 3 ′ ; ω ′ ) = v c,c ′ ϕ c (x 1 )ϕ * c (x 1 ′ )ϕ c ′ (x 2 )ϕ * c ′ (x 2 ′ )ϕ v (x 3 )ϕ * v (x 3 ′ ) ω -ϵ c -(ϵ c ′ -ϵ v ) + iη + v,v ′ c ϕ v (x 1 )ϕ * v (x 1 ′ )ϕ v ′ (x 2 )ϕ * v ′ (x 2 ′ )ϕ c (x 3 )ϕ * c (x 3 ′ ) ω -ϵ v + (ϵ c -ϵ v ′ ) -iη . (3.39)
From the above expression, we see that the poles of this contribution correspond to the sum of an eigenenergy and an eigenenergy difference of a conduction and valence state. This shows that G e+h 3 already contains information about satellites in the noninteracting limit. Therefore, even with only a static 3-body self-energy the resulting 1-GF (obtained from Eq. (3.32)) will, in general, include satellites. The main task of a static 3-body self-energy is to modify the position (and spectral weight) of the poles, both due to quasiparticles and satellites, and bring them closer to the exact removal and addition energies. For these reasons, in the following we will focus on a static 3-body self-energy.

The 3-GF in a basis set

For practical calculations we project the 3-GF in the basis that diagonalizes G 0 1 . With the transformation (2.22), it becomes

G e+h ijl;mok (ω) = dx 1 dx 2 dx 3 dx 1 ′ dx 2 ′ dx 3 ′ ϕ * i (x 1 )ϕ * j (x 2 )ϕ l (x 3 ′ )ϕ m (x 1 ′ )ϕ o (x 2 ′ )ϕ * k (x 3 )G e+h (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ , x 3 ′ ; ω). (3.40)
Its spectral representation is

G e+h 3(ijl;mok) (ω) = n X ijl n X † mok n ω -(E N +1 n -E N 0 ) + iη + n Z ijl n Z † mok n ω -(E N 0 -E N -1 n ) -iη , (3.41)
Theory 53 in which

X ijl n = ⟨Ψ N 0 |ĉ † l ĉj ĉi |Ψ N +1 n ⟩ X † mok n = ⟨Ψ N +1 n |c † m ĉ † o ĉk |Ψ N 0 ⟩ Z ijl n = ⟨Ψ N -1 n |ĉ † l ĉj ĉi |Ψ N 0 ⟩ Z † mok n = ⟨Ψ N 0 |c † m ĉ † o ĉk |Ψ N -1 n ⟩. (3.42)
In particular, the non-interacting part is equal to

G 0(e+h) ijl,mok (ω) = [G 0 im G 0 jo G 0 lk ](ω) + G 0 io (ω)G 0 jl G 0 mk + G 0 jm (ω)G 0 il G 0 ok -G 0 im (ω)G 0 jl G 0 ok -G 0 jo (ω)G 0 il G 0 mk -[G 0 io G 0 jm G 0 lk ](ω) (3.43)
where G 0 ij (ω) is defined in equation (2.47), and

G 0 im = iδ im f i (3.44) [G 0 im G 0 jo G 0 lk ](ω) = δ im δ jo δ lk (f i -f l )(f j -f l ) ω -ϵ 0 i -(ϵ 0 j -ϵ 0 l ) + iηsign(ϵ 0 i -µ) . (3.45)
Here the occupation numbers (f i -f l )(f j -f l ) restrict G 0 3 to its 2h1e and 2e1h contributions. This means that the indices i, j, m, o refer to conduction (valence) states while l and k refer to valence (conduction) states. Therefore, the indices i, j, m, o can be associated to a pp process, and the same symmetry relation of Eq. (2.133) holds

G e+h ijl,mok = -G e+h jil,mok = -G e+h ijl,omk = G e+h jil,omk . (3.46)
The conditions i > j and m > o are applied to avoid double counting.

Let us analyze separately the quasi-particle and the satellite contributions. We have four quasi-particle terms, which yield

G 0 io (ω)G 0 jl G 0 mk + G 0 jm (ω)G 0 il G 0 ok -G 0 im (ω)G 0 jl G 0 ok -G 0 jo (ω)G 0 il G 0 mk = δ jl f j f k δ im δ ok -δ io δ mk ω -ϵ 0 i + iηsign(ϵ 0 i -µ) + δ il f i f k δ jo δ mk -δ jm δ ok ω -ϵ 0 j + iηsign(ϵ 0 j -µ) . (3.47)
This expression contains redundant information, i.e. the poles of G 0 1 repeated multiple times, which makes its matrix representation non-invertible. To avoid this problem we can reduce the space of definition of the quasi-particle contribution to the one of G 0 im (ω) only, without losing information. The sum of the two satellite terms is

[G 0 im G 0 jo G 0 lk ](ω) -[G 0 io G 0 jm G 0 lk ](ω) = [δ im δ jo -δ io δ jm ]δ lk (f i -f l )(f j -f l ) ω -ϵ 0 i -(ϵ 0 j -ϵ 0 l ) + iηsign(ϵ 0 i -µ) . (3.48)
Thanks to the conditions i > j and m > o, the deltas of the second term δ io δ jm are never fulfilled. Therefore, the satellite contribution is completely described by the first term.

The 3-GF and the 4-GF namely

dx 1 dx 4 dx 1 ′ dx 4 ′ δ(x 1 -x 1 ′ )δ(x 4 -x 4 ′ )G eh 4 (x 1 , x 2 , x 3 , x 4 , x 1 ′ , x 2 ′ , x 3 ′ , x 4 ′ ; ω) = dx 1 dx 4 G eh 4 (x 1 , x 2 , x 3 , x 4 , x 1 , x 2 ′ , x 3 ′ , x 4 ; ω) = dx 1 dx 4 i n ⟨Ψ N 0 | ψ(x 1 ) ψ(x 2 ) ψ † (x 2 ′ ) ψ † (x 1 )|Ψ N n ⟩⟨Ψ N n | ψ(x 4 ) ψ(x 3 ) ψ † (x 3 ′ ) ψ † (x 4 )|Ψ N 0 ⟩ ω -(E N n -E N 0 ) + iη -i n ⟨Ψ N 0 | ψ(x 4 ) ψ(x 3 ) ψ † (x 3 ′ ) ψ † (x 4 )|Ψ N n ⟩⟨Ψ N n | ψ(x 1 ) ψ(x 2 ) ψ † (x 2 ′ ) ψ † (x 1 )|Ψ N 0 ⟩ ω + (E N n -E N 0 ) -iη =(1 -N ) 2 i n ⟨Ψ N 0 | ψ(x 2 ) ψ † (x 2 ′ )|Ψ N n ⟩⟨Ψ N n | ψ(x 3 ) ψ † (x 3 ′ )|Ψ N 0 ⟩ ω -(E N n -E N 0 ) + iη -i n ⟨Ψ N 0 | ψ(x 3 ) ψ † (x 3 ′ )|Ψ N n ⟩⟨Ψ N n | ψ(x 2 ) ψ † (x 2 ′ )|Ψ N 0 ⟩ ω + (E N n -E N 0 ) -iη =(N -1) 2 G eh 2 (x 2 , x 3 , x 2 ′ , x 3 ′ ; ω). (3.58)
Therefore

G eh 2 (x 2 , x 3 , x 2 ′ , x 3 ′ , ω) = 1 (N -1) 2 dx 1 dx 4 G eh 4 (x 1 , x 2 , x 3 , x 4 , x 1 , x 2 ′ , x 3 ′ , x 4 , ω). (3.59)
Let us now analyze the non-interacting contribution with the choice of times in Eq. (3.53). Due to the large amount of terms created by equation (3.52), we will give only some examples. Each term is composed of four non-interacting 1-GF. They can depend on a time difference or be instantaneous. Depending on the number of time-dependent 1-GF, each term describes the propagation of 4, 2, or 0 particles.

When each G 0 1 in the group of four depends on a time difference, the term describes the propagation of four particles. An example is the term

G 0 1 (1, 4 ′ )G 0 1 (2, 3 ′ )G 0 1 (3, 2 ′ )G 0 1 (4, 1 ′ ) = = G 0 1 (x 1 , x 4 ′ ; τ )G 0 1 (x 2 , x 3 ′ ; τ )G 0 1 (x 3 , x 2 ′ ; -τ )G 0 1 (x 4 , x 1 ′ ; -τ ), (3.60)
Theory 57 which, in its spectral representation, becomes

[G 0 1 (x 1 , x 4 ′ )G 0 1 (x 2 , x 3 ′ )G 0 1 (x 3 , x 2 ′ )G 0 1 (x 4 , x 1 ′ )](ω) ≡ ≡G 0(eh) 4 (x 1 , x 2 , x 3 , x 4 , x 1 ′ , x 2 ′ , x 3 ′ , x 4 ′ ; ω) = = -i c,c ′ ,v,v ′ ϕ c (x 1 )ϕ * c (x 4 ′ )ϕ c ′ (x 2 )ϕ * c ′ (x 3 ′ )ϕ v ′ (x 3 )ϕ * v ′ (x 2 ′ )ϕ v (x 4 )ϕ * v (x 1 ′ ) ω -(ϵ 0 c -ϵ 0 v ) -(ϵ 0 c ′ -ϵ 0 v ′ ) + iη - c,c ′ ,v,v ′ ϕ v (x 1 )ϕ * v (x 4 ′ )ϕ v ′ (x 2 )ϕ * v ′ (x 3 ′ )ϕ c ′ (x 3 )ϕ * c ′ (x 2 ′ )ϕ c (x 4 )ϕ * c (x 1 ′ ) ω + (ϵ 0 c -ϵ 0 v ) + (ϵ 0 c ′ -ϵ 0 v ′ ) -iη . (3.61)
From this expression, we clearly see that the poles of this term correspond to a sum of two eigenenergy differences of a conduction and valence state. This shows the description of double neutral excitations already at the non-interacting level.

Terms in which two G 0 1 depend on time, while the other two are static, describe the propagation of two particles. For example in the term

G 0 1 (1, 1 ′ )G 0 1 (2, 3 ′ )G 0 1 (3, 2 ′ )G 0 1 (4, 4 ′ ) = = G 0 1 (x 1 , x 1 ′ ; 0 + )G 0 1 (x 2 , x 3 ′ ; τ )G 0 1 (x 3 , x 2 ′ ; -τ )G 0 1 (x 4 , x 4 ′ ; 0 + ), (3.62)
the two time dependent G 0 1 define the poles of the 4-GF in accordance with the eh channel of the 2-GF (see Eq. (2.79)).

The terms where all the G 0 1 are static, like for example the term

G 0 1 (1, 1 ′ )G 0 1 (2, 2 ′ )G 0 1 (3, 3 ′ )G 0 1 (4, 4 ′ ) = = G 0 1 (x 1 , x 1 ′ ; 0 + )G 0 1 (x 2 , x 2 ′ ; 0 + )G 0 1 (x 3 , x 3 ′ ; 0 + )G 0 1 (x 4 , x 4 ′ ; 0 + ), (3.63)
do not describe the propagation of particles. This terms do not have poles and, therefore, they can be excluded.

The 4-GF in a basis set

Analogously to the 3-GF case, projecting the 4-GF on the basis that diagonalizes G 0 1 is very important to keep only the essential information and avoid double counting.

Using the following transformation

G eh ijln;mokp (ω) = dx 1 dx 2 dx 3 dx 4 dx 1 ′ dx 2 ′ dx 3 ′ dx 4 ′ ϕ * i (x 1 )ϕ * j (x 2 )ϕ l (x 2 ′ )ϕ n (x 1 ′ ) × G eh 4 (x 1 , x 2 , x 3 , x 4 , x 1 ′ , x 2 ′ , x 3 ′ , x 4 ′ ; ω)ϕ m (x 4 ′ )ϕ o (x 3 ′ )ϕ * k (x 3 )ϕ * p (x 4 ), (3.64)
where

G 0 n (ω) =      G 0 i1,j1 (ω) 0 ... 0 0 [G 0 i1,j1 G 0 i2,j2 G 0 i3,j3 ](ω) 0 ... ... ... 0 ... 0 0 ... 0 [G 0 i1,j1 ...G 0 in,jn ](ω)      , (4.2)
if n ≥ 3 is odd, while

G 0 n (ω) =     [G 0 i 1 ,j 1 G 0 i 2 ,j 2 ](ω) 0 ... 0 0 [G 0 i 1 ,j 1 G 0 i 2 ,j 2 G 0 i 3 ,j 3 G 0 i 4 ,j 4 ](ω) 0 ... ... ... 0 ... 0 0 ... 0 [G 0 i 1 ,j 1 ...G 0 in,jn ](ω)     , (4.3)
if n ≥ 4 is even. We notice that these matrices have to be defined in a proper space to avoid double counting and to make the matrices invertible. Moreover, each term must have the same poles in its spectral representation, like in the equations (2.45) and (3.24) or equations (2.79) and (3.55).

The multi-channel self-energy is defined as

Σ n =          Σ 1p Σ 1p/
(n-2)p Σ (n-2)p/np 0 ... ... 0 Σnp/(n-2)p Σ np          (4.4) if n is odd, while Σ n =          Σ 2p Σ 2p/4p 0 ... ... 0 Σ4p/2p Σ 4p Σ 4p/6p 0 ... ... 0 Σ6p/4p Σ 6p Σ 6p/8p 0 ... ... 0 ... ... ... 0 ... ... 0 Σ(n-2)p/(n-4)p Σ (n-2)p Σ (n-2)p/np 0 ... ... 0 Σnp/(n-2)p Σ np          (4.5)
if n is even.

Here, for notational convenience, we do not report the indices of each term. In general, they depend on the number of particles they are related to, i.e., Σ (n-2)p/np ≡ Σ (n-2)p/np i1...in-2;j1...jn , Σnp/(n-2)p ≡ Σ np/(n-2)p i1...in;j1...jn-2 , and Σ np ≡ Σ np i1...in;j1...jn

Multi-channel 61

Each term in the multi-channel self-energy correlates particle pairs at the RPA+exchange (RPAx) level. This is an approximation for the multi-channel selfenergy. We will explicitly show the multi-channel self-energy for the 1-3-GF in section 4.2 and the one for the 2-4-GF in section 4.3. Knowing how to construct these terms allows the reader to construct subsequent multi-channel self-energy terms.

This strategy follows a well-known idea: by starting from a large space, a static self-energy is able to capture the dynamical properties of a smaller space, cite [START_REF] Berger | Potential energy surfaces without unphysical discontinuities: The coulomb hole plus screened exchange approach[END_REF][START_REF] Schirmer | New approach to the one-particle green's function for finite fermi systems[END_REF][START_REF] Von Niessen | Computational methods for the one-particle green's function[END_REF].

Multi-channel Dyson equation for the 1-and 3-GF

In this section, we will show how to construct the multi-channel self-energy for the 1-3-GF channel. We will analyze which diagrams of the 1-GF we obtain from this procedure.

Our starting point is the multi-channel Dyson equation (4.1) with n = 3

G 3 (ω) = G 0 3 (ω) + G 0 3 (ω)Σ 3 G 3 (ω), (4.6)
where G 0 3 (ω) is the noninteracting 3-GF (3.49) and Σ 3 is the multi-channel self-energy. To have all the information needed to understand this equation in this section, we report the noninteracting term already discussed at the end of section 3.1. It reads

G 0 3 (ω) = G 0 1 (ω) 0 0 G 0(3p) 3 (ω) , (4.7)
where

G 0 im (ω) = δ im ω -ϵ 0 i + iηsign(ϵ 0 i -µ) , (4.8)
G 0(3p) i>jl;m>ok (ω) = δ im δ jo δ lk (f i -f l )(f j -f l ) ω -ϵ 0 i -(ϵ 0 j -ϵ 0 l ) + iηsign(ϵ 0 i -µ) . (4.9)
It becomes clear that the upper left block of Eq. (4.7) is in the electron-hole space and the lower right block is in the 2e1h and 2h1h space. Without loss of generality, the non-interacting 3-GF can be replaced with an independent-particle 3-GF. It is convenient to use a Hartree-Fock 3-GF since it already contains Hartree and exchange contributions. In this way, the three-body self-energy in Eq. (4.6) is reduced to its correlation part. Therefore, in the following, G 0 3 will denote the Hartree-Fock 3-GF. The multi-channel self-energy is defined as

Σ 3 = Σ 1p Σ c Σc Σ 3p . (4.10)

Approximation to Σ 3

For practical calculations, we need an approximation to Σ 3 . To achieve this we correlate only pairs of particles in Σ 3p . We let each pair interact via a direct and an exchange interaction, i.e., at the RPA+exchange (RPAx) level, as in the electron-hole and particle-particle channels of the 2-GF [START_REF] Leo | Quantum statistical mechanics[END_REF][START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF][START_REF] Strinati | Application of the green's functions method to the study of the optical properties of semiconductors[END_REF]. The four-point couplings Σ c and Σc correspond to two-particle channels which, for consistency, are also treated at the RPAx level. Since all the correlation is included in Σ c , Σc , and Σ 3p , the head, Σ 1p , vanishes.

We thus arrive at the following static approximation,

Σ 3p ijl;mok = [(1-f i )(1-f j )f l -f i f j (1-f l )][δ lk vijom +δ mj viklo +δ io vjklm -δ oj viklm -δ im vjklo ], (4.11)
Σ c i;mok = vikom , (4.12) Σc ijl;m = vijlm , (4.13)
Σ 1p i;m = 0, (4.14)
where vikom = v ikom -v ikmo , with the bare Coulomb potential defined in Eq. (2.136).

The first term on the right-hand side of Eq. (4.11), i.e, the one involving δ lk vijom , accounts for all the two-particle interactions in Σ 3p , while the last four terms account for all the one-electron-one-hole interactions. The occupation numbers in Eq. (4.11) ensure that Σ 3p is block diagonal and it has opposite signs for the 2e1h channel and the 2h1e channel. We note that our approximate three-body self-energy is hermitian.

It is interesting to mention that in Eqs. (4.11)-(4.13) can also be obtained using a technique similar to the adiabatic diagrammatic construction [START_REF] Bintrim | Full-frequency GW without frequency[END_REF][START_REF] Schirmer | The two-particle-hole tamm-dancoff approximation (2ph-tda) equations for closed-shell atoms and molecules[END_REF][START_REF] Schirmer | New approach to the one-particle green's function for finite fermi systems[END_REF][START_REF] Von Niessen | Computational methods for the one-particle green's function[END_REF].

Diagrammatic analysis

It is instructive to represent the multi-channel Dyson equation in Eq. (4.6) diagrammatically according to (4.15) From Eq. (4.15) we see that the 3-GF contains the 1-GF ( ) plus an explicit three-body part G 3p 3 ( ) and the couplings between the 1-GF and G 3p 3 , namely G c ( ) and Gc 3 ( ). The G 3p 3 , G c 3 and Gc 3 terms add correlation to the non-interacting 1-GF ( ) and 3-GF ( ). We represent Σ c ijl,m and Σc i,mok by right triangles to reflect their dimensions.

To analyze the diagrammatic structure of the multi-channel Dyson equation in Eq. (4.15), it is convenient to represent it in real space, since it is the space where Feynman diagrams were discussed in chapter 2.2.6 and 2.3.4. Applying the change of basis of Eq. (3.40) to the multi-channel self-energy, we obtain

Σ 3p ijl;mok = dx 1 dx 2 dx 3 dx 1 ′ dx 2 ′ dx 3 ′ Σ 3p (x 1 , x 2 , x 3 ′ , x 1 ′ , x 2 ′ , x 3 ) × ϕ * i (x 1 )ϕ * j (x 2 )ϕ l (x 3 ′ )ϕ m (x 1 ′ )ϕ o (x 2 ′ )ϕ * k (x 3 ), (4.16) Σ c i;mok = dx 1 dx 3 dx 1 ′ dx 2 ′ Σ c (x 1 , x 1 ′ , x 2 ′ , x 3 )ϕ * i (x 1 )ϕ m (x 1 ′ )ϕ o (x 2 ′ )ϕ * k (x 3 ), (4.17) Σc ijl;m = dx 1 dx 2 dx 1 ′ dx 3 ′ Σc (x 1 , x 2 , x 3 ′ , x 1 ′ )ϕ * i (x 1 )ϕ * j (x 2 )ϕ l (x 3 ′ )ϕ m (x 1 ′ ). (4.18)
The multi-channel self-energy becomes

Σ 2e1h (x 1 , x 2 , x 3 ′ , x 1 ′ , x 2 ′ , x 3 ) = -Σ 2h1e (x 1 , x 2 , x 3 ′ , x 1 ′ , x 2 ′ , x 3 ) = δ(x 3 , x 3 ′ )[δ(x 1 , x 1 ′ )δ(x 2 , x 2 ′ ) -δ(x 1 , x 2 ′ )δ(x 2 , x 1 ′ )]v(r 1 , r 2 ) + δ(x 1 ′ , x 2 )[δ(x 1 , x 2 ′ )δ(x 3 , x 3 ′ ) -δ(x 1 , x 3 ′ )δ(x 3 , x 2 ′ )]v(r 1 , r 3 ) + δ(x 1 , x 2 ′ )[δ(x 2 , x 1 ′ )δ(x 3 , x 3 ′ ) -δ(x 2 , x 3 ′ )δ(x 3 , x 1 ′ )]v(r 2 , r 3 ) -δ(x 2 , x 2 ′ )[δ(x 1 , x 1 ′ )δ(x 3 , x 3 ′ ) -δ(x 1 , x 3 ′ )δ(x 3 , x 1 ′ )]v(r 1 , r 3 ) -δ(x 1 , x 1 ′ )[δ(x 2 , x 2 ′ )δ(x 3 , x 3 ′ ) -δ(x 2 , x 3 ′ )δ(x 3 , x 2 ′ )]v(r 2 , r 3 ) (4.19) Σ c (x 1 , x 1 ′ , x 2 ′ , x 3 ) = [δ(x 1 , x 1 ′ )δ(x 2 ′ , x 3 ) -δ(x 1 , x 2 ′ )δ(x 1 ′ , x 3 )]v(r 1 , r 3 ) (4.20) Σc (x 1 , x 2 , x 3 ′ , x 1 ′ ) = [δ(x 1 , x 1 ′ )δ(x 2 , x 3 ′ ) -δ(x 1 , x 3 ′ )δ(x 2 , x 1 ′ )]v(r 1 , r 2 ). (4.21)
From these expressions in real space, it becomes easier to understand the diagrammatic structure. The 2e1h body part (4. [START_REF] Viktor | Application of quantum field theory methods to the many body problem[END_REF]) is given by .

(4.22) The 2h1e diagrams are equal to the 2e1h diagrams in Eq. (4.22) but with an overall minus sign. In the diagrams given above, each dotted line represents a Dirac delta, To understand which diagrams of G 1 are included at each order in our approach, it suffices to iterate the multi-channel Dyson equation in Eq. (4.15) and inspect the head of the matrix. One iteration does not change the head, i.e., no correlation is added to G 0 1 . A second and a third iteration yield the two self-energy insertions depicted in Fig. 4.1.

In general, the n-th order proper skeleton diagrams can be obtained from the (n -1)th-order diagram by inserting between Σ c and Σc a Σ 3p rectangle linked to a Multi-channel Dyson equation and a ladder diagram as (4.28) Notice that, due to the condition i > j and m > o in Eq. (4.9) on each G 0 3 line, the sum of four different multi-channel diagrams, each composed of different multi-channel self-energy terms, are needed to obtain a third-order diagram. Here, for simplicity, only one of them is shown. A more detailed analysis is reported in appendix C.

The bubble diagram in Eq. (4.27) is also contained in the GW approximation which contains bubble diagrams up to infinite order. In our approach, higher-order diagrams are obtained by iterating the multi-channel Dyson equation (4.15) further. By doing so, it is possible to check that all GW diagrams are included in our approach, and screening effects are thus accounted for. Moreover, our approach goes beyond the RPA screening included in GW . As an example, we report the following fourth-order diagram (4.29) where a first-order vertex correction to the RPA bubble is shown. Although the self-energy only involves the bare Coulomb potential, screening beyond the RPA is naturally included in the multi-channel Dyson equation. Therefore, there is no need to use a screened interaction in our theory since this would lead to the double counting of diagrams. This principle is the same as in the exchange term of the eh BSE, reported in Figs. 2.14 and 2.16.

Finally, to understand the number of proper skeleton diagrams produced at each order by our theory, one has to consider the number of Σ 3 terms, i.e. two terms for each coupling term and ten terms for Σ 3p . Furthermore, due to the conditions i > j and m > o, for each G 0(3p) 3 present in a three-body diagram, one needs two combinations of Σ 3 to find one G 1 diagram. Therefore, the total number of proper

c G 0(3p) 3 Σc G 0 1 Σ c G 0(3p) 3
Σc (see Fig. 4.2). Therefore, each time a G 0 1 enters between two proper self-energy terms, repetitions of proper self-energy terms appear, exactly as in the Dyson equation (2.56). The non-skeleton diagrams are present due to G 0(3p) in Eq. (3.49). To improve them, one can dress G 0(3p) beyond HF using, e.g., second Born, GW or the T-matrix within the quasiparticle approximation. In this case, one should assume that the correlated G 0(3p) 3 is diagonal in the same basis as G 0 3 . This approximation is also used, for example, to derive the cumulant approximation [START_REF] Guzzo | Valence electron photoemission spectrum of semiconductors: Ab initio description of multiple satellites[END_REF][START_REF] Martin | Interacting Electrons: Theory and Computational Approaches[END_REF].

Effective three particle Hamiltonian

To obtain an equation that can be solved in practice using standard numerical tools, we use a strategy similar to the one used for the Bethe-Salpeter equation [START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF], i.e., we map Eq. (4.1) onto an effective three-particle Hamiltonian according to (4.30) in which the effective Hamiltonian H eff is given by

G 3(ijl;mok) = ωI -H eff -1 ijl;mok
H eff = H 1p Σ c Σc H 3p , (4.31)
where

H 1p i;m = ϵ 0 i δ im , (4.32)
H 3p ijl;mok = (ϵ 0 i -(ϵ 0 l -ϵ 0 j ))δ im δ jo δ lk + Σ 3p ijl;mok . (4.33)
Details of this derivation are shown in appendix D.

Multi-channel Dyson equation

Since our final goal is to calculate photoemission spectra, which are linked to the one-body-channel of G 3 , it suffices to extract the head from G 3 in Eq. (4.30).

A straightforward diagonalisation of Eq. (4.31) would scale as N 9 , where N is the number of electrons in the system. However, the scaling can be significantly reduced using standard iterative methods such as the Haydock-Lanczos solver [START_REF] Haydock | Electronic structure based on the local atomic environment for tight-binding bands[END_REF][START_REF] Hernandez | Slepc: A scalable and flexible toolkit for the solution of eigenvalue problems[END_REF][START_REF] Schmidt | Efficient O(N 2 ) method to solve the bethe-salpeter equation[END_REF]. With it, the scaling reduces to N 6 since the scaling is then determined by the construction of the effective Hamiltonian. We also note that methods that describe trions have a similar scaling as our approach, and they have successfully been applied to real systems [START_REF] Deilmann | Three-particle correlation from a many-body perspective: Trions in a carbon nanotube[END_REF][START_REF] Torche | First-principles many-body theory for charged and neutral excitations: Trion fine structure splitting in transition metal dichalcogenides[END_REF].

Multi-channel Dyson equation for the 2-and 4-GF

In this section, we introduce the multi-channel Dyson equation, and the corresponding self-energy, for the 2-4-body channel. By using the Feynman diagram technique, we will analyze which diagramsof the electron-hole 2-GF we obtain.

By restricting the multi-channel Dyson equation (4.1) to n = 4, we obtain

L(ω) = L0 (ω) + L0 (ω)Σ 4 L(ω) (4.34)
where L0 is the non-interacting correlation function (3.67), and Σ 4 is the multi-channel self-energy. For completeness, we report the non-interacting term

L0 (ω) = L0 2 (ω) 0 0 L0 4 (ω) , (4.35)
where

L0 jl,ok (ω) = δ jo δ lk (f j -f l ) ϵ 0 j -ϵ 0 l -ω + iηsign(f j -f l ) , (4.36) L0 i>jl>n,m>ok>p (ω) = δ im δ jo δ lk δ np (f i -f n )(f i -f l )(f j -f n ) (ϵ 0 i -ϵ 0 n ) + (ϵ 0 j -ϵ 0 l ) -ω + iηsign(f j -f l ) . ( 4

.37)

As for the 1-3-body case, the above non-interacting terms can be replaced with an independent-particle L0 . As we will show in the diagrammatic analysis, it is convenient to use Hartree-Fock orbitals and energies.

Approximation to Σ 4

The multi-channel self-energy is defined as

Σ 4 = Ξ Σ c Σc Σ 4p . (4.38)
The approximation that we use for it is analogous to the approximation we have used for the 1-3-body channel. We let each pair of particles interact at the RPAx level. where vijkm = v ijkm -v ijmk with v ijkm defined in Eq. (2.136). The first line of the multi-channel self-energy body (4.39), accounts for all the two-particle interactions, while all the other terms describe all the one-electron-one-hole interactions. We note that the multi-channel self-energy is hermitian.

Diagrammatic analysis

The

Σ c (x 2 , x 3 , x 4 , x 2 ′ , x 3 ′ , x 4 ′ ) = δ(x 2 , x 3 ′ )[δ(x 2 ′, x 3 )δ(x 4 ′, x 4 )-δ(x 2 ′, x 4 )δ(x 4 ′, x 3 )]v(r 2 ′,r 4 ′) +δ(x 2 ′, x 3 )[δ(x 2 , x 4 ′ )δ(x 4 , x 3 ′ )-δ(x 2 , x 3 ′ )δ(x 4 , x 4 ′ )]v(r 2 , r 4 ) +δ(x 2 ′, x 4 )[δ(x 2 , x 3 ′ )δ(x 3 , x 4 ′ )-δ(x 2 , x 4 ′ )δ(x 3 , x 3 ′ )]v(r 2 , r 3 ) +δ(x 2 , x 4 ′ )[δ(x 2 ′, x 4 )δ(x 3 ′, x 3 )-δ(x 2 ′, x 3 )δ(x 3 ′, x 4 )]v(r 2 ′,r 3 ′), (4.43) Σc (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ x 3 ′ ) = δ(x 2 , x 3 ′ )[δ(x 1 , x 1 ′ )δ(x 3 , x 2 ′ )-δ(x 1 , x 2 ′ )δ(x 3 , x 1 ′ )]v(r 1 , r 3 ) +δ(x 2 ′ , x 3 )[δ(x 1 ′, x 2 )δ(x 3 ′, x 1 )-δ(x 1 ′, x 1 )δ(x 3 ′, x 2 )]v(r 1 ′ , r 3 ′ ) +δ(x 1 ′ , x 3 )[δ(x 2 ′, x 1 )δ(x 3 ′, x 2 )-δ(x 2 ′, x 2 )δ(x 3 ′, x 1 )]v(r 2 ′ , r 3 ′ ) +δ(x 1 , x 3 ′ )[δ(x 2 , x 2 ′ )δ(x 3 , x 1 ′ )-δ(x 2 , x 1 ′ )δ(x 3 , x 2 ′ )]v(r 2 , r 3 ), (4.44) Σ 4p (x 1 , x 2 , x 3 , x 4 , x 1 ′ , x 2 ′ , x 3 ′ , x 4 ′ ) = = δ(x 1 , x 4 ′ )δ(x 2 , x 3 ′ )[δ(x 2 ′ , x 4 )δ(x 1 ′ , x 3 ) -δ(x 2 ′ , x 3 )δ(x 1 ′ , x 4 )]v(r 2 ′ , r 1 ′ ) + δ(x 2 ′ , x 3 )δ(x 1 ′ , x 4 )[δ(x 1 , x 3 ′ )δ(x 2 , x 4 ′ ) -δ(x 1 , x 4 ′ )δ(x 2 , x 3 ′ )]v(r 1 , r 2 ) + δ(x 2 ′ , x 4 )δ(x 1 ′ , x 3 )[δ(x 1 , x 4 ′ )δ(x 2 , x 3 ′ ) -δ(x 1 , x 3 ′ )δ(x 2 , x 4 ′ )]v(r 1 , r 2 ) + δ(x 1 , x 3 ′ )δ(x 2 , x 4 ′ )[δ(x 2 ′ , x 3 )δ(x 1 ′ , x 4 ) -δ(x 2 ′ , x 4 )δ(x 1 ′ , x 3 )]v(r 1 ′ , r 2 ′ ) + δ(x 1 , x 4 ′ )δ(x 1 ′ , x 4 )[δ(x 2 , x 3 ′ )δ(x 3 , x 2 ′ ) -δ(x 2 , x 2 ′ )δ(x 3 , x 3 ′ )]v(r 2 , r 3 ) + δ(x 1 , x 4 ′ )δ(x 2 ′ , x 3 )[δ(x 2 , x 3 ′ )δ(x 4 , x 1 ′ ) -δ(x 2 , x 1 ′ )δ(x 4 , x 3 ′ )]v(r 2 , r 4 ) + δ(x 2 , x 3 ′ )δ(x 1 ′ , x 4 )[δ(x 1 , x 4 ′ )δ(x 3 , x 2 ′ ) -δ(x 1 , x 2 ′ )δ(x 3 , x 4 ′ )]v(r 1 , r 3 ) + δ(x 2 , x 3 ′ )δ(x 2 ′ , x 3 )[δ(x 1 , x 4 ′ )δ(x 4 , x 1 ′ ) -δ(x 1 , x 1 ′ )δ(x 4 , x 4 ′ )]v(r 1 , r 4 ) -δ(x 1 , x 4 ′ )δ(x 2 ′ , x 4 )[δ(x 2 , x 3 ′ )δ(x 3 , x 1 ′ ) -δ(x 2 , x 1 ′ )δ(x 3 , x 3 ′ )]v(r 2 , r 3 ) -δ(x 1 , x 4 ′ )δ(x 1 ′ , x 3 )[δ(x 2 , x 3 ′ )δ(x 4 , x 2 ′ ) -δ(x 2 , x 2 ′ )δ(x 4 , x 3 ′ )]v(r 2 , r 4 ) -δ(x 2 , x 3 ′ )δ(x 2 ′ , x 4 )[δ(x 1 , x 4 ′ )δ(x 3 , x 1 ′ ) -δ(x 1 , x 1 ′ )δ(x 3 , x 4 ′ )]v(r 1 , r 3 ) -δ(x 2 , x 3 ′ )δ(x 1 ′ , x 3 )[δ(x 1 , x 4 ′ )δ(x 4 , x 2 ′ ) -δ(x 1 , x 2 ′ )δ(x 4 , x 4 ′ )]v(r 1 , r 4 ) -δ(x 1 , x 3 ′ )δ(x 1 ′ , x 4 )[δ(x 2 , x 4 ′ )δ(x 3 , x 2 ′ ) -δ(x 2 , x 2 ′ )δ(x 3 , x 4 ′ )]v(r 2 , r 3 ) -δ(x 1 , x 3 ′ )δ(x 2 ′ , x 3 )[δ(x 2 , x 4 ′ )δ(x 4 , x 1 ′ ) -δ(x 2 , x 1 ′ )δ(x 4 , x 4 ′ )]v(r 2 , r 4 ) -δ(x 2 , x 4 ′ )δ(x 1 ′ , x 4 )[δ(x 1 , x 3 ′ )δ(x 3 , x 2 ′ ) -δ(x 1 , x 2 ′ )δ(x 3 , x 3 ′ )]v(r 1 , r 3 ) -δ(x 2 , x 4 ′ )δ(x 2 ′ , x 3 )[δ(x 1 , x 3 ′ )δ(x 4 , x 1 ′ ) -δ(x 1 , x 1 ′ )δ(x 4 , x 3 ′ )]v(r 1 , r 4 ) + δ(x 1 , x 3 ′ )δ(x 2 ′ , x 4 )[δ(x 2 , x 4 ′ )δ(x 3 , x 1 ′ ) -δ(x 2 , x 1 ′ )δ(x 3 , x 4 ′ )]v(r 2 , r 3 ) + δ(x 1 , x 3 ′ )δ(x 1 ′ , x 3 )[δ(x 2 , x 4 ′ )δ(x 4 , x 2 ′ ) -δ(x 2 , x 2 ′ )δ(x 4 , x 4 ′ )]v(r 2 , r 4 ) + δ(x 2 , x 4 ′ )δ(x 2 ′ , x 4 )[δ(x 1 , x 3 ′ )δ(x 3 , x 1 ′ ) -δ(x 1 , x 1 ′ )δ(x 3 , x 3 ′ )]v(r 1 , r 3 ) + δ(x 2 , x 4 ′ )δ(x 1 ′ , x 3 )[δ(x 1 , x 3 ′ )δ(x 4 , x 2 ′ ) -δ(x 1 , x 2 ′ )δ(x 4 , x 3 ′ )]v(r 1 , r 4 ). (4.45)
The above expressions of the multi-channel self-energy in real space makes it easier to understand the diagrammatic structure. The coupling terms in Eqs. (4.43) 2.68)). Therefore, as for the 1-3-body case, the screening is naturally included in the theory although the self-energy only involves the bare Coulomb potential. Analyzing subsequent orders, one can realize that the screening goes beyond the RPA.

Multi-channel Dyson equation

In the third order in the interaction, also the body part of the multi-channel selfenergy (4.48) enters in the head of the matrix of G 4 . There are hence several terms which are created. We report an example in Fig. 4.8.

To improve non-skeleton diagrams, one can go beyond HF in the non-interacting 4-body term (4.37), i.e., using single-particle eigenvalues ϵ 0 obtained from second Born, GW or T-matrix quasiparticle calculations. Notice that, diagrams that add An example of a combination of multi-channel self-energy coupling terms that creates a mixing term. Here, the eh attraction is between a hole and an electron that is already dressed. correlation to a single G 0 1 line cannot be removed by improving the non-interacting two-body term (4.36) because they are needed to create also other terms. For example, the diagrams in Figs. 4.3 and 4.5 are created with the same coupling term

Multi-channel Dyson equation

Effective four-particle Hamiltonian

We recast the multi-channel Dyson equation (4.34) in an eigenvalue problem. Following the same strategy of the BSE [START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF], we rewrite the equation as Lijln;mokp = H eff -ωI -1 ijln;mokp (4.49) in which the effective Hamiltonian H eff is given by

H eff = H 2p Σ c Σc H 4p , (4.50)
where

H 2p jl;ok = (ϵ 0 j -ϵ 0 l )δ jo δ lk + Ξ jl;ok (f l -f j ), (4.51)
H 4p ijln;mokp = [(ϵ 0 i -ϵ 0 n )+(ϵ 0 j -ϵ 0 l )]δ im δ jo δ lk δ np +Σ 4p ijln;mokp (f n -f i )(f l -f i )(f n -f j ). (4.52)
To obtain the information related to the two-particle correlation function, after the inversion of the matrix in Eq. (4.49) one has to extract the head from L.

CHAPTER 5

Results: the symmetric Hubbard dimer

In order to illustrate the strategy discussed in the previous chapter we consider the symmetric Hubbard dimer. It consists of two degenerate sites each containing one orbital. Moreover, only electrons on the same site interact with each other. In this chapter we will test the multi-channel Dyson-equation for the 1-3-GF. We will show that the approximate multi-channel self-energy (see Eqs. The Hamiltonian corresponding to the symmetric Hubbard dimer is given by

H = -t i̸ =j=1,2 σ ĉ † iσ ĉjσ + U 2 i=1,2 σσ ′ ĉ † iσ ĉ † iσ ′ ĉiσ ′ ĉiσ + ϵ 0 i=1,2 σ niσ , (5.1)
in which -t, U and ϵ 0 represent the hopping kinetic energy, the (spin-independent) on-site interaction and the orbital energy, respectively, and n iσ = c † iσ ĉiσ is the number operator. We made explicit the spin σ index in the above equation. We note that the amount of electron correlation in the system is proportional to the ratio U/t. We consider weak correlation when U/t = 1, while when U/t = 4 we are in the strong correlation regime. The model is exactly solvable and, therefore, allows us to test the accuracy of our results in both the weakly and strongly correlated regimes. In appendix E we report the eigenvalues and eigenfunctions of the above Hamiltonian for the cases of one, two and three electrons in the bonding/antibonding basis. From these ingredients we can calculate the exact 1-and 3-GF from their definitions, see Eqs. (2.45) and (3.24), respectively. From the 1-GF we will plot the spectral function

The Hubbard dimer

A(ω) = i A ii (ω) = i 1 π sign(µ -ω)ImG ii (ω). (5.2)
More details of the exact 1-GF are given in, e.g., Refs. [START_REF] Di | Reduced density-matrix functional theory: Correlation and spectroscopy[END_REF][START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF].

To test the accuracy of the multi-channel Dyson equation (4.1), we will calculate the multi-channel self-energy in two ways: 1) from the exact results by inverting Eq. (4.1), namely

Σ 3 = [G 0 3 ] -1 -G -1 3 ;
2) from the approximations in Eqs. (4.11)-(4.14). We will show that the two methods give the same results. Therefore, the multichannel Dyson equation with the approximate multi-channel self-energy in Eqs. (4.11)-(4.14) is exact for the symmetric Hubbard dimer.

1/4 filling

In the bonding/antibonding basis the ground state of the 1/4 filling regime is the state |b ↑⟩ (we could have equivalently chosen |b ↓⟩ as ground state). For this reason the removal part of the multi-channel 3-GF in Eq. (4.1) is composed of only one element: the 1-GF matrix element where the b ↑ electron is removed. The addition part is composed of three 1-GF and three 3-GF matrix elements. The 1-GF elements describe the addition of an electron to the state b ↓, a ↓, or a ↑. Due to the restriction of the three-particle space defined in Eq. (4.9), namely G 3p i>jl;m>ok , the only allowed

Exact results

To construct the exact multi-channel 3-GF we use the eigenvalues and eigenvectors of the symmetric Hubbard dimer (see App. E). We obtain

G 3 =                b ↑ b ↓ a ↓ a ↑ a ↓ a ↑ b ↑ a ↓ b ↓ b ↑ a ↑ b ↓ b ↑ b ↑ G b↑;b↑ 0 0 0 0 0 0 b ↓ 0 G b↓;b↓ 0 0 G b↓;a↓a↑b↑ 0 0 a ↓ 0 0 G a↓;a↓ 0 0 0 G a↓;a↑b↓b↑ a ↑ 0 0 0 G a↑;a↑ 0 0 0 a ↓ a ↑ b ↑ 0 G a↓a↑b↑;b↓ 0 0 G a↓a↑b↑;a↓a↑b↑ 0 0 a ↓ b ↓ b ↑ 0 0 0 0 0 G a↓b↓b↑;a↓b↓b↑ 0 a ↑ b ↓ b ↑ 0 0 G a↑b↓b↑;a↓ 0 0 0 G a↑b↓b↑;a↑b↓b↑                , (5.3)
where the diagonal elements are

G b↑;b↑ (ω) = 1 ω -(ϵ 0 -t) -iη , G b↓;b↓ (ω) = (1 -A) 2 /a 2 ω -(ϵ 0 + (U -c)/2 + t) + iη + (1 -B) 2 /b 2 ω -(ϵ 0 + (U + c)/2 + t) + iη , G a↑;a↑ (ω) = 1 ω -(ϵ 0 + t) + iη , G a↓;a↓ (ω) = 1 2 1 ω -(ϵ 0 + t) + iη + 1 ω -(ϵ 0 + U + t) + iη , G a↓a↑b↑;a↓a↑b↑ (ω) = (1 + A) 2 /a 2 ω -(ϵ 0 + (U -c)/2 + t) + iη + (1 + B) 2 /b 2 ω -(ϵ 0 + (U + c)/2 + t) + iη , G a↓b↓b↑;a↓b↓b↑ (ω) = 1 ω -(ϵ 0 + t) + iη , G a↑b↓b↑;a↑b↓b↑ (ω) = 1 2 1 ω -(ϵ 0 + t) + iη + 1 ω -(ϵ 0 + U + t) + iη , (5.4)
and the off-diagonal elements are

G b↓;a↓a↑b↑ (ω) = G a↓a↑b↑;b↓ (ω) = = (1 -A 2 )/a 2 ω -(ϵ 0 + (U -c)/2 + t) + iη + (1 -B 2 )/b 2 ω -(ϵ 0 + (U + c)/2 + t) + iη , G a↓;a↑b↓b↑ (ω) = G a↑b↓b↑;a↓ (ω) = 1 2 1 ω -(ϵ 0 + t) + iη - 1 ω -(ϵ 0 + U + t) + iη , (5.5)
The Hubbard dimer

in which c = √ 16t 2 + U 2 , a = √ 2 c-U 16t 2 + (c -U ) 2 , A = 4t U -c , b = √ 2 c+U 16t 2 + (c + U ) 2 and B = 4t U +c .
We observe that the addition part of G 3 contains four distinct poles for U > 0. Three of them are quasi-particle peaks, while the one at energy ϵ 0 + (U + c)/2 + t is a satellite peak since it collapses to ϵ 0 + 3t at U = 0, which deviates from the antibonding peak ϵ 0 + t by 2t, which is the pole of the RP A polarizability. Moreover, this same pole has zero intensity at U = 0 in the 1-body channel (B = 1 when the interaction U is zero). Notice that this satellite peak has a non-vanishing amplitude for U = 0 in the 3-body channel, while all the off-diagonal elements go to zero. The exact non-interacting G 3 hence reads

G 0 3 (ω) = diag(G 0 b↑;b↑ (ω), G 0 b↓;b↓ (ω), G 0 a↓;a↓ (ω), G 0 a↑;a↑ (ω), G 0 a↓a↑b↑;a↓a↑b↑ (ω), G 0 a↓b↓b↑;a↓b↓b↑ (ω), G 0 a↑b↓b↑;a↑b↓b↑ (ω)); (5.6)
with

G 0 b↑;b↑ (ω) = 1 ω -(ϵ 0 -t) -iη , G 0 b↓;b↓ (ω) = 1 ω -(ϵ 0 -t) + iη , G 0 a↓;a↓ (ω) = G 0 a↑;a↑ (ω) = G 0 a↓b↓b↑;a↓b↓b↑ (ω) = G 0 a↑b↓b↑;a↑b↓b↑ (ω) = 1 ω -(ϵ 0 + t) + iη , G 0 a↓a↑b↑;a↓a↑b↑ (ω) = 1 ω -(ϵ 0 + 3t) + iη . (5.7)
In Fig. 5.2 we show a comparison between the spectral function corresponding to the non-interacting 1-GF and the spectral function corresponding to the noninteracting multi-channel 3-GF function G 0 3 , defined in Eqs. (5.6) and (5.7). One can see that, in the non-interacting limit, the peak at the satellite position is present with a nonvanishing spectral weight only in the spectral function obtained from the G 0 3 . The non-interacting 1-GF is retrieved from G 0 3 by taking the head of the matrix in Eq. ( 5.3).

To compare with the multi-channel self-energy in Eqs.

The Hubbard dimer

given by

G HF b↑;b↑ (ω) = 1 ω -(ϵ 0 -t) -iη , G HF b↓;b↓ (ω) = 1 ω -(ϵ 0 -t + U/2) + iη , G HF a↓;a↓ (ω) = G HF a↑b↓b↑;a↑b↓b↑ (ω) = 1 ω -(ϵ 0 + t + U/2) + iη , G HF a↑;a↑ (ω) = 1 ω -(ϵ 0 + t) + iη , G HF a↓b↓b↑;a↓b↓b↑ (ω) = 1 ω -(ϵ 0 + t + U ) + iη , G HF a↓a↑b↑;a↓a↑b↑ (ω) = 1 ω -(ϵ 0 + 3t + U/2) + iη . (5.8)
The exact multi-channel self-energy calculated from the inverse of the multi-channel Dyson equation (

= [G HF 3 ] -1 -G -1 3 is Σ exact 3 =          0 0 0 U 2 0 0 0 0 0 0 0 -U 2 0 0 0 0 0 0 U 2 0 0 0 0 0 0 0 0 0 -U 0 0 -U 2 0 0 0 0          . 4.1), i.e., Σ exact 3
We notice that the exact multi-channel self-energy is static (as it should be since one cannot have more than three particles in the system, i.e. the added electron and the electron-hole pair which it creates) unlike the exact self-energy for the 1-GF [START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF] which is dynamical.

Approximate results

To construct the approximate multi-channel self-energy in Eqs. where U is the on-site interaction, i = {σ i ĩ} specifies the orbital ĩ with spin σ i , and ĩ is equal to bonding or anti-bonding. With these ingredients the approximate multi-channel self-energy is equal to the exact one, which is given in Eq. (5.9). For completeness, we report the matrix representation of the effective Hamiltonian (4.31). It is given by

H eff =          ϵ 0 -t + U 2 0 0 U 2 0 0 ϵ 0 + t + U 2 0 0 0 -U 0 0 ϵ 0 + t 0 0 U 2 0 0 ϵ 0 + 3t + U 2 0 0 0 0 0 ϵ 0 + t 0 -U 2 0 0 0 ϵ 0 + t + U 2          .
(5.11) Solving [H eff -ω] -1 and taking the head (the upper left 3 × 3 block) of the resulting G 3 , we obtain the exact result for the 1-GF. Notice that, the only non-zero term in the self-energy body, the one with -U , does not affect the result of the 1-GF, but it is important to obtain the exact 3-GF.

In Fig. 5.3 we compare the spectral functions for U/t = 1 obtained with Σ 3 and

The Hubbard dimer

U / t = 4 A(ω) Exact, Σ 3 (ω=0) U / t = 4 A(ω) Σ 1 GW U / t = 4 A(ω) Σ 1 (ω=0) -4 -3 -2 -1 0 1 2 3 4 5 6 7 A 3 (ω) ω Figure 5.4:
The addition part of the 1-and 3-body spectral functions for the Hubbard dimer at 1/4 filling at strong interaction (U/t = 4). Top panel: the spectral function A(ω) obtained with various levels of theory: the exact 1-GF obtained from the multichannel Dyson equation (black solid line); the GW approximation (blue dotted line); the exact static 1-GF (green dashed line). The rightmost peak is a satellite. Bottom panel: the exact 3-body spectral function A 3 (ω). All spectra correspond to ϵ 0 = 1.

Σ 1 (ω = 0), i.e., the exact 1-body self-energy in the static approximation. In the former case we consider only the head of the matrix in Eq. (5.11) to retrieve G 1 . For completeness, we also report the corresponding 3-body spectral function in the bottom panel of Fig. 5.3. At 1/4 filling the spectral function obtained from Σ 3 (ω = 0) is exact. Instead, the spectral function corresponding to Σ 1 (ω = 0) misses the (small) satellite peak, as was expected, and greatly underestimates the position of the highest-energy quasiparticle peak. We also report the spectral function obtained from a dynamical 1-body self-energy, namely the popular GW approximation to the 1-body self-energy.

The analytical result for the GW approximation can be found in Ref. [START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF]. We see that Σ GW 1 yields a very good spectral function at 1/4 filling and weak interaction. Notice the presence of a very small unphysical peak at negative energy caused by the self-screening problem GW suffers from [START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF].

When we increase the interaction strength to U/t = 4 (see Fig. 5.4) we observe that the spectral weight of the satellite in the spectral function has also increased. Instead, in the 3-body spectral function the spectral weight related to the satellite is not influenced by the interaction strength, only its position depends on it. Again, after takeing the head of the matrix in Eq. (5.11) we retrieve G 1 from G 3 which, as mentioned above, leads to the exact spectral function in the case of 1/4 filling. From the spectral function obtained from Σ 1 (ω = 0) we observe that the underestimation of the position of the highest-energy quasiparticle peak is even larger than in the case at weak interaction strength. Moreover, the position of the lowest-lying quasiparticle peak is overestimated. Finally, at strong correlation the energies of the quasiparticles and the satellite in the GW spectral functions are either substantially overestimated or underestimated. Moreover, the extra non-physical peak caused by the self-screening is much more visible.

1/2 filling

We now study the multi-channel 3-GF for the Hubbard dimer at 1/2 filling. The 1-GF channel is a 4 × 4 matrix composed by the four possibilities a ↓, a ↑, b ↓, and b ↑. The

3-GF channel is composed of ijl = (a ↓ a ↑ b ↑), (a ↓ a ↑ b ↓), (b ↓ b ↑ a ↑), (b ↓ b ↑ a ↓).
The Hubbard dimer

Exact results

To construct the exact multi-channel 3-GF we use the eigenvalues and eigenvectors of the symmetric Hubbard dimer (see App. E). We obtain

G 3 =                  a ↓ a ↑ b ↓ b ↑ a ↓ a ↑ b ↑ a ↓ a ↑ b ↓ b ↓ b ↑ a ↑ b ↓ b ↑ a ↓ a ↓
                 , (5.12)
where the diagonal elements are

G a↓;a↓ (ω) = G a↑;a↑ (ω) = = 1 a 2 (1 -A) 2 ω -(ϵ 0 + (c + U )/2 -t) + iη + (1 + A) 2 ω -(ϵ 0 + (U -c)/2 -t) -iη G b↓;b↓ (ω) = G b↑;b↑ (ω) = = 1 a 2 (1 + A) 2 ω -(ϵ 0 + (c + U )/2 + t) + iη + (1 -A) 2 ω -(ϵ 0 + (U -c)/2 + t) -iη G a↓a↑b↑;a↓a↑b↑ (ω) = G a↓a↑b↓;a↓a↑b↓ (ω) = = 1 a 2 (1 -A) 2 ω -(ϵ 0 + (c + U )/2 + t) + iη + (1 + A) 2 ω -(ϵ 0 + (U -c)/2 + t) -iη G b↓b↑a↑;b↓b↑a↑ (ω) = G b↓b↑a↓;b↓b↑a↓ (ω) = = 1 a 2 (1 + A) 2 ω -(ϵ 0 + (c + U )/2 -t) + iη + (1 -A) 2 ω -(ϵ 0 + (U -c)/2 -t) -iη , (5.13)
and the off-diagonal elements are

G a↓;b↓b↑a↑ (ω) = G b↓b↑a↑;a↓ (ω) = = 1 a 2 - (1 -A 2 ) ω -(ϵ 0 + (c + U )/2 -t) + iη + (1 -A 2 ) ω -(ϵ 0 + (U -c)/2 -t) -iη G a↑;b↓b↑a↓ (ω) = G b↓b↑a↓;a↑ (ω) = = 1 a 2 (1 -A 2 ) ω -(ϵ 0 + (c + U )/2 -t) + iη - (1 -A 2 ) ω -(ϵ 0 + (U -c)/2 -t) -iη G b↓;a↓a↑b↑ (ω) = G a↓a↑b↑;b↓ (ω) = = 1 a 2 - (1 -A 2 ) ω -(ϵ 0 + (c + U )/2 + t) + iη + (1 -A 2 ) ω -(ϵ 0 + (U -c)/2 + t) -iη G b↑;a↓a↑b↓ (ω) = G a↓a↑b↓;b↑ (ω) = = 1 a 2 (1 -A 2 ) ω -(ϵ 0 + (c + U )/2 + t) + iη - (1 -A 2 ) ω -(ϵ 0 + (U -c)/2 + t) -iη .
(5.14)

We observe that, in the 1-body channel, two amplitudes go to zero in the noninteracting limit (A = -1 when U = 0). We refer to those contributions as satellite peaks. In the 3-GF channel, satellite peaks have non-zero amplitude in the noninteracting limit, while all the off-diagonal elements go to zero. The non-interacting

, 1 ω -(ϵ 0 -3t + U 2 ) -iη . (5.16)
From the above expressions, it is possible to calculate the exact multi-channel self energy. It is equal to

Σ 3 =               0 0 0 0 0 0 U 2 0 0 0 0 0 0 0 0 -U 2 0 0 0 0 U 2 0 0 0 0 0 0 0 0 -U 2 0 0 0 0 U 2 0 0 0 0 0 0 0 0 -U 2 0 0 0 0 U 2 0 0 0 0 0 0 0 0 -U 2 0 0 0 0 0 0               .
(5.17)

We notice that the exact multi-channel self-energy is static unlike the exact 1-GF self-energy [START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF].

Approximate results

The approximate multi-channel self-energy is obtained using Eqs. (4.10)-(4.11) and (5.10).

As for the 1/4 filling case, it is equal to the exact multi-channel self-energy above. The effective Hamiltonian (Eq. (4.31)) reads

H eff =               ϵ 0 + t + U 2 0 0 0 0 0 U 2 0 0 ϵ 0 + t + U 2 0 0 0 0 0 -U 2 0 0 ϵ 0 -t + U 2 0 U 2 0 0 0 0 0 0 ϵ 0 -t + U 2 0 -U 2 0 0 0 0 U 2 0 ϵ 0 + 3t + U 2 0 0 0 0 0 0 -U 2 0 ϵ 0 + 3t + U 2 0 0 U 2 0 0 0 0 0 ϵ 0 -3t + U 2 0 0 -U 2 0 0 0 0 0 ϵ 0 -3t + U 2               .
(5.18)

The Hubbard dimer Solving [H eff -ω] -1 and taking the head (the upper left 4 × 4 block), we obtain the exact 1-GF. In Figures 5.5 and 5.6 we compare the exact spectral function for the Hubbard dimer at 1/2 filling for U/t = 1 and U/t = 4, respectively, to the spectral functions obtained with the following two approximations: 1) the calculation of the 1-GF using the exact 1-body self-energy in the static approximation (Σ 1 (ω = 0)) and 2) the GW approximation [START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF]. At weak interaction (Figure 5.5) the quasiparticle peaks are very well described by the two approximations considered. On the contrary, the agreement with the exact result for the satellites is not very good. The GW approximation overestimates the energy difference with the nearest quasiparticle energy, while in the spectral function obtained from the static approximation to the 1-body self-energy the satellite amplitudes are not present at all. At strong interaction (Fig. 5.6) the quasiparticle energies are still well described by the static approximations (Σ 1 (ω = 0)), although the gap between the two quasiparticle energies is slightly overestimated. Instead, the GW approximation significantly underestimates this gap. The satellites are absent in the spectral function obtained from the static approximation to the The Hubbard dimer 1-body self-energy, while the GW approximation completely fails to reproduce the positions of the satellites and severely underestimates its amplitudes.

In conclusion, the multi-channel Dyson equation with the approximated static multi-channel self-energy (Eqs. (4.11)-(4.14)) produces the exact results for the symmetric Hubbard dimer at both 1/4 and 1/2 filling unlike the widely used GW approximation.

CHAPTER 6

Conclusions and outlook

In this thesis we addressed the theoretical description of spectroscopy, in particle photoemission and absorption, using Green's function theory. These two spectral probes give information about the electronic structure, which is the cornerstone of our understanding of material properties. From the theory point of view, many-body perturbation theory (MBPT) based on Green's functions is the method of choice to obtain photoemission spectra in an ab-initio way, i.e. without empirical or adjustable parameters. Indeed the photoemission spectrum (PES) can be obtained in a straightforward way through the one-body Green's function (1-GF). Within MBPT the GW approach (where G is the 1-GF and W the screened Coulomb interaction) has been by far the method of choice for most applications. Instead, to describe optical spectra within MBPT one can use the two-body Green's function. Each Green's function can be obtained from a Dyson equation which links the non-interacting Green's function to the interacting one through a kernel that contains the many-body effects. For example, the non-interacting 1-GF is linked to the interacting one through the self-energy. Similar Dyson-like equations exist for higher-order propagators, e.g., the Bethe-Salpeter equation (BSE), which is related to the two-body Green's function. Each of these Dyson equations involves just a single type of many-body Green's func-93

Conclusion

tion.

In this thesis, we extended the concept of the Dyson equation to a multi-channel version which allows us to couple two or more Green's functions. Thus, for a given physical process, the pertinent physics of each of the Green's functions involved can be exploited. In this thesis, we focused on the coupling of the one-body Green's function and the three-body Green's function to calculate photoemission spectra, and on the coupling of the two-body Green's function and the four-body Green's function to calculate optical spectra, but our method can be easily generalized to the coupling of other many-body Green's functions. We introduced an approximate self-energy based on the well-established "Random Phase approximation +exchange" (RPAx) approximation to the particle-particle and electron-hole BSE kernel. By reformulating the multi-channel Dyson equation in terms of an effective Hamiltonian, along the same line as for the BSE, the solution of the original equation reduces to an eigenvalue problem. Finally, we tested the validity of our approach for the 1and 3-GF channel using the symmetric Hubbard dimer. Despite the simplicity of this model, it is challenging for standard approximations to the self-energy, such as GW and T-matrix, to capture the correct physics of the dimer at both 1/2 and 1/4 filling. In this respect, our approach is very promising since it gives the exact result at both fillings. Our next step is to apply the multi-channel Dyson equation to real systems, starting from simple molecules to simple semiconductors (e.g. silicon) to strong correlated systems, such as transition metal oxides, which are a challenge for ab initio theories. As far as the part on the 2-4-GF coupling is concerned, the next natural step is to test it on the Hubbard dimer. This system at 1/2 filling has a neutral transition with a double excitation character which cannot be described by the BSE with a static kernel (which is the standard protocol used in real systems). Therefore, it will be a useful test for our theory. Successively, we can move on to real systems. Although higher-order excitations are usually dark in the experiments, they can play an important role for the description of ground and optically active excites states. We will apply our approach to closed-shell molecules such as polyenes, for example, where the lowest-lying singlet state is known to have a HOMO 2 -LUMO 2 double-excitation character.

CHAPTER 7

Résumé en français

Le monde physique est composé de particules interagissantes. En physique de la matière condensée, les électrons déterminent de nombreuses propriétés matérielles telles que la conductivité électronique et thermique, la magnétisation, etc. Par conséquent, élaborer une théorie ab initio qui explique et prédit également le comportement électronique des matériaux est très important. En général, il est impossible de résoudre l'équation de Schrödinger en raison de la croissance exponentielle des dimensions de la fonction d'onde à plusieurs corps avec le nombre croissant d'électrons [START_REF] Montina | Exponential complexity and ontological theories of quantum mechanics[END_REF]. Des approximations sont donc nécessaires. Certaines des plus réussies sont la théorie de la fonctionnelle de la densité (DFT) et la théorie des perturbations à plusieurs corps (MBPT). Ces théories évitent l'utilisation de la fonction d'onde et utilisent des quantités plus simples, telles que la densité électronique et la fonction de Green à un corps (1-GF). En utilisant un cadre simplifié, il est nécessaire d'introduire des potentiels effectifs contenant la complexité du problème initial.

En DFT [START_REF] Hohenberg | Inhomogeneous electron gas[END_REF][START_REF] Kohn | Self-consistent equations including exchange and correlation effects[END_REF], on remplace le système interagissant par un système fictif de particules non interagissant, le système de Kohn-Sham (KS), dans lequel, sous l'influence d'un potentiel effectif (le potentiel de Kohn-Sham), la densité électronique du système interagissant est reproduite. On doit alors résoudre un ensemble d'équations à une 95 résumé en français particule. Ces équations doivent être résolues de manière auto-cohérente car le potentiel de KS dépend de la densité. Cependant, les approximations au potentiel de Kohn-Sham ne sont pas intuitives, et, de plus, il n'est pas toujours connu comment extraire l'observable d'intérêt de la densité.

La MBPT basée sur les fonctions de Green [START_REF] Bechstedt | Many-body approach to electronic excitations[END_REF][START_REF] Nakanishi | A general survey of the theory of the bethe-salpeter equation[END_REF][START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF] est la méthode d'intérêt dans cette thèse. Les approches standard sont basées sur la 1-GF ou la 2-GF. La 1-GF est très utile pour calculer les spectres de photoémission car le spectre de photoémission peut être obtenu de manière directe à partir de la partie imaginaire de la 1-GF. De plus, elle donne accès à de nombreuses autres propriétés, telles que les surfaces d'énergie potentielle, les densités électroniques, les distributions de moment, etc. Dans l'approximation dite GW [START_REF] Hedin | New method for calculating the one-particle green's function with application to the electron-gas problem[END_REF], la MBPT est devenue la méthode de choix pour le calcul des structures de bande de quasi-particules [1,2,[START_REF] Ia Nechaev | Evidence for a direct band gap in the topological insulator bi 2 se 3 from theory and experiment[END_REF][START_REF] Vidal | Strong interplay between structure and electronic properties in CuIn(S, Se) 2 : A first-principles study[END_REF][START_REF] Vidal | Falsepositive and false-negative assignments of topological insulators in density functional theory and hybrids[END_REF][START_REF] Waroquiers | Band widths and gaps from the tran-blaha functional: Comparison with many-body perturbation theory[END_REF] et des spectres de photo-émission directe et inverse [START_REF] Athanasios | Quasiparticle self-consistent gw method applied to localized 4f electron systems[END_REF][START_REF] Sergey | All-electron selfconsistent gw approximation: Application to si, mno, and nio[END_REF][START_REF] Gatti | Understanding correlations in vanadium dioxide from first principles[END_REF][START_REF] Lim | Angle-resolved photoemission and quasiparticle calculation of zno: The need for d band shift in oxide semiconductors[END_REF][START_REF] Papalazarou | Valence-band electronic structure of v 2 o 3 : Identification of v and o bands[END_REF] de nombreux matériaux, améliorant considérablement les résultats fournis par les méthodes statiques de structure électronique en champ moyen. Cependant, GW présente certaines lacunes fondamentales telles que l'auto-écrantage [START_REF] Nelson | Self-interaction in green'sfunction theory of the hydrogen atom[END_REF][START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF], le manque de complète auto-cohérence (dans les systèmes étendus) [START_REF] Caruso | Unified description of ground and excited states of finite systems: The self-consistent gw approach[END_REF][START_REF] Stan | Levels of selfconsistency in the GW approximation[END_REF][START_REF] Van Schilfgaarde | Quasiparticle self-consistent gw theory[END_REF], le manque de taille-consistance [START_REF] Erik Dahlen | Variational energy functionals of the green function and of the density tested on molecules[END_REF] et, en particulier, il ne décrit pas bien la forte corrélation (sans imposer un ordre magnétique) [START_REF] Romaniello | Double excitations in finite systems[END_REF]. Par exemple, l'oxyde de nickel paramagnétique, qui est un isolant, est décrit comme un métal par GW . C'est parce que GW , à part l'échange, donne une description classique du système et de sa réponse à un électron ou un trou supplémentaire [START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF]. De plus, il ne fonctionne pas bien pour les caractéristiques spectrales qui sont liées à d'autres excitations du système, telles que les plasmons [START_REF] Guzzo | Valence electron photoemission spectrum of semiconductors: Ab initio description of multiple satellites[END_REF]. La 2-GF est quant à elle liée aux excitations neutres. Ces énergies d'excitation sont obtenues en résolvant l'équation de Bethe-Salpeter (BSE) [START_REF] Martin | Interacting Electrons: Theory and Computational Approaches[END_REF][START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF]. La BSE est une équation de type Dyson et repose sur la 2-GF, spécifiquement une partie appelée la fonction de corrélation à deux corps. La description d'une paire électron-trou est bien décrite par l'approximation la plus courante et largement utilisée du noyau de la BSE, qui est basée sur l'approximation GW . Cette approximation du noyau de la BSE est dynamique, c'est-à-dire qu'elle dépend de la fréquence. Cependant, dans les calculs pratiques, elle est souvent considérée comme statique, c'est-à-dire qu'elle devient indépendante de la fréquence. Par conséquent, ce traitement exclut la possibilité de décrire les excitations doubles, c'est-à-dire la formation simultanée de deux paires électron-trou [START_REF] Romaniello | Double excitations in finite systems[END_REF][START_REF] Sangalli | Double excitations in correlated systems: A many-body approach[END_REF], et en général les excitations d'ordre supérieur.

Dans cette thèse, nous nous concentrerons sur la description théorique de la photoémission et de la spectroscopie optique. Ces deux techniques sont des outils précieux pour obtenir des informations sur la structure électronique, qui est la pierre angulaire de notre compréhension des propriétés des matériaux. french 97

La spectroscopie de photoémission

La spectroscopie de photoémission directe et inverse sont les principales techniques pour étudier la structure électronique des solides. Dans la spectroscopie de photoémission directe (PES) [START_REF] Hüfner | Photoelectron Spectroscopy: Principles and Applications[END_REF][START_REF] Nordling | Precision method for obtaining absolute values of atomic binding energies[END_REF], des photons de haute énergie, généralement des rayons X, sont dirigés sur la surface du matériau. Ces photons ont une énergie suffisante pour ioniser les électrons des atomes dans le matériau. Lorsqu'un électron est éjecté, il crée un trou, c'est-à-dire un niveau d'énergie non occupé dans la structure électronique. En mesurant l'énergie cinétique des électrons émis, des informations sur la structure électronique sont obtenues, telles que la structure des bandes de valence. La spectroscopie de photoémission inverse (IPES) [START_REF] Binnig | Tunneling spectroscopy and inverse photoemission: Image and field states[END_REF][START_REF] Cardona | Photoemission in Solids[END_REF] exploite le processus opposé. Ici, un électron est absorbé dans le matériau, provoquant l'émission d'un photon. En mesurant l'énergie du photon, des informations sur la structure des bandes de conduction sont obtenues.

Pour les matériaux qui ne sont pas fortement corrélés, il est possible de distinguer deux types de pics dans le spectre : les quasiparticules et les satellites. Les pics de quasiparticules sont dus aux excitations habillées à une particule, c'est-à-dire qu'il est possible de les décrire comme des particules non interagissantes avec une durée de vie finie et une masse effective qui contient l'effet de l'interaction avec toutes les particules du système, telles que les noyaux et les autres électrons. Les pics de satellites sont dus à l'interaction résiduelle et entraînent des excitations supplémentaires. Ils ne peuvent pas être traités dans un modèle de particules indépendantes. Ces excitations supplémentaires des pics de satellites sont généralement dues au couplage des électrons avec des paires d'électrons ou de (quasi)bosons, tels que les plasmons (voir par exemple dans le sodium [START_REF] Steiner | Simple metals[END_REF]), les phonons et les magnons. Dans le cas d'une forte corrélation, cette image s'effondre car les quasiparticules et les satellites se mélangent. Ce régime est particulièrement intéressant car des propriétés électroniques non triviales ou "exotiques" des matériaux émergent, représentant ainsi un défi pour les méthodes théoriques. La figure 7.1 montre une représentation schématique des deux techniques.

Absorption

L'absorption se produit lorsqu'un système subit une transition d'un état d'énergie inférieure à un état d'énergie supérieure en absorbant de l'énergie, normalement transportée par des photons. Ces photons n'ont pas assez d'énergie pour ioniser les électrons, mais ils sont suffisamment énergétiques pour exciter un électron d'une bande de valence à une bande de conduction. Ce processus correspond à la création d'une paire électron-trou. La figure 7.2 représente de manière schématique le processus d'absorption. D'autres excitations, ou satellites, peuvent également apparaître dans les spectres. Par exemple, il est bien connu que la forte localisation spatiale peut conduire à la formation simultanée de deux paires électron-trou [START_REF] Sangalli | Double excitations in correlated systems: A many-body approach[END_REF].

Pour une analyse complète de la spectroscopie de photoémission et de l'absorption, de nombreux effets différents doivent être pris en compte. Par exemple, le photon entrant (électron) subit des effets de surface avant de pénétrer dans le matériau; dans le matériau, tous les électrons interagissent les uns avec les autres et avec les noyaux, générant des contributions électron-phonon. De plus, les effets thermiques, les défauts et les impuretés modifieront les spectres finaux. Dans cette thèse, nous concentrerons notre attention sur l'interaction électron-électron. Par conséquent, nous considérerons des systèmes à température nulle sans défauts ni impuretés. De plus, nous limiterons l'analyse à l'approximation de Born-Oppenheimer, où les électrons interagissent avec les noyaux uniquement par l'intermédiaire d'un potentiel de champ moyen. Les interactions électron-phonon ne sont pas prises en compte.

La fonction de Green à trois corps

Pour étudier la spectroscopie de photoémission directe et inverse, la fonction de Green à un corps est la quantité la plus naturelle, car sa fonction spectrale est directement liée à l'elevement et à l'ajout d'électrons. Cependant, les approximations nécessaires pour le calculer donnent des spectres qui ne sont pas toujours précis. En particulier, la description de la structure des satellites est très difficile. Dans ce chapitre, nous montrons qu'il est possible de calculer les spectres en partant également de la fonction de Green à trois corps. Nous montrerons qu'au niveau de la fonction de Green à trois corps non interagissante, il existe déjà des informations liées aux satellites. Par conséquent, comme nous le montrerons, une self-énergie statique (selfénergie à trois corps) est suffisante pour obtenir à la fois les quasi-particules et les satellites dans les spectres de photoémission. Nous notons que c'est une stratégie générale : plus d'informations contenues dans la quantité fondamentale signifient moins d'informations requises dans le potentiel effectif, c'est-à-dire l'self-énergie dans notre cas, pour décrire les effets many-body pertinents. La fonction de Green à trois corps est définie par

G 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = i⟨Ψ N 0 |T [ ψ(1) ψ(2) ψ(3) ψ † (3 ′ ) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N 0 ⟩. (7.1)
La fonction de Green à trois corps dépend de six différences de temps ou cinq différences de temps lorsque l'hamiltonien est indépendant du temps, et le nombre total de permutations des opérateurs de champ dans l'équation (7.1) en raison de l'opérateur T est de 6! = 720. En fonction de l'ordre des opérateurs de champ (et donc des temps), la fonction de Green à trois corps fournit des informations différentes.

résumé en français

En général, elle décrit la propagation de trois particules (électrons ou trous), et la fonction de Green à trois corps peut donc être divisée en quatre composantes :

G eee 3 , G hhh 3 , G eeh 3 et G hhe 3 .
Afin de rendre cette séparation explicite, on peut réécrire les six opérateurs de champ ordonnés dans le temps dans l'équation (7.1) comme une somme de produits de deux termes, chacun contenant trois opérateurs de champ ordonnés dans le temps Tous les termes dans autres permutations sont non singuliers aux fréquences égales aux énergies de enlevement ou d'ajout d'électrons. Par conséquent, leur contribution aux spectres finaux est nulle. À partir de l'équation (3.22), il est clair que les pôles des parties eeh et hhe de la 3-GF sont les mêmes que les pôles de la 1-GF, voir l'équation (2.53). En analogie avec la différence de temps de la 1-GF, il est possible d'interpréter τ comme le temps de propagation combinée de la partic-french 101 ule ajoutée (électron ou trou) et de la paire électron-trou. Les quatre différences de temps restantes correspondent aux processus physiques suivants : 1) le temps entre la particule ajoutée et la création de la paire électron-trou ; 2) le temps nécessaire pour créer la paire électron-trou ; 3) le temps nécessaire pour recombiner la paire électron-trou ; 4) le temps entre la recombinaison et le enlevement de la particule. La correspondance entre la différence de temps et chaque processus dépend de l'ordre des temps. Pour la description de la spectroscopie (inverse) photoélectronique, les quatre processus peuvent être considérés comme instantanés. Par conséquent, nous pouvons prendre la limite τ ij → 0 pour chacune des quatre différences de temps. Il est pratique de choisir un ordre temporel qui est cohérent avec la chronologie du processus de (inverse) photoémission. Par exemple, dans la spectroscopie de photoémission directe, d'abord un électron est émis du système conduisant à la création de paires électron-trou. Après un temps τ , les paires électron-trou recombinent et finalement un électron est ajouté. Cela correspond à l'ordre suivant des opérateurs de champ ψ † ψ † ψ ψ † ψ ψ agissant sur |Ψ N 0 ⟩. Cet ordre des opérateurs de champ peut être obtenu avec le choix suivant pour les différences de temps,

G 3 (1, 2, 3, 1 ′ , 2 ′ , 3 ′ ) = = i n ⟨Ψ N 0 |T [ ψ(1) ψ(2) ψ(3)]|Ψ N +3 n ⟩⟨Ψ N +3 n |T [ ψ † (3 ′ ) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N 0 ⟩ -i n ⟨Ψ N 0 |T [ ψ † (3 ′ ) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N -3 n ⟩⟨Ψ N -3 n |T [ ψ(1) ψ(2) ψ(3)]|Ψ N 0 ⟩ -i n ⟨Ψ N 0 |T [ ψ(1) ψ(2) ψ † (3 ′ )]|Ψ N +1 n ⟩⟨Ψ N +1 n |T [ ψ(3) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N 0 ⟩ + i n ⟨Ψ N 0 |T [ ψ(3) ψ † (2 ′ ) ψ † (1 ′ )]|Ψ N -1 n ⟩⟨Ψ N -1 n |T [ ψ(1) ψ(2) ψ † (3 ′ )]|Ψ N 0 ⟩ + autres permutations. ( 7
τ 12 = 0 -, τ 23 ′ = 0 -, τ 1 ′ 2 ′ = 0 + , τ 2 ′ 3 = 0 + . (7.3)
Nous remarquons qu'en raison de la présence des fonctions de Heaviside dans les équations (3.9), (3.17

G e+h 3 (x 1 , x 2 , x 3 , x 1 ′ , x 2 ′ , x 3 ′ ; ω) = = n X n (x 1 , x 2 , x 3 ′ )X * n (x 1 ′ , x 2 ′ , x 3 ) ω -(E N +1 n -E N 0 ) + iη + n Z * n (x 1 ′ , x 2 ′ , x 3 )Z n (x 1 , x 2 , x 3 ′ ) ω -(E N 0 -E N -1 n ) -iη , (7.4)
où les amplitudes électron-électron-trou et trou-trou-électron, X n et Z n respectivement, sont définies comme

X n (x 1 , x 2 , x 3 ′ ) = ⟨Ψ N 0 | ψ † (x 3 ′ ) ψ(x 2 ) ψ(x 1 )|Ψ N +1 n ⟩, (7.5)
Z n (x 1 , x 2 , x 3 ′ ) = ⟨Ψ N -1 n | ψ † (x 3 ′ ) ψ(x 2 ) ψ(x 1 )|Ψ N 0 ⟩. (7.6)
Pour plus de complétude, nous donnons également ici les expressions explicites des conjugués complexes de ces amplitudes,

X * n (x 1 ′ , x 2 ′ , x 3 ) = ⟨Ψ N +1 n | ψ † (x 1 ′ ) ψ † (x 2 ′ ) ψ(x 3 )|Ψ N 0 ⟩, (7.7)
Z * n (x 1 ′ , x 2 ′ , x 3 ) = ⟨Ψ N 0 | ψ † (x 1 ′ ) ψ † (x 2 ′ ) ψ(x 3 )|Ψ N -1 n ⟩. (7.8)
résumé en français La fonction de Green à trois corps non interagissante peut être écrite dans la représentation matricielle suivante

G 0 3 (ω) = G 0 1 (ω) 0 0 G 0(3p) (ω) , (7.9
) avec G 0 1 (ω) définie comme G 0 ij (ω) = δ ij ω -ϵ 0 i + sign(ϵ 0 i -µ)iη . (7.10) et G 0(3p) (ω) définie comme G 0(3p) i>jl;m>ok (ω) = δ im δ jo δ lk (f i -f l )(f j -f l ) ω -ϵ 0 i -(ϵ 0 j -ϵ 0 l ) + iηsign(ϵ 0 i -µ) . (7.11)
Pour calculer la fonction de Green à trois corps, introduisons l'équation de Dyson multicanal (7.12) La self-énergie multi-canal est définie comme

Σ 3 = Σ 1p Σ c Σc Σ 3p . (7.13)
Nous aboutissons à l'approximation statique suivante, Remarquez que, en raison de la condition i > j et m > o dans l'équation (7.11) sur chaque ligne G 0 3 , la somme de quatre diagrammes multi-canaux différents, chacun composé de termes d'self-énergie multi-canal différents, est nécessaire pour obtenir un diagramme du troisième ordre. Ici, pour simplifier, seulement l'un d'entre eux est montré.

Σ 3p ijl;mok = [(1-f i )(1-f j )f l -f i f j (1-f l )][δ lk vijom +δ mj viklo +δ io vjklm -δ oj viklm -δ im vjklo ], ( 7
Σ 2e1h (x 1 , x 2 , x 3 ′ , x 1 ′ , x 2 ′ , x 3 ) = -Σ 2h1e (x 1 , x 2 , x 3 ′ , x 1 ′ , x 2 ′ , x 3 ) = δ(x 3 , x 3 ′ )[δ(x 1 , x 1 ′ )δ(x 2 , x 2 ′ ) -δ(x 1 , x 2 ′ )δ(x 2 , x 1 ′ )]v(r 1 , r 2 ) + δ(x 1 ′ , x 2 )[δ(x 1 , x 2 ′ )δ(x 3 , x 3 ′ ) -δ(x 1 , x 3 ′ )δ(x 3 , x 2 ′ )]v(r 1 , r 3 ) + δ(x 1 , x 2 ′ )[δ(x 2 , x 1 ′ )δ(x 3 , x 3 ′ ) -δ(x 2 , x 3 ′ )δ(x 3 , x 1 ′ )]v(r 2 , r 3 ) -δ(x 2 , x 2 ′ )[δ(x 1 , x 1 ′ )δ(x 3 , x 3 ′ ) -δ(x 1 , x 3 ′ )δ(x 3 , x 1 ′ )]v(r 1 , r 3 ) -δ(x 1 , x 1 ′ )[δ(x 2 , x 2 ′ )δ(x 3 , x 3 ′ ) -δ(x 2 , x 3 ′ )δ(x 3 , x 2 ′ )]v(r 2 , r 3 ) (7.18) Σ c (x 1 , x 1 ′ , x 2 ′ , x 3 ) = [δ(x 1 , x 1 ′ )δ(x 2 ′ , x 3 ) -δ(x 1 , x 2 ′ )δ(x 1 ′ , x 3 )]v(r 1 , r 3 ) (7.19) Σc (x 1 , x 2 , x 3 ′ , x 1 ′ ) = [δ(x 1 , x 1 ′ )δ(x 2 , x 3 ′ ) -δ(x 1 , x 3 ′ )δ(x 2 , x 1 ′ )]v(r 1 , r
Le diagramme buble dans l'équation (7.26) est également contenu dans l'approximation GW qui inclut des diagrammes bubles jusqu'à un ordre infini. Dans notre approche, les diagrammes d'ordre supérieur sont obtenus en itérant davantage l'équation de Dyson multi-canal (7.12). En le faisant, il est possible de vérifier que tous les diagrammes GW sont inclus dans notre approche, et les effets de écrantage sont ainsi pris en compte. De plus, notre approche va au-delà du écrantage RPA inclus dans GW . À titre d'exemple, nous rapportons le diagramme du quatrième ordre suivant (7.28) où une correction de premier ordre au vertex de la buble RPA est montrée. Bien que l'self-énergie n'implique que le potentiel de Coulomb nu, le écrantage au-delà de la RPA est naturellement inclus dans l'équation de Dyson multi-canal. Par conséquent, il n'est pas nécessaire d'utiliser une interaction écrantage dans notre théorie, car cela conduirait au double comptage des diagrammes.

Enfin, pour comprendre le nombre de diagrammes proper skeleton produits à chaque ordre par notre théorie, il faut tenir compte du nombre de termes Σ 3 , c'est-àdire deux termes pour chaque terme de couplage et dix termes pour Σ G 1 . Par conséquent, le nombre total de diagrammes proper skeleton inclus dans l'équation (7.12), dans l'self-énergie multi-canal approximative des équations (7.14)-(7.17), à l'ordre n est donné par 2×2×10 n-2 2 n-1 = 10 n-2

2 n-3 pour n ≥ 2. L'itération de l'équation de Dyson multi-canal (7.12) crée également des termes tels que Σ c G 0(3p) 3

Σc G 0 1 Σ c G 0(3p) 3
Σc (voir Fig. 7.4). Par conséquent, chaque fois qu'un G 0 1 entre entre deux termes d'self-énergie proper, des répétitions de termes d'selfénergie proper apparaissent, exactement comme dans l'équation de Dyson (2.56). Les diagrammes non-skeleton sont présents en raison de G 0(3p) dans l'équation (3.49). Pour les améliorer, on peut habiller G 0(3p) au-delà de HF en utilisant, par exemple, le second Born, GW ou la matrice T dans l'approximation des quasi-particules. Dans ce cas, on devrait supposer que le G 0(3p) corrélé est diagonal dans la même base que G 0 3 . Cette approximation est également utilisée, par exemple, pour dériver l'approximation des cumulants [START_REF] Guzzo | Valence electron photoemission spectrum of semiconductors: Ab initio description of multiple satellites[END_REF][START_REF] Martin | Interacting Electrons: Theory and Computational Approaches[END_REF].

Pour obtenir une équation qui peut être résolue en pratique à l'aide d'outils numériques standard, nous utilisons une stratégie similaire à celle utilisée pour l'équation de Bethe-Salpeter [START_REF] Onida | Electronic excitations: density-functional versus many-body green's-function approaches[END_REF], c'est-à-dire, nous transformons l'équation (7.12) en un hamiltonien effectif à trois particules selon G 3(ijl;mok) = ωI -H eff -1 ijl;mok (7.29) dans lequel l'hamiltonien effectif H eff est donné par

H eff = H 1p Σ c Σc H 3p , (7.30)
où

H 1p i;m = ϵ 0 i δ im , (7.31)
H 3p ijl;mok = (ϵ 0 i -(ϵ 0 l -ϵ 0 j ))δ im δ jo δ lk + Σ 3p ijl;mok . (7.32)

résumé en français

Puisque notre objectif final est de calculer les spectres de photoémission, qui sont liés au canal à une particule de G 3 , il suffit d'extraire la tête de G 3 dans l'équation (7.29).

Une diagonalisation directe de l'équation (7.30) aurait une complexité en temps de l'ordre de N 9 , où N est le nombre d'électrons dans le système. Cependant, cette complexité peut être significativement réduite en utilisant des méthodes itératives standard telles que le solveur Haydock-Lanczos [START_REF] Haydock | Electronic structure based on the local atomic environment for tight-binding bands[END_REF][START_REF] Hernandez | Slepc: A scalable and flexible toolkit for the solution of eigenvalue problems[END_REF][START_REF] Schmidt | Efficient O(N 2 ) method to solve the bethe-salpeter equation[END_REF]. Avec cette méthode, la complexité est réduite à N 6 car elle est alors déterminée par la construction de l'hamiltonien effectif. Nous notons également que les méthodes qui décrivent les trions ont une complexité similaire à notre approche, et elles ont été appliquées avec succès à des systèmes réels [START_REF] Deilmann | Three-particle correlation from a many-body perspective: Trions in a carbon nanotube[END_REF][START_REF] Torche | First-principles many-body theory for charged and neutral excitations: Trion fine structure splitting in transition metal dichalcogenides[END_REF]

.9: La fonction spectrale du dimère de Hubbard avec un remplissage de 1/2 et une forte interaction (U/t = 4) obtenue avec différents niveaux de théorie. Résultat exact obtenu à partir de l'équation de Dyson multi-canal (ligne pleine noire) ; l'approximation GW (ligne pointillée bleue) ; la 1-GF statique exacte (ligne en pointillés verts). Les pics externes sont les satellites. Les spectres correspondent à ϵ 0 = -U/2, ce qui garantit la symétrie particule-trou. french 113

Conclusion

Dans cette thèse, nous avons abordé la description théorique de la spectroscopie, en particulier la photoémission de particules et l'absorption, en utilisant la théorie des fonctions de Green. Ces deux sondes spectrales fournissent des informations sur la structure électronique, qui est la pierre angulaire de notre compréhension des propriétés des matériaux.

Du point de vue théorique, la théorie des perturbations à plusieurs corps (MBPT) basée sur les fonctions de Green est la méthode privilégiée pour obtenir des spectres de photoémission de manière ab initio, c'est-à-dire sans paramètres empiriques ou ajustables. En effet, le spectre de photoémission (PES) peut être obtenu de manière directe grâce à la fonction de Green à une particule (1-GF). Dans le cadre de la MBPT, l'approche GW (où G est la 1-GF et W l'interaction de Coulomb écrantée) a largement été la méthode de choix pour la plupart des applications. En revanche, pour décrire les spectres optiques dans le cadre de la MBPT, on peut utiliser la fonction de Green à deux corps.

Chaque fonction de Green peut être obtenue à partir d'une équation de Dyson qui relie la fonction de Green non interactive à celle interactive à travers un noyau qui contient les effets à plusieurs corps. Par exemple, la 1-GF non interactive est liée à la 1-GF interactive par l'intermédiaire de l'self-énergie. Des équations de Dyson similaires existent pour des propageurs d'ordre supérieur, par exemple l'équation de Bethe-Salpeter (BSE), qui est liée à la fonction de Green à deux corps. Chacune de ces équations de Dyson implique un seul type de fonction de Green à plusieurs corps.

Dans cette thèse, nous avons étendu le concept de l'équation de Dyson à une version multi-canal qui nous permet de coupler deux ou plusieurs fonctions de Green. Ainsi, pour un processus physique donné, la physique pertinente de chacune des fonctions de Green impliquées peut être exploitée. Dans cette thèse, nous nous sommes concentrés sur le couplage de la fonction de Green à une particule et de la fonction de Green à trois corps pour calculer les spectres de photoémission, et sur le couplage de la fonction de Green à deux corps et de la fonction de Green à quatre corps pour calculer les spectres optiques, mais notre méthode peut être facilement généralisée au couplage d'autres fonctions de Green à plusieurs corps. Nous avons introduit une self-énergie approximative basée sur l'approximation bien établie "Random Phase approximation + exchange" (RPAx) pour le noyau de l'équation de Bethe-Salpeter des trous électroniques et des particules. En reformulant l'équation de Dyson multi-canal en termes d'un hamiltonien effectif, sur la même ligne que pour la BSE, la solution de l'équation originale se réduit à un problème de valeurs propres.

Enfin, nous avons testé la validité de notre approche pour les canaux 1 et 3-GF résumé en français en utilisant le dimère de Hubbard symétrique. Malgré la simplicité de ce modèle, il est difficile pour les approximations standard de l'self-énergie, telles que GW et T-matrix, de capturer correctement la physique du dimère aux remplissages de 1/2 et 1/4. À cet égard, notre approche est très prometteuse car elle donne le résultat exact aux deux remplissages. Notre prochaine étape est d'appliquer l'équation de Dyson multi-canal à des systèmes réels, en commençant par des molécules simples jusqu'à des semi-conducteurs simples (par exemple, le silicium) et des systèmes fortement corrélés, tels que les oxydes de métaux de transition, qui représentent un défi pour les théories ab initio.

En ce qui concerne la partie sur le couplage 2-4-GF, la prochaine étape naturelle est de le tester sur le dimère de Hubbard. Ce système à un remplissage de 1/2 a une transition neutre avec un caractère d'excitation double qui ne peut pas être décrit par la BSE avec un noyau.

2ϵ 0 + U 0 0 -1 √ 2 1 √ 2 0 0 2ϵ 0 0 0 0 0 -1 0 2ϵ 0 0 -1 0 0 0 0 2ϵ 0 0 0 -1 √ 2 -1 √ 2 0 0 Table E.1:
Eigenvalues and coefficients of the expansion of the eigenstates in the bounding/antibonding basis for the two electrons system (1/2 filling).

3/4 filling

The eigenstate of the Hamiltonian at 1/4 filling are |b ↑ b ↓ a ↑⟩, |b ↑ b ↓ a ↓⟩, |b ↑ a ↑ a ↓⟩, |b ↓ a ↑ a ↓⟩ with energies 3ϵ 0 + U -t, 3ϵ 0 + U -t, 3ϵ 0 + U + t, 3ϵ 0 + U + t, respectively.

Figure 1 . 1 :

 11 Figure 1.1: A schematic representation of direct (PES) and inverse (IPES) photoemission spectroscopy. The removal of an electron (hole) from the system leaves a photohole (photoelectron), which excites the system and leads to the creation of electron-hole pairs.

Figure 1 . 2 :

 12 Figure 1.2: A schematic representation of the creation of an electron-hole pair due to the absorption of a photon.

Figure 2 . 1 :

 21 Figure 2.1: Singularities of G(ω) in the complex plane. The removal energies E N 0 -E N -1 m are located above the real axis, whereas the addition energies E N +1 m -E N 0 are located below. The chemical potential µ lies inside the band gap ∆E g .

20 BackgroundFigure 2 . 2 :

 2022 Figure 2.2: Illustrative process of the Hedin pentagon.On the left there is the selfconsistent procedure that is difficult to do in practical calculation. On the right there is the Hedin pentagon in the GW approximation. In this later case, the equations to solve reduce to four, and accordingly the pentagon reduces to a quadrangle.

Background

21

 21

22 BackgroundFigure 2 . 3 :

 2223 Figure 2.3: Diagrammatic representation of the three combinations described in Eq. (2.75). Each arrow line represents a non-interacting 1-GF, and the wiggly lines represent the Coulomb potential v. The first term from the left represents the firstorder Hartree term, the second is the first-order Fock contribution and the last one corresponds to a vacuum diagram.

Figure 2 . 4 :

 24 Figure 2.4: Three examples of second-order self-energy diagrams. The first from the left is an improper diagram: when we remove the G 0 line in the middle, it disconnects the diagram into two separate first-order exchange self-energy contributions. The diagram in the middle is an example of a non-skeleton diagram since it has an exchange self-energy insertion. The last diagram is both proper and skeleton, representing the self-energy at first-order in the RP A irreducible polarizability (see Eq. (2.68)).

Figure 2 .

 2 Figure 2.4 shows some examples of self-energy diagrams.Having a proper self-energy is crucial when solving the Dyson equation. Using an improper self-energy leads to a double-counting problem, where certain diagrams enter the 1-GF more times than they should. This issue becomes evident when examining the Dyson equation order by order, namelyG = G 0 + G 0 ΣG 0 + G 0 ΣG 0 ΣG 0 + ....The concept of a skeleton self-energy is significant in simplifying the understanding of diagrams, as it allows us to focus on the essential contributions while neglecting some higher-order effects. To illustrate this point, let us consider the GW self-energy Σ = iGW . It relies on the 1-GF both directly and through the irreducible polarizability present in W , as given in equations (2.67) and (2.68). When considering this 1-GF as non-interacting, only skeleton diagrams are present. Non-skeleton diagrams

Figure 2 . 5 :Figure 2 . 6 :

 2526 Figure 2.5: Diagrammatic representation of the Hedin equations. (a) The Dyson equation (2.65), where the double line represents the full interacting 1-GF. (b) The equation for the self-energy (Eq. (2.66)): the first term corresponds to the Hartree contribution and the second term represents the exchange-correlation part. (c) The Dyson-like equation for the screened Coulomb interaction (Eq. (2.67)), where the double wiggly line represents the screened Coulomb interaction W . (d) The equation for the irreducible polarizability (Eq. (2.68)). (e) The equation for the irreducible vertex (2.69).

Figure 2 . 7 :

 27 Figure 2.7: First and second order contribution in the irreducible polarizability χ of the screened Coulomb interaction in the RP A approximation.

Figure 2 . 8 :

 28 Figure 2.8: Comparison between the two first order topologically distinct connected diagrams of the pp channel. The diagram on the left has zero interchanges between lines, i.e., l = 0. The diagram on the right has one interchange, i.e., l = 1.

Figure 2 . 9 :

 29 Figure 2.9: Two examples of diagrams that are at second order in the interaction. On the left there is a diagram where the two interactions cross each other. On the right there is a first-order screening diagram in the interaction between the two particle lines.

Figure 2 . 10 :

 210 Figure 2.10: Comparison between the two first order topologically distinct connected diagrams of the eh channel. The diagram on the left has zero interchanges between lines, i.e., l = 0. The diagram on the right represents the interchange of two lines in the eh channel, i.e., l = 1.

Figure 2 . 11 :

 211 Figure 2.11: The two second-order diagrams that are created by repeating the same first-order term.

Figure 2 . 13 :

 213 Figure 2.13: Diagrammatic structure of the kernel of the BSE in the HF approximation. The dotted line represents points collapsing in a Dirac delta and the wiggly line represent the bare Coulomb interaction.

Figure 2 .

 2 [START_REF] Erik Dahlen | Variational energy functionals of the green function and of the density tested on molecules[END_REF] shows diagrams not present in the BSE with a GW kernel.

Figure 2 . 14 :

 214 Figure 2.14: The BSE with only the exchange kernel. By iterating the Dyson-like equation the RPA screening naturally appears.

Figure 2 .

 2 Figure 2.15: The BSE with only the attraction kernel. By iterating the Dyson-like equation, the same interaction is repeated. No screening effects are included.

Figure 2 .

 2 Figure 2.16: A third-order term of the BSE with the GW kernel. This diagram shows a term of the screening beyond the RPA.

Figure 2 . 18 :

 218 Figure 2.18: Diagrammatic representation of the pp BSE kernel (2.137).

. 5 ) 3 , G eeh 3 and G hhe 3 , 3 , nine to G hhe 3 and nine to G eeh 3 .

 5333333 These terms describe G eee 3 , G hhh respectively. For each term, we have used the closure relation in Fock space to conserve the particle number. In total 6!/(3! 3!) = 20 different couples of three time-ordered operators can be formed, one that corresponds to G eee 3 , one to G hhh The term other orderings refers to all the other possibilities to create G hhe 3 or G eeh 3 terms.

3 and G eeh 3 ,

 3 see Fig 1.1. In order to have a more compact notation, we use here and in the following G h 3 and G e 3 , for G hhe 3 and G eeh 3 , respectively, i.e., the presence of the electron-hole pair is implied. Let us analyze the G e

. 11 )

 11 Using Eqs. (3.9) and (3.10) in Eq. (3.6) yields

3 in

 3 Eq. (3.24) involves the N -body ground state wave function as well as excited-state wave functions of the corresponding N + 1 and N -1 electron systems. It is therefore convenient to introduce an effective potential that links G e+h 3 to G e+h 03 , which is the non-interacting G e+h 3

Figure 4 . 1 :

 41 Figure 4.1: Diagrammatic representation of a general second-and third-order 1body self-energy obtained by iterating the multi-channel Dyson equation (4.15). By inserting multi-channel self-energy diagrams, Eqs. (4.22)-(4.24), all the second-and third-order proper skeleton 1-body self-energy diagrams are obtained.

Figure 4 . 2 : 2 2 n- 1 =

 4221 Figure 4.2: A fourth-order term created by the repetition of two second-order terms.

  diagrammatic representation of the 2-4-body multi-channel Dyson equation (4.34) is (4.42) To represent diagrammatically the multi-channel self-energy terms, we move to real space. Following the change of basis in Eq. (3.64), the multi-channel self-energy Multi-channel Dyson equation in Eqs. (4.39)-(4.41) becomes

Figure 4 . 3 :

 43 Figure 4.3: An example of a combination of multi-channel self-energy coupling terms that improve a single G 0 1 line.

Figure 4 . 4 :

 44 Figure 4.4: An example of a combination of multi-channel self-energy coupling terms that improve the interaction between the two particles. Here they create a diagram where the two interaction lines cross each other.

Figure 4 .

 4 Figure 4.6 shows the eh attraction potential screened at first-order in the irreducible polarizability χ (see Eq. (2.68)). Therefore, as for the 1-3-body case, the screening is naturally included in the theory although the self-energy only involves the bare Coulomb potential. Analyzing subsequent orders, one can realize that the screening goes beyond the RPA.

Figure 4 . 5 :

 45 Figure 4.5: An example of a combination of multi-channel self-energy coupling terms that improve the interaction between the two particles. Here they create a ladder diagram.

Figure 4 . 6 :

 46 Figure 4.6: An example of a combination of multi-channel self-energy coupling terms that improve the interaction between the two particles. Here they screened an eh attraction.

Figure 4 . 7 :

 47 Figure 4.7: An example of a combination of multi-channel self-energy coupling terms that creates a mixing term. Here, the eh attraction is between a hole and an electron that is already dressed.

Figure 4 . 8 :

 48 Figure 4.8: An example of a third-order term in the bare Coulomb interaction.

  (4.11)-(4.14)) is enough to obtain exact results for both the 1/4 and the 1/2 filling regime. The test of the 2-4-GF multi-channel Dyson equation is left as future work.

Figure 5 . 1 :

 51 Figure 5.1: Schematic representation of the Hubbard dimer. Each site is at the same energy that we fix to one through all the chapter, with -t the off-diagonal kinetic term (hopping term), and U the onsite Coulomb interaction.

  indices (i, j, l) are (a ↓ a ↑ b ↑), (a ↓ b ↓ b ↑), and (a ↑ b ↓ b ↑).

Figure 5 . 2 :

 52 Figure5.2: The addition part of the 1-and 3-body spectral function for the Hubbard dimer at 1/4 filling in the non-interacting limit (U = 0). Top panel: the exact spectral function A(ω). Bottom panel: the exact three-body spectral function A 3 (ω). The peak at ω = 4, which is present only in the three-body spectral function, it is related to a satellite. The spectra correspond to ϵ 0 = 1.

( 4 . 2 δ

 42 [START_REF] Csanak | Green's function technique in atomic and molecular physics[END_REF])-(4.14) we need an expression of the potential in the b/a basis. In this basis the Coulomb potential is given by v σiσm δ σj σ l (δ ĩ mδ jl + δ ĩj δ ml + δ ĩl δ j m -2δ ĩ mδ jl δ ĩj ).(5.10)

Figure 5 . 3 :

 53 Figure 5.3: The addition part of the 1-and 3-body spectral functions for the Hubbard dimer at 1/4 filling at weak interaction (U/t = 1). Top panel: the spectral function A(ω) obtained with various levels of theory: the exact 1-GF obtained from the multichannel Dyson equation (black solid line); the GW approximation (blue dotted line); the exact static 1-GF (green dashed line). Inset (a): zoom of the satellite peak. Bottom panel: the exact 3-body spectral function A 3 (ω). All spectra correspond to ϵ 0 = 1.

Figure 5 . 5 :

 55 Figure 5.5: The spectral function of the Hubbard dimer at 1/2 filling at weak interaction (U/t = 1) obtained with various levels of theory. Exact result obtained from the multi-channel Dyson equation (black solid line); the GW approximation (blue dotted line); the exact static 1-GF (green dashed line). The outer peaks are the satellites. (a) zoom of the removal satellite; (b) zoom of the addition satellite. The spectra correspond to ϵ 0 = -U/2 which guarantees the particle-hole symmetry.

Figure 5 . 6 :

 56 Figure 5.6: The spectral function of the Hubbard dimer at 1/2 filling at strong interaction (U/t = 4) obtained with various levels of theory. Exact result obtained from the multi-channel Dyson equation (black solid line); the GW approximation (blue dotted line); the exact static 1-GF (green dashed line). The outer peaks are the satellites. The spectra correspond to ϵ 0 = -U/2 which guarantees the particle-hole symmetry.

Figure 7 . 1 :

 71 Figure 7.1: Une représentation schématique de la spectroscopie de photoémission directe (PES) et inverse (IPES). L'élimination d'un électron (trou) du système laisse un phototrou (photoélectron), ce qui excite le système et conduit à la création de paires électron-trou.

Figure 7 . 2 :

 72 Figure 7.2: Une représentation schématique de la création d'une paire électron-trou due à l'absorption d'un photon.

. 2 )G hhh 3 , G eeh 3 et G hhe 3 . 3 , neuf à G hhe 3 et neuf à G eeh 3 .h 3 et G e 3 , pour G hhe 3 et G eeh 3 ,

 2333333333 Ces termes décrivent respectivement G eee 3 , Pour chaque terme, nous avons utilisé la relation de fermeture dans l'espace de Fock pour conserver le nombre de particules. Au total, 6!/(3!; 3!) = 20 couples différents de trois opérateurs ordonnés dans le temps peuvent être formés, un correspondant à G eee 3 , un à G hhh Le terme autres permutations fait référence à toutes les autres possibilités de créer des termes G hhe 3 ou G eeh 3 . Pour étudier les spectres de photoémission, nous ne sommes pas intéressés par les parties eee et hhh. Par conséquent, nous concentrons notre attention sur G hhe 3 et G eeh 3 . Afin d'avoir une notation plus compacte, nous utilisons ici et dans ce qui suit G respectivement, c'est-à-dire que la présence de la paire électron-trou est implicite. Analysons les termes G e 3 et G h 3 considérés dans l'équation (7.2) pour montrer qu'ils contiennent des informations sur les énergies de enlevement et d'ajout. Nous suivons une procédure similaire à celle que Csanak et al.[START_REF] Csanak | Green's function technique in atomic and molecular physics[END_REF] ont utilisée pour trouver la partie eh/he de la fonction de Green à deux particules.

  ),(3.10), et(3.18), un seul terme dans la somme reste après fixation des différences de temps. Nous notons que d'autres choix pour les différences de temps sont possibles afin d'obtenir le même ordre des opérateurs de création et d'annihilation mentionnés ci-dessus. Le résultat final ne dépend pas de ce choix.Avec les différences de temps données dans l'équation (3.23), nous obtenons l'expression suivante pour G e+h

3

 3

Figure 7 . 3 :

 73 Figure 7.3: Représentation diagrammatique d'une self-énergie générale du premier ordre obtenue par itération de l'équation de Dyson multi-canal (7.12). En insérant des diagrammes d'self-énergie multi-canal, les équations (7.21)-(7.16) permettent d'obtenir tous les diagrammes d'self-énergie proper skeleton du premier ordre, du deuxième et du troisième ordre.

Figure 7 . 4 :

 74 Figure 7.4: Un terme du quatrième ordre créé par la répétition de deux termes du deuxième ordre.

  .

7. 4

 4 Résultats : le dimère de Hubbard symétrique Afin d'illustrer la stratégie précédemment décrite, nous considérons le dimère de Hubbard symétrique. Il se compose de deux sites dégénérés, chacun contenant une orbitale. De plus, seuls les électrons sur le même site interagissent entre eux. Nous montrerons que l'self-énergie approchée multi-canal (voir les équations (4.11)-(4.14)) est suffisante pour obtenir des résultats exacts à la fois pour le régime de remplissage 1/4 et 1/2.L'hamiltonien correspondant au dimère de Hubbard symétrique est donné parH = -t i̸ =j=1,2 σ iσ ′ ĉiσ ′ ĉiσ + ϵ 0 i=1,2 σ niσ ,(7.33)dans lequel -t, U et ϵ 0 représentent respectivement l'énergie cinétique de hopping, l'interaction sur site (indépendante du spin) et l'énergie orbitale, et n iσ = c † iσĉiσ est l'opérateur de nombre. Nous avons explicité l'indice de spin σ dans l'équation ci-dessus. Nous notons que le degré de corrélation électronique dans le système est proportionnel au rapport U/t. Nous considérons une faible corrélation lorsque U/t = 1, tandis que lorsque U/t = 4, nous sommes dans le régime de forte corrélation. Le modèle est exactement soluble et nous permet donc de tester la précision de nos résultats dans les régimes de corrélation faible et forte.Les résultats pour le remplissage 1/4. Fig.7.6 montre les résultats pour U/t = 1 et Fig.7.7 montre les résultats pour U/t = 4Les résultats pour le remplissage 1/2. Fig.7.8 montre les résultats pour U/t = 1 et Fig.7.9 montre les résultats pour U/t = 4

Figure 7 . 5 :Figure 7 . 6 :

 7576 Figure 7.5: Représentation schématique du dimère de Hubbard. Chaque site a la même énergie que nous fixons à un tout au long du chapitre, avec -t terme cinétique hors-diagonale (terme de hopping) et U l'interaction de Coulomb sur site.

Figure 7 . 7 :Figure 7 . 8 :

 7778 Figure 7.7: La partie d'addition des fonctions spectrales à un et trois corps pour le dimère de Hubbard avec un remplissage de 1/4 et une interaction forte (U/t = 4). Panneau supérieur : la fonction spectrale A(ω) obtenue avec différents niveaux de théorie : la 1-GF exacte obtenue à partir de l'équation de Dyson multi-canal (ligne pleine noire) ; l'approximation GW (ligne pointillée bleue) ; la 1-GF statique exacte (ligne en pointillés verts). Le pic le plus à droite est un satellite. Panneau inférieur : la fonction spectrale exacte à trois corps A 3 (ω). Tous les spectres correspondent à ϵ 0 = 1.

1/ 2 filling

 2 In E.1 we report the eigenvalues and the coefficients of the eigenvectors for the symmetric Hubbard dimer at 1/2 filling. To have a more compact notation we definec = √ 16t 2 + U 2 , a = √ 2 c-U 16t 2 + (c -U ) 2 , A = 4t U -c , b = √ 2 c+U 16t 2 + (c + U ) 2 and B = 4t U +c . E i |b ↑ b ↓⟩ |b ↑ a ↑⟩ |b ↑ a ↓⟩ |b ↓ a ↓⟩ |b ↓ a ↓⟩ |a ↑ a ↓⟩ 2ϵ 0 + (U -c)/2

  . Different times choices might have different normalization prefactors. The normalization prefactor can be 1 N 2 or

	1 (N +1) 2 for Eq. (3.34), while it can be 1 N 2 or	1 (N -1) 2 for Eq. (3.35).

3.1.2 Dyson equation for G e+h

3

  Therefore, the head of the matrix (4.38) corresponds to the kernel in Eq. (2.106), and, for the other components, we obtain the following static approximation Σ 4p ijln;mokp = δ im δ jo vlnkp + δ lk δ np vijmo + δ lp δ nk vijom + δ io δ jm vlnpk + δ im δ np vjklo + δ im δ lk vjpno + δ jo δ np viklm + δ jo δ lk vipnm -δ im δ nk vjplo -δ im δ lp vjkno -δ jo δ nk viplm -δ jo δ lp viknm (4.39)

-δ io δ np vjklm -δ jm δ np viklo -δ io δ lk vjpnm -δ jm δ lk vipno + δ io δ lp vjknm + δ io δ nk vjplm + δ jm δ lp vikno + δ jm δ nk viplo Σ c jl,

mokp = δ jo vlmpk + δ lk vjpom + δ lp vjkmo + δ jm vlokp (4.40) Σc ijln,ok = δ jo vikln + δ lk vnoij + δ nk vloji + δ io vjknl (4.41)

The 3-GF and the 4-GF After removing all the redundant information, the non-interacting 3-GF can be written in the following matrix representation

with G 0 1 (ω) defined in equation (2.47), and G 0(3p) (ω) defined in equation (3.45). Equation (3.49) sets the space for the 1-3-GF multi-channel Dyson equation, it will be discussed in detail in section 4.2. The main idea is to use the three-body channel G 0(3p) 3

(ω) as a reservoir that is coupled with the 1-body channel through a multichannel self-energy. This strategy allows us to directly obtain the 1-GF without using equation (3.32).

The four-body Green's function

The 4-GF is defined by

As for the 3-GF, it can be decomposed into terms that depend on the number of interacting particles

where the superscript 0, 1, 2 refers to the number of interactions between particles, and δG 4 is the term where all four particles interact with each other. Let us analyse the non-interacting term G 0 4 which is important in order to construct a Dyson-like equation. Thanks to Wick's theorem, it can be written as the following determinant

which generates 4! = 24 terms, each composed of a product of four 1-GF. The 4-GF depends on eight times or seven time differences when the Hamiltonian is time-independent. The total number of permutations of the field operators in Eq. (3.50) due to the T operator is 8! = 40320. Different choices of time ordering bring different orders of the field operators and, therefore, different physical information. The 4-GF describes the propagation of four particles (electrons or holes) and it can be split into five components: G eeee

The eehh 4-GF

In this thesis, we aim to use the 4-GF to improve the eh BSE by treating single and double excitation on equal footing. For this reason, we analyze only the eehh component that, for simplicity, we call G eh 4 . With the time definition of Eq. (3.13), we define the G eh 4 with the following time differences

obtaining

Following the same ideas we used for the 3-GF in the previous section we obtain the following spectral representation

with

Equation (3.55) has the same poles as the eh channel of the 2-GF (see Eq. (2.79)). For this reason, the eehh part of the 4-GF contains the eh part of the 2-GF. To show that, we follow the same strategy adopted to obtain the 1-GF from the 3-GF (Eq. (3.32)),

The 3-GF and the 4-GF the spectral representation in Eq. (3.55) can be rewritten as

in which

We notice that the 4-GF remains the same, except for a sign, under the permutation of the following couples of indices: (i, j), (l, n), (m, o), and (k, p). To avoid double counting and to make the 4-GF invertible the conditions i > j, l > n, m > o, and k > p are used.

With these restrictions, the non-interacting 4-GF can be reduced as follows: among all the G 1 G 1 G 1 G 1 combinations that describe the propagation of four particles, only the term defined in Eq. (3.61) is non-zero and, therefore, is considered. All the G 1 G 1 G 1 G 1 combinations that describe the propagation of two particles contains redundant information, i.e. the poles of the eh channel of the 2-GF repeated multiple times, which makes its matrix representation non-invertible. To avoid this problem we can reduce the space of definition of the eh channel to the one of [G 0 jo G 0 lk ](ω) only, without losing information. Combinations where all the G 0 1 are static are not considered since they do not have poles.

With all these restrictions, and the definitions of the two-and four-body correlation function as L0 2 = iG 0(eh) 2

, and L0 4 = iG 0(eh) 4

, the non-interacting 4-body correlation function simply reads

where

.69)

The above equations hold for any independent-particle picture. This is the starting point for the multi-channel Dyson equation for the 2-and 4-body channels. We will discuss it in the next chapter.

CHAPTER 4

The multi-channel Dyson equation

In this chapter we present the main result of this thesis: the multi-channel Dyson equation. It is a Dyson equation that couples different GF channels through a multichannel self-energy that has the advantage of being static and containing only the bare Coulomb interaction. This is an approximation, since in general the self-energy can be frequency dependent. This equation couples n-GF channels that have the same poles in the spectral representation. Previously we have shown that this is the case for the 1-GF and 3-GF (equations (2.45) and (3.24)), and, the for the 2-GF and 4-GF (equations (2.79)) and (3.55)) cases. In general, every n-GF can be coupled with the (n±2)-GF.

General case

The multi-channel Dyson equation is defined as We note that each diagram on the right-hand side of Eqs. (4.25) and (4.26) is the sum of two diagrams. Because of the restriction of the space mentioned before (see Eq. (4.9)), both diagrams are needed, and there is no double counting.

Appendix C shows a detailed calculation on the equality between diagrams calculated from the multi-channel Dyson equation (4.1) and from the Feynman rules discussed in section 2.2.6.

It can be verified that inserting Eqs.(4.22)-(4.24) into the diagram on the right in Fig. 4.1 yields all ten third-order proper skeleton diagrams, which include both bubble and ladder diagrams. For example, a bubble diagram is obtained as follows, To study which diagrams of the eh BSE appear from the multi-channel Dyson equation (4.34), we iterate it and inspect the head of the matrix. The first iteration is equivalent to the first iteration of the BSE with the RPAx kernel (see Eqs. (2.104), (2.106)). At the second iteration, the coupling terms enter in the head of G 4 . They have three different tasks: 1) they add correlation to a single G 0 1 line (see Fig. 4.3). Due to this, the non-interacting two-body term (Eq. (4.36)) has to be considered at the Hartree-Fock level. Adding more correlation to L0

2 will create double-counting problems; 2) they improve the interaction between the two particles (see Figs.

(5.15)

As was the case at 1/4 filling, the non-interacting 3-body channel contains information about satellites, i.e., the terms corresponding to the poles at ω = ϵ 0 ± 3t, have nonvanishing amplitudes.

The HF 3-GF is given by

Useful Relations

In this appendix, we demonstrate some useful relations. Let us start from the identity

Using the Martin-Schwinger hierarchy Eq. (2.54)

which can be substituted in Eq. (A.1). We hence arrive at

We are now going to analyse the inverse for a 2-GF. For these functions the inverse is defined starting from the non-interacting case and extended at the interacting case. Let us start with the non-interacting 2-GF

We define as inverse the following expression

First of all, let us to emphasize that the term "inverse" is improper here because the two combination of deltas are not the identity. This definition can be extended to the interacting part

To understand how to use it we exploit the symmetry property of G n

where every permutation to near indices brings a minus. Notice that we cannot permute indices with and without ′ because they depend on different field operators (construction or annihilation). With this symmetry we have

Therefore, if we want to apply G 2 G -1 2 we have two rules 1. we can apply it to indices related to the same field operator;

2. we have to add 1 2 .

δG 3

In our journey to find a Dyson equation for the 3-GF we have followed also another strategy, which is based on the work of H. Hübener [START_REF] Huebener | Second-order response bethe-salpeter equation[END_REF] for the core of G 3 . We report it here for completeness

in the last equality Eqs. (3.2), (3.3), and (3.4) have been used. The other equation obtained from the functional derivative is

where

In Eq. (B.2), it is important to notice that Ξ 3 appears in the term without L 3 , while the term with L 3 interacts through a BSE kernel. Moreover, neglecting the variation of the potential with respect to the 1-GF, as is usually done, and considering the GW self-energy, the three-body kernel is zero, i.e., Ξ GW 3 = 0. In general, the lowest order in the Coulomb interaction present in Ξ 3 is the second order. This is due to the two-body nature of the potential.

Let us continue with the exact equations. To find an equation for δG 3 , we make use of Eqs (B.1), and (B.2), obtaining

Here it is clear that the lower order in the interaction is the second. Early, a lot of work has been done to try to find approximations to reduce Eqs. (B.2) and (B.4) to a Dyson-like equation, but these approximations are not very practical.

Diagram calculation

In this appendix, we show how to formally compare diagrams obtained from the multi-channel Dyson equation 4.15 and the ones obtained from the Wick theorem.

The second order in the potential 1-body self-energy obtained from the multichannel Dyson equation is shown on the left of Fig. 4.1. Making use of eqs. (4.9), (4.12), and (4.13), we obtain

The last equality is determined by the condition (

Where in the first equality, the second term is obtained using c ↔ c ′ for the second term in the numerator, and to obtain the last term the relation

The other second-order contribution, Fig.C.2, is obtained when one potential is direct (inverse) and the other inverse (direct). Following the same strategy as in Eq.(C.2), this self-energy term is

To obtain the equation related to each diagram directly from the diagram, we start from the equation in position space. Following the Feynman rules explain in section 2.2.6, the diagram in Fig. C.1 is expressed by the following equation moving to the frequency domain and using Eq.(3.39)

finally, projecting in a basis and using Eq. (2.136)

The equation for the diagram in Fig.

following the same procedure as for the other terms, we obtain

Appendix C

The calculation of all terms is space demanding and repetitive. We will calculate two third-order proper skeleton diagrams out of the total 10, i.e., those represented in

The other three terms are calculated with the same expression but with different potential indices. Accordingly, to the diagrams in Fig. C.3, the other terms are obtained with the substitution

using the same strategy as in Eq.(C.2) these differences create Therefore, the sum of the four diagrams gives

like before, the other terms are calculating changing the self-energy indices according to the diagrams [START_REF] Sergey | All-electron selfconsistent gw approximation: Application to si, mno, and nio[END_REF] using the same strategy as Eq.(C.2) these differences create

summing the four components we obtain We are now going to do the same calculation but starting from the diagrammatic definition. To calculate the diagrams we will follow the same idea as in the second order, namely, we start in position space and we project in a basis to compare it with the multi-channel results, eqs. [START_REF] Viktor | Application of quantum field theory methods to the many body problem[END_REF] the above equation can be split into two cases, t > t ′ and t < t ′ . For construction, t ′′ has to be between the two. Defining Σ + for the case t > t ′ , making use of Eq. (2.31), we obtain

since t ′′′ is a generic time between t and t ′ , it is possible to consider it equals to t ′′ . In this way, two terms with the following time dependence G 0(e) (τ )G 0(e) (τ )G 0(h) (-τ ) appear. In the frequency space, they can be written accordingly to Eq. (3.45), and the self-energy becomes

By projecting in a basis and exploiting the definition of v ijkl , Eq.(2.136), we obtain

that is the second part of Eq.(C.14). To obtain the other term, we have to follow the same procedure but with the condition t ′ < t. The difference in sign of Eq.(C.

As before, to study it we dived in the t > t ′ and t ′ > t cases. We analyze the case with t > t ′ . Using the equation (2.30), we find

Expressing the self-energy in the frequency domain

127 projecting in a basis, and using the Eq.(2.136), we obtain

that is the second part of the right-hand side of Eq.(C.18). To obtain the other term we need to study the condition t < t ′ . Whit it, the 1-GF to split in two is a hole part, i.e., Eq. (2.31) has to be used. Doing all the calculations as before, we obtain the first term of Eq.(C.18).

APPENDIX D

From the multi-channel Dyson equation to an eigenvalue problem

In this appendix, we show how to obtain equation (4.30) starting from the multichannel Dyson equation (4.1). We start by defining a quantity that links the interacting and the non-interacting 3-GF as

Eigenvalues and eigenvectors of the symmetric Hubbard dimer

In Refs. [START_REF] Romaniello | The self-energy beyond GW: Local and nonlocal vertex corrections[END_REF] the eigensystem of the Hamiltonian in Eq. (5.1) is solved in the site basis, where the eigenstates of the system are linear combinations of Slater determinants composed by sites 1 or 2 with spin ↑, or ↓.

Here we report the results in the bounding/antibonding basis since in this basis the non-interacting 1-GF is diagonal. Therefore, it is the basis where the multi-channel Dyson is defined.

It is obtained with the change of basis and represent both the site and the spin.

1/4 filling

The eigenstate of the Hamiltonian at 1/4 filling are |b ↑⟩, |b ↓⟩, |a ↑⟩, |a ↓⟩ with energies ϵ 0 -t, ϵ 0 -t, ϵ 0 + t, ϵ 0 + t, respectively. The degeneracy of the ground state can be removed applying a small magnetic field, this allow us to choose freely the ground state. We chose the state |b ↑⟩ as ground state.