Close
About
FAQ
Home
Collections
Login
USC Login
Register
0
Selected
Invert selection
Deselect all
Deselect all
Click here to refresh results
Click here to refresh results
USC
/
Digital Library
/
University of Southern California Dissertations and Theses
/
Coulomb interactions and superconductivity in low dimensional materials
(USC Thesis Other)
Coulomb interactions and superconductivity in low dimensional materials
PDF
Download
Share
Open document
Flip pages
Contact Us
Contact Us
Copy asset link
Request this asset
Transcript (if available)
Content
University of Southern California Doctoral Dissertation Roelof E. Groenewald Coulomb Interactions and Superconductivity in Low Dimensional Materials COULOMB INTERACTIONS AND SUPERCONDUCTIVITY IN LOW DIMENSIONAL MATERIALS by Roelof E. Groenewald supervised by Dr. Stephan Haas A Dissertation Presented to the FACULTY OF THE USC GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulllment of the Requirements for the Degree DOCTOR OF PHILOSOPHY (PHYSICS & ASTRONOMY) December, 2018 ...kry my op die pad, die pad waar die dromers loop, son oor ons skouers ons oe vol hoop, kom saam, kom saam dis al wat dit kos, ek en jy is gemaak om spore te los. - Elvis Blue I Abstract Screening of the Coulomb interaction is generally suppressed in low dimensional materials and heavily aected by substrates. This in turn aects the interaction between electrons, which plays a major role in many-body eects. In the eld of superconductivity, the Coulomb interaction is usually reduced to an approximate static parameter (retarded Morel-Anderson Coulomb potential), which describes only the eective Coulomb repulsion, and therefore is responsible for reducing the superconducting transition temperature. Our calculations of reveal that the Coulomb interactions, renormalized by the re- duced layer thickness and the substrate properties, can shift the onset of the electron- phonon driven superconducting phase in monolayer MoS 2 , but do not signicantly aect the critical temperature at optimal doping. We overcome the inadequate handling of the Coulomb interaction as a static parameter and present an ab initio based dynamical Coulomb description. We apply this method to MoS 2 to account, in a realistic manner, for intrinsic material screening, substrate screening as well as dynamical screening processes. Using this method we can reliably describe the plasmonic excitations and their coupling to the electrons. Plasmonic excitations, brought on by the dynamical screening of the Coulomb interac- tion, behave fundamentally dierent in layered materials in comparison to bulk systems. They form gapless modes, which in turn couple at low energies to the electrons. Thereby they can strongly in uence superconducting instabilities. We show how these excitations can be controlled from the outside via changes in the dielectric environment or in the doping level, which allows for external tuning of the super- conducting transition temperature. By solving the gap equation for an eective system, we nd, generally, that the plasmonic in uence can both strongly enhance or reduce the transi- tion temperature, depending on the details of the plasmon-phonon interplay. We formulate simple experimental guidelines to nd plasmon-induced elevated transition temperatures in layered materials. Lastly, we provide details of how this approach can be used to explore the interplay of plasmons and phonons in real materials (such as MoS 2 ) and nd an optimal ratio between substrate screening and electron doping which maximizes the superconducting critical tem- perature of the resulting induced phononic and plasmonic superconducting state. II Contents 1 Introduction 1 2 General Theory Overview 6 2.1 Bloch Theorem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 2.2 Single Particle Properties . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.2.1 Tight-Binding Method . . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.2.2 Density Functional Theory . . . . . . . . . . . . . . . . . . . . . . . 11 2.2.3 Wannier Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 2.3 Many Particle Properties . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 2.3.1 Plasmons . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 2.3.2 Phonons . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 2.4 Particles Interacting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 2.4.1 Particle Propagators and Dyson's Equation . . . . . . . . . . . . . . 32 2.4.2 Coulomb Interaction . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 2.4.2.1 Electron-plasmon coupling . . . . . . . . . . . . . . . . . . 38 2.4.3 Electron-Phonon Interaction . . . . . . . . . . . . . . . . . . . . . . 40 2.4.3.1 Electron-phonon coupling . . . . . . . . . . . . . . . . . . . 40 2.4.4 Spectral Function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 3 Superconductivity Theory 43 3.1 BCS-theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 3.1.1 Cooper instability . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 3.1.2 BCS Gap Equation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 3.2 DFT for Superconductors . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 3.2.1 Approximation for the exchange-correlation kernel . . . . . . . . . . 56 3.2.2 Energy Averaged Kernels . . . . . . . . . . . . . . . . . . . . . . . . 60 3.2.3 Derivation of McMillan's equation in SC-DFT . . . . . . . . . . . . . 63 3.2.4 Plasmonic Eliashberg function in 2d . . . . . . . . . . . . . . . . . . 69 III CONTENTS CONTENTS 4 Model System - 2d Square Lattice 73 4.1 Model Description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73 4.1.1 Electronic Dispersion and DOS . . . . . . . . . . . . . . . . . . . . . 76 4.1.2 Phononic Properties . . . . . . . . . . . . . . . . . . . . . . . . . . . 77 4.2 Plasmonic Properties . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 4.2.1 Plasmon Dispersion . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 4.2.2 Plasmonic Eliashberg Function . . . . . . . . . . . . . . . . . . . . . 80 4.2.3 Self-energy Contribution . . . . . . . . . . . . . . . . . . . . . . . . . 84 4.3 Superconductivity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87 4.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93 5 Material Realistic Calculations for MoS 2 96 5.1 Ab initio Calculations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 5.1.1 DFT - electronic dispersion . . . . . . . . . . . . . . . . . . . . . . . 98 5.1.2 DFPT - phonon dispersion . . . . . . . . . . . . . . . . . . . . . . . 98 5.2 Model - Wannier Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . 100 5.3 Substrate Screening . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 5.4 Plasmonic Properties . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105 5.5 Superconductivity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 5.6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115 6 Conclusions and Outlook 116 6.1 Outlook . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 Appendices 123 A Theory Overview Sections 123 A.1 Hartree-Fock Equations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123 A.2 Attractive electron-electron interaction from phonon scattering . . . . . . . 124 B Computational Parameters 125 B.1 MoS 2 Ab initio Calculations . . . . . . . . . . . . . . . . . . . . . . . . . . . 125 B.2 MoS 2 Wannier . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 B.3 2d Square Lattice SC-DFT calculations . . . . . . . . . . . . . . . . . . . . 129 B.4 MoS 2 Plasmon calculations . . . . . . . . . . . . . . . . . . . . . . . . . . . 129 C Code Implementation Benchmarks 129 C.1 SC-DFT - Step-Like Density of States . . . . . . . . . . . . . . . . . . . . . 130 C.2 SC-DFT - MoS 2 Monolayer . . . . . . . . . . . . . . . . . . . . . . . . . . . 131 C.3 Three-band model for TMDC's . . . . . . . . . . . . . . . . . . . . . . . . . 132 Bibliography 134 IV CONTENTS CONTENTS Acknowledgments 150 V Chapter 1 Introduction The Nobel laureate, P. W. Anderson , wrote a highly celebrated article titled, \More is Dierent - Broken symmetry and the hierarchical nature of science" [2]. In the article he argues against the notion that the purpose of the physicist is to uncover the fundamental laws of nature. After the foundation of quantum mechanics was built, the brilliant P. Dirac y is alleged to have said `the rest is chemistry'. Similarly, many scholars even until today, believe biology is `just' applied chemistry, and psychology `just' applied biology etc. Anderson explains that this hierarchical view of science is misleading and the physical world cannot be fully understood from its' fundamental laws alone. Instead complex systems often exhibit behavior fundamentally dierent than the behavior of their constituent parts. As Pavarini and Koch noted in Ref. [4], gold is made up of gold atoms, but a gold atom does not conduct electricity or even appear golden in color. It's only when many gold atoms are placed together that we recover the properties of a gold sample. Anderson said of this emergence of new properties, \we can see how the whole becomes not only more than but very dierent from the sum of its parts." Perhaps the most compelling example against the constructionist's point of view lies in the history of superconductivity. All the fundamental laws necessary to explain super- conductivity have been known since the 1920's but it wasn't until 1957 that a theory was presented to explain it and even then, that theory was built phenomenologically and not based on rst principles. Today, we understand that a large system can spontaneously break symmetries present in the parts that formed it. For example, nature has a continuous translational invariance which is re ected in the positions of atoms or molecules in a liquid state z . But when a jointly received the 1977 Nobel prize in physics along with Sir N. Mott and J. van Vleck \for their fundamental theoretical investigations of the electronic structure of magnetic and disordered systems" [1] y shared the 1933 Nobel prize in physics with E. Schr odinger \for the discovery of new productive forms of atomic theory" [3] and also predicted the existence of antimatter z this examples comes from Ref. [5], it should be noted that liquids have short range order but here we 1 CHAPTER 1. INTRODUCTION large enough collection of such atoms or molecules is cooled down below a critical temper- ature (1064 o C in the case of gold) the atoms condense into a less symmetric state, forming a crystal lattice with discrete translational invariances. This reduction of symmetry goes hand in hand with the emergence of new properties not present in the original system. The rigidity of solid crystals is an obvious example. The exact nature of the broken symmetry is not always obvious, though, as it is the case with superconductivity , and even less so the mechanism that leads to this symmetry reduction. In my view the goal of the con- densed matter physicist is to determine the mechanisms that enable symmetry breaking and understand how these mechanism emerge from the underlying fundamental laws of physics. To this end it is not our goal to nd exact solutions of the fundamental equations but rather construct models that enable us to study the exact mechanism we are interested in. The work presented here is entirely built on this philosophy. We use rst principle calculations to construct model systems that inherit the crucial properties of nature, which we can then reduce to minimal systems that exhibit the mechanisms we are interested in studying. Specically, as the title suggests, this thesis addresses some questions con- cerning properties of the superconducting state of systems with reduced dimensionality. In particular, we consider two dimensional (2d) systems. In the general context of com- putational physics, researchers often exploit a continuous, translational invariance in a particular direction of a three dimensional (3d) system to reduce the computational space to two dimensions. This is not the situation in this thesis. We are interested in systems consisting of single atomic layer crystals and therefore only extend into 2d y . In this intro- ductory chapter we will describe why such systems are of interest. We'll also give a brief description of what superconductivity is and attempt to convince the reader that, due to the nature of 2d materials, their superconducting states are particularly interesting. Superconductivity was rst discovered by Kamerlingh Onnes in 1911 [9] after he suc- ceeded in liquifying helium. He experimented with materials at liquid helium tempera- ture and observed that a sample of mercury loses all its' electrical resistance below 4:1K (269:05 o C). This is the dening characteristic of a material in a superconducting state - zero resistivity. The consequence is that electrical current can ow through the material with no loss in the form of Joule heating. For all practical purposes a current in a super- conducting loop can ow perpetually. Another, characteristic of superconductors lies in the so called Meissner eect [10]. If a magnetic eld is applied to a superconductor cur- rent will ow in the conductor to perfectly expel the magnetic eld from the inside of the superconductor, similar to theoretical perfect diamagnetism. This perfect diamagnetism mean translational invariance of the thermodynamic state the broken symmetry in this case is a loss of gauge invariance, discussed in detail in chapter 3 y see Ref. [6] for a list of 2d crystals. The paper includes in its authors the Nobel laureates, Novoselov and Geim, who were the rst to obtain a monolayer of graphene [7,8]. 2 CHAPTER 1. INTRODUCTION can only be sustained up to a critical eld strength, H c , which is found to be temperature dependent. The temperature at which the critical eld vanishes is called the critical tem- perature of the superconductor, i.e. H c (T c ) = 0. At the critical temperature the system undergoes a phase transition from the normal metallic state to the superconducting state. Experiments have shown that H c (T ) is a continuous function of T indicating that the normal-to-superconductor transition is a second-order phase transition . In 1913, lead was found to become superconducting at 7:2K, and in 1930 the critical temperature of niobium was shown to be 9:2K [11]. More than 50 elements have been shown to superconduct in labs around the world [12], since then, and many more metal compounds as well. The current conrmed record for the highest temperature superconductor found, is from a 2015 experiment by Drozdov et al. in which they found sulfur hydride (H 2 S) under 90 gigapas- cals become superconducting at 203K [13]. Excitingly, a very recent report from Thapa et al. [14] claims that they witnessed silver nano-particles embedded in a gold lattice become superconductive under ambient pressure at 236K y . In 1957 Bardeen, Cooper and Schrieer [15] came up with an explanation for this ex- otic state of matter. They showed that the interaction between electrons via exchange of phonons becomes attractive at low temperatures, which leads to the formation of quasi- particles (called Cooper-pairs) consisting of pairs of electrons. These quasi-particles are eectively bosons and therefore can form Bose-Einstein condensates, which leads to the macroscopic current carrying state seen in superconductors (similar to the frictionless ow of a super uid). Superconductors that can be explained via this theory are now referred to as conventional superconductors. However, the electron-phonon interaction is not the only mechanism that can lead to a superconducting phase. Alternative, uncon- ventional mechanisms with electron-electron pairing occurring due to plasmons, excitons, spinons and vibrons have also been considered [16{18]. The necessary condition for strong electron-electron coupling by any bosonic mode is that the bosons exhibit low energy dis- persions. This is necessary since Cooper-pairs necessarily form between electrons close to the Fermi-level. In 2d the plasmonic dispersion exhibits a characteristic low-energy acous- tic mode!(q)/ p q originating from low-momentum electron scattering [19,20], which has been observed experimentally [21{23] and studied extensively from a theoretical point of view [24{28] in graphene. Furthermore, it has been predicted that additional linear plas- mons with!(q)/q arise due to high-momentum scattering processes between degenerated valleys such as K and K 0 in graphene [29]. Coupling of the electrons with such intrinsic gapless bosonic modes may lead to instabilities, such as charge density wave and supercon- ducting phases [30{33], similar to the eect of phonons. For this reason, 2d materials are of particular interest in investigating the possibility of plasmon induced superconductivity. Other authors have focused on plasmonic eects on superconductivity in 3d [31,32,34,35] but the plasmon dispersion in 3d is characterized by a high-energy optical modes making rst order transitions are characterized by a discontinuous order parameter while the superconducting order parameter follows Hc y this experiment is yet to be conrmed by another group 3 CHAPTER 1. INTRODUCTION their coupling to electrons weak. Ever since the successful exfoliation of singe layer graphene, interest in 2d materials has been high due to the promising application of these materials in next generation elec- tronic devices. Currently transistors (the key component to all processors) are made from silicon with metal-oxide deposits to allow for eld eect gating. These components rely on etching features with length-scales of around 14nm [36], and considering that the silicon lattice constant is about 0:54nm [37] it is clear the limit of how much smaller transistors can get, is close. On the other hand if a 2d transistor can be built, stacking them allows much denser computational power per unit volume. Unfortunately, graphene is not suited for this task since it doesn't have a band-gap and therefore can't reach the on-o current ratio needed for electronic applications [38]. For this reason monolayer transition metal dichalcogenides, such as MoS 2 , have been receiving increasing attention. These materials do possess a band-gap, and a transistor based on MoS 2 has already been demonstrated [39]. Furthermore, collective electronic excitations are of great interest in these low-dimensional materials which are characterized by reduced dielectric screening of Coulomb interactions. As a prominent example, plasmon modes in layered systems might form the basis to build optical devices, wave guides or so called plasmonic circuits [23,27,40,41]. Due to reduced screening, layered materials promote sizable electron-electron and electron-phonon inter- actions, leading to a plethora of non-trivial phenomena. In addition to a variety of many- body instabilities, such as charge or superconducting order [42{45], we also nd pronounced many-body excitations such as excitons with binding energies on the eV scale [46,47]. Based on changes to the electron-electron (Coulomb) interaction, these inherent many- body eects can be tuned from the outside by adjusting the doping level or the dielectric environment [48{51]. For example, it is well known that the instantaneous Coulomb inter- action is repulsive and thus decreases the transition temperature (T c ) in conventional su- perconductors [52{54], but unconventional paring mechanisms based on the instantaneous interaction have been shown to increase T c [55]. The experimental capability to design Van-der-Waals heterostructures in a Lego-like fashion [56] thus promises a fundamentally new opportunity to design material-properties on demand by tuning the interactions from the outside. Critical temperatures for the transition to superconducting or magnetic order may thereby be strongly increased. Outline This thesis is organized in order to guide the reader step-wise through it. Chapter 2 provides an overview of the ab initio techniques used to obtain electronic and phononic properties, as well as the general many-body theory employed in this thesis. We stress the emergence of collective excitations and how these eects can be calculated. An in-depth discussion of 4 CHAPTER 1. INTRODUCTION superconductivity and the superconductor-density functional theory, follows in chapter 3. Chapter 4 is dedicated to a model based study of the impact electron-plasmon coupling has on superconductivity. Then in chapter 5 we attempt to bind everything together with an in-depth study of monolayer MoS 2 , using ab initio calculations to determine the plasmonic and phononic properties and determining the interplay of their respective couplings with the electrons on the superconducting transition temperature of the material. Finally, in chapter 6 we provide a conclusion of the main results discussed throughout the thesis and an outlook for the future. 5 Chapter 2 General Theory Overview In this chapter we brie y discuss the crucial aspects of solid state physics on which this dissertation is built. Naturally, the reader is assumed to be familiar with the theory but for the sake of completeness we include a rather detailed review. First, we'll focus on the problem of solving the Schr odinger equation for an electron in a periodic potential (such as a crystal), discussing both simplied models and ab initio techniques. Then we'll proceed to discuss many-body interactions in such systems and their corresponding properties. 2.1 Bloch Theorem The subsequent discussion of Bloch's theorem follows closely the solid state physics textbook by Grosso and Parravicini [57] and lecture notes by Cao [58]. Consider an electron in a crystal, with a lattice translation vector ~ R, being acted upon with the Hamiltonian, H = ~ p 2m +V (~ r): (2.1) If we assume the electron's potential energy, V (~ r), comes entirely from the ion cores making up the crystal structure, it stands to reason that the function V (~ r) should have the same periodicity as the lattice, i.e. V (~ r) = V (~ r + ~ R). Now let T ~ R be the translation operator of vector ~ R, then for an arbitrary function f(~ r) we have, [T ~ R ;V (~ r)]f(~ r) =V (~ r + ~ R)f(~ r + ~ R)V (~ r)f(~ r + ~ R) = 0: (2.2) Remembering that the kinetic energy is invariant under translations, we see that T ~ R commutes with the Hamiltonian, which means both operators can be diagonalized simul- taneously. Let (~ r) be an eigenfunction of T ~ R expanded in a plane wave basis (that takes into account the boundary conditions of the problem) to give (~ r) = P q C q e iq~ r , then the corresponding eigenvalue is given by, 6 CHAPTER 2. GENERAL THEORY 2.1. BLOCH THEOREM T ~ R (~ r) =T ~ R X q C q e iq~ r = X q C q e iq~ r e iq ~ R =t ~ R X q C q e iq~ r : We see that t ~ R = e iq ~ R must be a constant which in turn means we need q ~ R = 2n + constant. Let q = k + G, where G is a reciprocal lattice vector i.e. G ~ R = 2n. Then the eigenvalue is, t ~ R =e ik ~ R , and the eigenvector associated with it is given by, k (~ r) = X G C k+G e i(k+G)~ r =e ik~ r X G C k+G e iG~ r =e ik~ r u k (~ r): (2.3) The eigenvalues and eigenstates can thus be labeled by the vector k. We also notice immediately thatu k =u k+G therefore the Bloch wavefuntions are periodic in the reciprocal space and we only need to study those vectors that can not be connected by G. The set of these vectors are called the rst Brillouin zone. Recall that T ~ R commutes with the Hamiltonian (as shown in Eq. 2.2) and therefore the k 's are also eigenfunctions of the Hamiltonian. This, nally, leads us to state Bloch's Theorem: Any physically acceptable solution of the Sch odinger equation in a periodic potential takes the form of a traveling plane wave modulated by an appropriate function with the lattice periodicity. We mentioned above that the expansion of the eigenfunctions needs to take into account the boundary conditions of the problem. In general we will employ periodic boundary conditions over a macroscopic crystal region, 0 < ~ r < N i ~ R i , where N i is the number of unit cells in the i'th direction in the crystal. Periodic boundary conditions amount to setting, (~ r +N i ~ R i ) = (~ r) for some eigenfunction, , i.e. the points ~ r and N i ~ R i are considered to be physically equivalent. Applying Bloch's theorem to this we see that, k (~ r +N i ~ R i ) =e ik(~ r+N i ~ R i ) u k (~ r +N i ~ R i ) =e iN i k ~ R i k (~ r); which means we need k ~ R i = 2 N i . This forces a quantization of the allowed values of k, and we see that this discretization is inversely related to the number of unit cells we consider in the crystal. This will become important later on when we consider the number of k-points needed to achieve converged results (results not aected by nite size eects) for numerical calculations. 7 2.2. SINGLE PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY 2.2 Single Particle Properties The subsequent discussion of the band theory of crystals follows closely from the solid state physics textbook by Grosso and Parravicini [57] chapters IV and V as well as lecture notes by Paxton [59]. The discussion of Wannier Functions follows the articles by Marzari and Vanderbilt [60] and Souza et. al. [61] as well as the Wannier90 user guide [62]. The cor- responding chapters in the PhD thesis of Dr. R osner [63] were also closely followed. The Hamiltonian for an interacting (non-relativistic) system of electrons (with positions r i and momenta p i ) and nuclei (with positions R I and momenta P I ), as can be found in atoms, molecules or solids (including crystals), is given by, H tot = X i p 2 i 2m + X I P 2 I 2M I + X i X I z I e 2 jr i R I j + 1 2 X i6=j e 2 jr i r j j + 1 2 X I6=J z I z J e 2 jR I R J j : (2.4) We are rstly interested in discussing the electronic properties of crystals and therefore we will neglect the kinetic energy of the nuclei (we will return to this when discussing phonons in Section 2.3.2). In this crude approximation the nuclei are taken to be xed in the equilibrium conguration of the crystal and their positions become parameters of the Hamiltonian. Furthermore, in this approximation the nuclear potential term is simply a constant that rigidly shifts the spectrum of the Hamiltonian, and can therefore also be ignored in considering electronic properties (where only dierences in energies matter). The electronic Hamiltonian can then be written as, H elec = X i p 2 i 2m + X i V nucl (r i ) + 1 2 X i6=j e 2 jr i r j j ; (2.5) where V nucl (r i ) is the potential energy given electron i due to the xed nuclei. We will now review two common methods, central to this thesis, used to nd the electronic dispersion of a crystal. Firstly, we'll employ the tight-binding method to show how we can construct model systems that will be used heavily in Chapter 4. Next we'll discuss how we can approximately solve the many electron case from rst principles using density functional theory and nally how we can combine the tight-binding method and DFT using Wannier constructions. 2.2.1 Tight-Binding Method The tight-binding method was rst suggested by Bloch in 1928 and Slater and Koster [64] published a general recipe for constructing tight-binding models. Sutton et al. [65] have also shown that one can arrive at the tight-binding method through a rigorous approximation of the density functional theory. The tight-binding method for us will consist in expanding 8 CHAPTER 2. GENERAL THEORY 2.2. SINGLE PARTICLE PROPERTIES the crystal eigenstates as a linear combination of atomic orbitals (LCAO). The specic orbitals are determined by the atoms making up the crystal. For the sake of simplicity we will derive the specics of the method for a homogeneous crystal with only one atom per unit cell. Let i (~ r ~ R m ) be an atomic orbital with quantum numbersi and energy i from the atom centered at unit cell ~ R m (of N total unit cells in the crystal). We construct the Bloch sum of wavevectors i (k;~ r) = 1 p N X ~ Rm e ik ~ Rm i (~ r ~ R m ); (2.6) We then take enough Bloch sums of vector k to describe the occupied (and perhaps lowest conduction) bands of the crystal we are considering to expand the crystal wavefunction, in the form, nk (~ r) = X i c n i (k) i (k;~ r): (2.7) With these eigenfunctions in hand we can replace the Hamiltonian eigenvalue problem with the linear eigenvalue problem, (H ij E n S ij )c n i = 0; (2.8) where E n is the n'th eigenenergy of the crystal, H ij are the Hamiltonian matrix elements in the LCAO basis and S ij the overlap matrix elements given by, H ij (k) =h j (k;~ r)jHj i (k;~ r)i; S ij (k) =h j (k;~ r)j i (k;~ r)i: (2.9) In the case that the atomic orbitals are heavily localized, the overlap of orbitals from one atom to the next becomes negligible and we can assume the atomic orbitals are orthonor- mal. In this case the Bloch sums are also orthonormal and the overlap matrixS(k) becomes the unit matrix. Next we consider the Hamiltonian matrix elements, H ij (k) = 1 N X ~ Rm; ~ Rn e ik( ~ Rn ~ Rm) D i (~ r ~ R m ) H j (~ r ~ R n ) E : Since the Hamiltonian has translational invariance we can choose ~ R m to be the reference unit cell i.e. ~ R m = 0, which allows us to drop the sum over it and cancel the factor of 1=N, to get H ij (k) = X ~ Rn e ik ~ Rn D i (~ r) H j (~ r ~ R n ) E : (2.10) We can now further simplify this expression if we assume the crystal potential can be expressed as the sum of spherically symmetric (s-like) atomic-like potentials, V a , centered 9 2.2. SINGLE PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY at the atom positions, so that the Hamiltonian becomes, H = p 2m +V a (~ r) + X ~ Rn6=0 V a (~ r ~ R n ) = p 2m +V a (~ r) +V 0 (~ r): (2.11) Recall that by denition the atomic orbitals are eigenfunctions of the kinetic energy oper- ator and atom potential from the center atom with energy i . Combining this, Eq. 2.10 and Eq. 2.11 we get, H ij (k) = X ~ Rn e ik ~ Rn Z i (~ r) h p 2m +V a (~ r) +V 0 (~ r) i j (~ r ~ R n )d~ r = i ij + X ~ Rn e ik ~ Rn Z i (~ r)V 0 (~ r) j (~ r ~ R n )d~ r (2.12) The ~ R n = 0 term in the above sum is often neglected on account of the assumption that the atomic potential should be roughly constant in the region where the orbitals overlap. Thus the whole term (called the crystal eld) contributes a constant to the energy and therefore doesn't aect the electronic dispersion. The number of ~ R n 's to consider is treated as a tuneable parameter when building a tight-binding model. In general only nearest neighbor couplings are considered, and the sum is cuto for translation extending further than one neighboring atom. The value of the hopping term can be determined semi-empirically, whereby the strength is tuned in order for the resulting band structure to match experimental band structures or ones calculated from ab initio techniques. Otherwise parameters are chosen to give bandwidths on the proper order for the situation being investigated. The general framework laid out above can easily be extended for crystals with more than one atom per unit cell by simply constructing the Bloch sums from atomic orbitals centered on each atom in the unit cell according to, i (k;~ r) = 1 p N X ~ Rm e ik ~ Rm i (~ r ~ d ~ R m ); where i (~ r ~ d ~ R m ) is an atomic orbital with energy i centered on the atom at ~ d in the unit cell ~ R m . The model is then built in the same way as before by constructing (or calculating) hopping terms between atomic orbitals. Tight-binding models can also be used in situations where the system lacks translational invariance. In such systems the Hamiltonian is constructed over the entire real-space system. This creates an easy framework, with relatively low computational cost, to study eects of disorder such as in the Anderson model. 10 CHAPTER 2. GENERAL THEORY 2.2. SINGLE PARTICLE PROPERTIES 2.2.2 Density Functional Theory The tight-binding method discussed above presents a recipe for constructing a model Hamil- tonian to describe the electronic properties of a crystal under harsh simplications. Fur- thermore, the resulting Hamiltonian has to be solved by full diagonalization to obtain the eigenenergies and eigenstates. In the following we will discuss a method to nd (ap- proximately) the ground state of the fully interacting many-electron system using density functional theory (DFT), for which W. Kohn received the Noble Prize in Chemistry in 1998. We are concerned with nding the ground state of the electronic Hamiltonian given by Eq. 2.5. The strength of DFT comes from the discovery that the ground state electronic density uniquely determines the ground state of the system. We will now justify this claim by relaying the essential building blocks of the theory. In order to do this we rearrange the electronic Hamiltonian into an internal part, consisting of the electronic kinetic energy and electron-electron interaction, and external part consisting of everything else (in this case the lattice potential V nucl ), H elec =H int +V ext : (2.13) In order to match standard notation we also rewrite V ext = X i Z v ext (~ r)(~ r~ r i )d~ r where v ext V nucl (~ r): (2.14) Now we assumeH int is xed i.e. constant electron number,N, and the electronic mass, charge and interaction are as usual. This leaves the external potential, v ext , as the only variable in the system. We will now present and prove the Hohenberg-Kohn theorem [66] which forms a central pillar of DFT. Hohenberg-Kohn Theorem There is a one-to-one correspondence between the electronic ground-state density of an N electron system and the external potential acting on it. In order to prove this theorem we rst dene the ground-state density, n 0 (~ r), as n 0 (~ r) =h 0 j X i (~ r~ r i )j 0 i; (2.15) wherej 0 i =j 0 (~ r 1 ;~ r 2 ;:::;~ r N ; )i is the crystal's electronic ground state. Next, combining Eq. 2.14 and 2.15 we get, h 0 jV ext j 0 i = Z n 0 (~ r)v ext (~ r)d~ r: (2.16) 11 2.2. SINGLE PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY Now, lets assume that there are two dierent external potentials, V ext 6=V 0 ext , that lead to the same ground state electronic density, n 0 . The dierent external potentials would lead to dierent Hamiltonians,H6=H 0 , which would therefore have dierent groundstates, j 0 i6=j 0 0 i. Consider now, 0 0 H 0 0 = 0 0 (H int +V ext +V 0 ext V 0 ext ) 0 0 =E 0 0 + Z n 0 (~ r)[v ext (~ r)v 0 ext (~ r)]d~ r; (2.17) whereE 0 0 is the ground state energy ofH 0 and we used Eq. 2.16. The ground state energy of H is strictly lower than the expectation value of H with any other state so we know that, E 0 <E 0 0 + Z n 0 (~ r)[v ext (~ r)v 0 ext (~ r)]d~ r (2.18a) Using the same argument with the opposite quantities we get, E 0 0 <E 0 + Z n 0 (~ r)[v 0 ext (~ r)v ext (~ r)]d~ r: (2.18b) Combining the above two expressions we nd, E 0 <E 0 0 + Z n 0 (~ r)[v ext (~ r)v 0 ext (~ r)]d~ r<E 0 ; (2.19) which is obviously a contradiction, thereby proving that dierent potentials V ext and V 0 ext cannot result in the same ground state electron density, n 0 . Thus n 0 uniquely determines the potential, V ext . ./ We have now proven the existence of a functional relationship, v ext (~ r) =F [n 0 (~ r)]; that would allow us to nd the external potential (and consequently the Hamiltonian) once the ground state density is known. Given that v ext denes the Hamiltonian and therefore also the ground state we can write, v ext (~ r)) 0 [v ext ] 8 > < > : )E[v ext ] )T [v ext ] )V ee [v ext ] (2.20) Since there is a one-to-one correspondence between v ext and n 0 , the expectation values in Eq. 2.20 are also functionals of n 0 . This allows us to establish a variational principle 12 CHAPTER 2. GENERAL THEORY 2.2. SINGLE PARTICLE PROPERTIES concerning the ground state density. Given a system of interacting electrons in a specic external potential, we consider the functional E HK [n(~ r)] =h 0 [n(~ r)]jH elec j 0 [n(~ r)]i =h 0 [n(~ r)]jT +V ee +V ext j 0 [n(~ r)]i; (2.21) whereV ext is xed butn(~ r) is allowed to vary and 0 [n(~ r)] is the ground state of the system with density n(~ r). The absolute minimum of the functional is achieved when 0 [n(~ r)] is the ground state of H elec which occurs when n(~ r) = n 0 (~ r) i.e. the exact ground state electron density of the system. Thus, the Hohenberg-Kohn energy functional (Eq. 2.21), which exists and is unique, has its minimal value at the exact ground state energy of the many-body, interacting electron system. Kohn-Sham Equations We will now carry out the variational procedure and thereby derive the so-called Kohn- Sham equations. The breakthrough by Kohn and Sham was to realize that if you can construct a non- interacting Hamiltonian that has the same ground state density as the interacting system you can solve the interacting problem simply by solving the simpler eective non-interacting Hamiltonian [67]. This procedure generally works as long asn(~ r) is well behaved [68]. The Kohn-Sham Hamiltonian, given by H KS = ~ 2 2m r 2 +v KS (~ r); (2.22) is thus constructed to have a ground state density equal to n 0 . The problem has now been casted into constructing the proper Kohn-Sham potential, v KS . Firstly, we decompose the ground state density to a sum of N orthonormal orbital contributions, n 0 (~ r) = N X i i (~ r) i (~ r): (2.23) Next we extract the non-interacting (single particle) mean kinetic energy, T s [n 0 ], and the Hartree potential, V H [n 0 ], T s [n 0 ] = X i h i (~ r)j ~ 2 2m r 2 j i (~ r)i (2.24a) V H [n 0 ] = 1 2 Z Z n 0 (~ r) e 2 j~ r~ r 0 j n 0 (~ r 0 )d~ rd~ r 0 1 2 X i;j h i j j e 2 r ij j i j i: (2.24b) The non-interacting form of the Hohenberg-Kohn functional from Eq. 2.21 can then be written as E HK [n(~ r)] =T s [n(~ r)] +V H [n(~ r)] + Z n(~ r)v ext (~ r)d~ r +E xc [n(~ r)]; (2.25) 13 2.2. SINGLE PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY where the exchange correlation functional E xc [n] is dened as E xc [n(~ r)] =T [n(~ r)]T s [n(~ r)] +V ee [n(~ r)]V H [n(~ r)]: (2.26) Here it is important to note that T [n] andV ee [n] is evaluated with the ground state wave- function of the full interacting system (see Eq. 2.20). Rewriting Eq. 2.25 in a clearer fashion, E HK [n(~ r)] = N X i h i j ~ 2 2m r 2 +v ext j i i + 1 2 N X i;j h i j j e 2 r ij j i j i +E xc [n(~ r)]; (2.27) we see that according to the standard variational procedure we vary 1 ; 2 ;:::; N to nd the minimum of the functional under the constraint of orthonormalization of the orbital wavefunctions. The variational expression for the rst two terms of Eq. 2.27 leads to the well known Hartree-Fock equations, discussed in Appendix A.1. The variation of E xc [n(~ r)] is given by E xc [n] = Z V xc (~ r)n(~ r)d~ r; (2.28) where V xc (~ r) E xc [n] n(~ r) : Combining the Hartree-Fock equations with Eq. 2.22 and Eq. 2.28 we arrive at the Kohn-Sham equations ~ 2 2m r 2 +v KS (~ r) i (~ r) = i i (~ r); (2.29) where the i 's are Lagrange multipliers associated with the orbital energies i [69] although this has to be justied case by case. At this point we have determined the Kohn-Sham potential of the non-interacting system that leads to the same ground state density as the interacting system to be v KS (~ r) =v nucl (~ r) +v Hartree (~ r) +v xc (~ r); (2.30) where v nucl is the external or ionic potential and v Hartree comes from the Hartree-Fock equations. The exact ground state energy of the interacting system is given by E 0 = X i i 1 2 X ij h i j j e 2 r ij j i j i +E xc [n 0 ] Z V xc (~ r)n 0 (~ r)d~ r (2.31) If the exchange correlation functional was known we would be able to exactly solve the ground state energy of the interacting system, but unfortunately it is not known and we 14 CHAPTER 2. GENERAL THEORY 2.2. SINGLE PARTICLE PROPERTIES have to rely on appropriate approximations. Once an approximate exchange correlation functional is chosen the Kohn-Sham potential (Eq. 2.30) is calculated with a trial density and the Kohn-Sham equations (Eq. 2.29) are used to solve for the resulting wave functions which can then be used in Eq. 2.23 to calculate a new density. This procedure is repeated until self consistency is reached, at which pointn 0 (~ r) has been found (up to the errors from the approximation of E xc [n 0 ]). Local Density Approximation One of the most common approximation for the exchange correlation functional, E xc [n], is what is known as the local density approximation (LDA). Here it is assumed that the exchange correlation potential is purely local (valid for systems with slowly varying density n(~ r)) and given by: E LDA xc [n(~ r)] = Z n(~ r)" hom xc (n(~ r))d~ r: (2.32) The homogeneous electron gas exchange correlation energy density, " hom xc can be sepa- rated into an exchange and a correlation part " hom xc (r s ) =" hom x (r s ) +" hom c (r s ); where r s is a measure of the electron density dened by 1=n = (4=3)(r s a B ) 3 , with a B = ~ 2 =me 2 being the Bohr radius. Expressions for the two terms have been obtained by parameterizing results obtained from Monte Carlo calculations for the homogeneous elec- tron gas. Specically, Perdew and Zunger (1998) [70] reported the following expressions from calculations on the unpolarized electron gas: " hom x (r s ) = 3 4 3 2 2=3 e 2 r s a B (2.33a) " hom c (r s ) = ( 0:1423=(1 + 1:0529 p r s + 0:3334r s ) for r s 1 0:0480 + 0:0311 lnr s 0:0116r s + 0:0020r s lnr s for r s 1 (2.33b) where the energies expressed in " hom c (r s ) is given in Hartrees (1 Hartree = 2 Rydberg). Crystal Potential We already stressed that the only theoretical obstacle to DFT is in a proper approximation for the exchange correlation functional. However, a technical diculty in applying the theory comes in when dealing with the crystal potential. In order to demonstrate this we will attempt to solve the Kohn-Sham equations given by Eq. 2.29 using the orthogonalized plane wave (OPW) method suggested by Sommerfeld and Bethe (1933) [71]. 15 2.2. SINGLE PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY We start by expanding the ground state wavefunction (with energy E 0 ) of the Kohn- Sham Hamiltonian (Eq. 2.22) in terms of plane waves k (~ r) = 1 p X G c k+G e i(k+G)~ r : (2.34) We dene a new variable k i = k + G i and use Eq. 2.34 to evaluate Eq. 2.29, thereby ending up with the innite set of equations, labeled by i. X j [(k 2 i E 0 ) ij +v KS (k i k j )] = 0: (2.35) Obviously in practice we would only take a nite number of k i vectors into account, establishing a cut-o energy such that ~ 2 =(2m)k 2 i <E cut , and here is where we run into a problem. Recall from Eq. 2.30 that v KS includes the crystal potential v nucl and a simple argument shows that a pure expansion of crystal states into plane waves is impractical due to the so called variational collapse problem. Consider a crystal with lattice constant a and fcc structure therefore giving a unit volume of =a 3 =4. The core state 1s of an atom with atomic numberZ can be estimated to have a radius,a 1s a B =Z, so in order to properly describe it with plane waves, requires a cut-ok max 2=a 1s = 2Z=a B . The number of reciprocal lattice vectors that t inside a sphere of radiusk max is given by (4=3)k 3 max =n(2) 3 = which givesnZ 3 (a=a B ) 3 . Now, for example, a silicon crystal (Si has Z = 14) has lattice constant a = 5:43 A = 10:26a B , showing that we would need n 10 6 number of k i vectors in the plane wave expansion to reliable resolve the crystal potential. This in turn means in order to solve the Kohn-Sham Hamiltonian (Eq. 2.22) would require diagonalizing a matrix of dimension 10 6 10 6 , which is not feasible seeing as matrix diagonalization has complexity O(n 3 ). There have been several advances to overcome the variational collapse problem in- cluding the orthogonalization procedure suggested by Herring [72], the cellular method of Wigner and Seitz [73] and the augmented plane waves (APW) method of Slater [74]. It is this last method we will focus on. The APW method relies on the assumption that, to a reasonable level, the crystal potential is spherically symmetric within non-overlapping spheres centered on each atom in the crystal and constant elsewhere. This assumption is generally referred to as the mun- tin approximation. Within the non-overlapping spheres (of radiusr s ), spherical harmonics Y L (~ r) are used to construct a solution, R l (;r), to the radial Schr odinger equation @ 2 (rR l (;r)) @r 2 = V (r) + l(l + 1) r 2 (rR l (;r)): (2.36) An augmented plane wave A(k i ;~ r;) is dened as a plane wave outside the sphere r s and a spherically symmetric function inside. At the boundaryj~ rj = r s the two parts of the 16 CHAPTER 2. GENERAL THEORY 2.2. SINGLE PARTICLE PROPERTIES APW join continuously. Using the expression for a plane wave in spherical harmonics the augmented plane wave can be expressed as: A APW (k i ;~ r;) = 8 > > > > < > > > > : A APW in (k i ;~ r;) = 4 p X lm i l j l (k i r s ) R l (;r s ) Y lm (k i )Y lm (~ r)R l (;r) A out (k i ;~ r;) = e ik i ~ r p 4 p X lm i l j l (k i r s )Y lm (k i )Y lm (~ r) (2.37) where A APW in (A out ) represents the augmented plane wave inside (outside) the mun spheres, is the unit cell volume and j l is the spherical Bessel functions of order l. In the above the energy dependence on the boundary has been circumvented by shifting the whole energy scale so that the interstitial regions have zero energy. This decreases the variational freedom in the method. A better approach is to expand R l (;r s ) in a Taylor series around some reference energy r . The inclusion of the energy derivative _ R l restores the variational freedom giving the linearized augmented plane waves (LAPW) A LAPW (k i ;~ r;) = 8 > > < > > : A LAPW in (k i ;~ r;) = X lm [a lm (k)R l (;r s ) +b lm (k) _ R l (;r s )]Y lm (~ r) A out (k i ;~ r;) = e ik i ~ r p (2.38) where in and out have the same meanings as before and a lm and b lm must be calculated with standard diagonalization [75]. With both APW and LAPW the spherical harmonics expansion has to be cut-o for some angular momentum l and magnetic quantum number m. These cut-os are determined by the nature of the specic atoms in the crystal. This general idea of building atomic potentials from dierent functions smoothly con- nected has been taken further with the full potential linearized augmented plane wave (FLAPW) basis [76] that relaxes the constraint of spherically symmetric potentials inside the mun spheres. Subsequently people have also introduced the idea of using precalcu- lated pseudopotentials [77{79] to mimic the crystal potential contributions from various atoms. The precalculated functions remove the nodal structure of the core potentials mak- ing them computationally easier to work with, while preserving the important aspects that contribute to material properties. In conclusion, we have outlined density functional theory, showing its powerful ap- proach to nding the groundstate of many-electron systems. We'd like to stress that the theory is exact but due to a lack of knowledge of the exchange correlation functional and computational diculties in evaluating the crystal potential it is necessary to make several approximations. Due to the nature of these approximations DFT shows its best results for systems where electron correlation is weak such as in simple metals and semi-conductors. 17 2.2. SINGLE PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY 2.2.3 Wannier Functions The preceding discussion concerning DFT showed that we can calculate the electronic states of a crystal in terms of Bloch states which we denote by nk (~ r), with n being the band index. These states are highly delocalized (which is why the reciprocal space is so useful), but the purpose of this thesis is to study many-body physics which can be more readily formulated in real space using a localized basis. Furthermore, DFT calculations are carried out for a relatively small number of k-vectors while we often require high energy resolution to reliable calculate many-body eects (as is the case with the superconduct- ing critical temperature). The Wannier construction, we will now discuss, allows us to transform the delocalised ab initio results in reciprocal space to highly localized real-space quantities. This gives us proper parameters to use in tight binding models that retain material-realistic electronic properties calculated from ab initio, but are much simpler to work with since we can isolate and detangle a subspace of the full energy window (ac- tive bands around the Fermi-level for example) thereby decreasing the dimensionality of the Hilbert space in which our tight-binding models live. Furthermore, the sophisticated Fourier-like transforms involved in the Wannier construction allows us to interpolate from a coarse k-grid to a highly dense one and thereby increasing our energy resolution. The Wannier function, w n ~ R (~ r), centered on lattice site ~ R is given by w n ~ R (~ r) = V (2) 3 Z BZ " X m U (k) mn mk (~ r) # e ik ~ R d ~ k; (2.39) whereV is the unit cell volume,U (k) is a unitary matrix and the integral is performed over the 1st Brillouin zone. The WF w n ~ R (~ r) is localized around ~ R but due to a gauge freedom in the Bloch functions, U (k) is not unique. Specically, a transformation of the form ~ nk E =e in(k) j nk i; will not change any observables in the system as long as n (k) is real and periodic, but it will change U (k) . Marzari and Vanderbilt [60] therefore suggested using this freedom to pick a gauge in which the Wannier functions are maximally localized (hence maximally localized Wannier functions). To this end they dene the spread of a Wannier function as = X n [h r 2 i n h~ r n i 2 ] = X n [ w n ~ 0 (~ r) r 2 w n ~ 0 (~ r) j w n ~ 0 (~ r) ~ r w n ~ 0 (~ r) j 2 ]: (2.40) We can split into a gauge invariant term I and a term ~ that is dependent on the gauge choice U (k) , as follows: I = X n [ w n ~ 0 (~ r) r 2 w n ~ 0 (~ r) X m ~ R j w m ~ R (~ r) ~ r w n ~ 0 (~ r) j 2 ] (2.41a) 18 CHAPTER 2. GENERAL THEORY 2.2. SINGLE PARTICLE PROPERTIES ~ = X n X m ~ R6=n ~ 0 j w m ~ R (~ r) ~ r w n ~ 0 (~ r) j 2 : (2.41b) Interestingly, the value of I coming from the valence bands of an insulator can be exper- imentally measured and is known as the mean-square quantum uctuation of the macro- scopic polarization [80]. Further value of this decomposition will become clear shortly. The details of calculating in reciprocal space are discussed in depth in Ref. [60] but due to the length of the derivations we will only quote the nal results, namely, I = 1 N X k;b w b J X mn jM k;b mn j 2 ! ; (2.42a) and ~ = 1 N X k;b w b X m6=n jM k;b mn j 2 + 1 N X k;b w b X n ( Im lnM k;b mn bh~ r n i); (2.42b) whereN is the number of k-vectors considered,J is the number of Bloch states counted, b is a vector connecting k to its neighboring points with w b being the corresponding weight needed to calculater k (see Appendix B of Ref. [60]) and h~ r n i = 1 N X k;b w b b Im lnM k;b mn : (2.43) The matrixM k;b measures overlap between the periodic components of Bloch states from neighboring k's M k;b mn =hu mk ju nk+b i: (2.44) Now suppose we would like to nd a localized, tight binding description of theN b bands surrounding the Fermi energy (E F ) of some crystal. The strategy to do this was laid out in Ref. [61], and works roughly as follows. Firstly, we dene an energy window in which the N b bands we are interested in falls. For every k-vector there will be a number N k of bands that fall inside the energy window. By construction we will have N k N b . This gives us an N k -dimensional Hilbert spaceF(k) spanned by the states within the window. For k-vectors where N k > N b we need to nd the N b -dimensional subspaceS(k)F(k) that leads to the smallest I (Eq. 2.42a). In order to understand why we need to minimize I we can rewrite Eq. 2.42a as I = 1 N X k;b w b T k;b ; (2.45) with T k;b =N b X mn jM (k;b) mn j 2 = Tr h ^ P k ^ Q k+b i ; (2.46) 19 2.2. SINGLE PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY where ^ P k = P n ju nk ihu nk j is the projector ontoS(k) and ^ Q k = 1 ^ P k . The quantityT k;b is called the `spillage' from spaceS(k) toS(k + b) since it measures mismatch between them. When T k;b is zero the two spaces are the same and I vanishes. In order to nd the minimum of I a variational procedure is used which leads to the following eigenvalue problem that has to be solved self-consistently: X b w b [ ^ P (i) k+b ] in ! u (i) mk E = (i) mk u (i) mk E ; (2.47) where [ ^ P (i) k+b ] in = ^ P (i1) k+b + (1)[ ^ P (i1) k+b ] in and is a hyperparameter tuned to achieve the fastest convergence (normally setting 0:5). Once a convergent solution is reached and the optimal subspaceS(k) is known we can diagonalize the Hamiltonian in the subspace to obtainN k eigenvalues (~ " nk ) and Bloch-like eigenfunctions ( ~ nk = e ik~ r ~ u nk ). The next step is to nd U (k) that can be used in Eq. 2.39 to construct the maximally localized Wannier functions. This is done (as explained in Ref. [60]) by iteratively updating the overlap matrix M k;b (Eq. 2.44) using U (k) and then updating U (k) with M k;b and so forth until convergence is reached. Tight Binding Hopping Parameters Finally, with the knowledge ofU (k) we can calculate the Hamiltonian matrix for the rotated states, H rot (k) = (U (k) ) y ~ H(k)U (k) ; (2.48) where ~ H mn (k) = ~ " mk m;n is the diagonalized Hamiltonian in the subspaceS. Using a standard Fourier transform we can obtain the real-space Hamiltonian for the subspace, H rot ( ~ R) = 1 N X k e ik ~ R H rot (k) =t ( ~ R); (2.49) where we use the subscripts and to match standard notation for tight-binding models with hopping parameters between chemical orbitals. With the knowledge of the hopping parameters we can construct H rot (k) at arbitrary k 0 -points using H rot (k 0 ) = X ~ R e ik 0 ~ R t ( ~ R): (2.50) Diagonalizing this matrix yields the band structure (eigenenergies) and eigenstates in the so-called band basis, H rot (k)jnki =" nk jnki; (2.51) as well as the unitary to transform an arbitrary operator from the orbital basis to the band basis. 20 CHAPTER 2. GENERAL THEORY 2.3. MANY PARTICLE PROPERTIES 2.3 Many Particle Properties The subsequent discussion of many-body eects follows closely from the solid state physics textbook by Grosso and Parravicini [57] chapters VII - IX (regarding plasmonic proper- ties) as well as the textbook by Fetter and Walecka [81] chapter 12 (regarding phononic properties). The lecture notes by Heid [82] were used in the discussion of DFPT. The corresponding chapters in the PhD thesis of Dr. R osner [63] were also closely followed. In the preceding section we focused on the energy eigenstates of a crystal based on the one-electron approximation with the lattice ions xed at their equilibrium positions. In this section we relax both these approximations and discuss the resulting many-body eects. These eects are manifested as global excitations (or collective excitations) that can in turn be regarded as a non-interacting (or very slightly interacting) system of pseudo-particles. Specically, in the case of a collective excitation of many electrons we call these pseudo- particles plasmons, described in section 2.3.1, while the interaction of electrons with the lattice ions create excitations called phonons, as discussed in section 2.3.2. 2.3.1 Plasmons Plasmons form a central theme of this thesis. We discuss their theoretical aspects here. Classical Derivation of Plasmon Dispersion Here we present a classical derivation of the plasmon dispersion valid for systems of any dimensionality, following the method showed by Santoyo and Castillo-Mussot [83]. We consider a system where electrons eectively move freely between a xed background of positive charge so that the system is neutral altogether. We use a hydrodynamic model to describe the electrons (meaning the 'sea' of electrons is seen as a uid) with a linearized Navier-Stokes equation of motion, @ ~ j(~ r;t) @t = n 0 e 2 m ~ E(~ r;t) 2 r(~ r;t); (2.52) where ~ j is the current density, is the induced density of charge (due to electrons moving) and is a parameter like a spring constant. In the above equation we attribute electron motion to both electrical elds ( ~ E) and \pressure" due to a density gradient. The total electric eld can be found in the normal way (assuming non-relativistic speeds) using ~ E =r tot , where, tot (~ r;t) = ext (~ r;t) + Z (~ r 0 ;t) j~ r~ r 0 j dr 0 : (2.53) Here we note that the total charge density is split into 2 terms, (~ r;t) = 0 +(~ r;t). This allows us to relate the induced charge density to the current density using the continuity 21 2.3. MANY PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY equation @ @t +r ~ j = 0: We can solve for normal modes of the electron density oscillations (plasmons) by dier- entiating the continuity equation above with respect to time (we can swap the order of dierentiation on ~ j since time and space coordinates are independent). This gives @ 2 @t 2 +r @ ~ j @t = 0: Next we take the divergence of Eq. 2.52 and substitute the above into it, to get ( @ 2 @t 2 + 2 r 2 )(~ r;t) = n 0 e 2 m r ~ E(~ r;t) = n 0 e 2 m r 2 Z (~ r 0 ;t) j~ r~ r 0 j dr 0 ; (2.54) where we let ext = 0. We then Fourier transform to momentum space according to (~ r;t) = X q;! (~ q;!)e i!t e i~ q~ r ; (2.55) (~ r;t) = X q;! (q)(~ q;!)e i!t e i~ q~ r : (2.56) Here we use(q) as a placeholder for the Fourier transform ofj~ r~ r 0 j 1 . Ind> 1 dimensions (q)/ q 1d . Specically for 2 and 3 dimensional systems (which are of course of most interest in condensed matter) we would have, (q) = 2 j~ qj 2d; (q) = 4 j~ qj 2 3d: (2.57) Plugging the Eq. 2.55 and Eq. 2.56 into Eq. 2.54 gives ! 2 2 q 2 n 0 e 2 m q 2 (q) (~ q;!) = 0; from where we can read o the normal plasmons modes as: (assuming is very small) ! 2 = n 0 e 2 m (q) + 2 q 2 n 0 e 2 m (q) q 2 (2.58) 22 CHAPTER 2. GENERAL THEORY 2.3. MANY PARTICLE PROPERTIES Thus in d dimensions we have !(q)/q (3d)=2 . Plugging in (q) from Eq. 2.57 we get the plasmon dispersions for these cases: !(q) = 2n 0 e 2 m p q 2d; !(q) = 4n 0 e 2 m 3d: (2.59) This simple derivation shows the well known 2d acoustic plasmon dispersion with !(q) p q. We also see 3d plasmons are dispersionless (constant). The dierence in plasmonic behavior in 2d and 3d system is a central topic of this thesis which we will come back to often. The above derivation is useful for understanding the emergence of \stable" collective oscillations of the Fermi sea, but in order to calculate material realistic plasmonic proper- ties we need a quantum mechanical description of this many-body eect. In linear response theory we assume the induction electric eld is related to the total electric eld by ~ D =" ~ E; (2.60) where we dene a function " called the dielectric function. We can frame this statement in terms of electric potential functions, in order to make things clearer, by considering the eect of an external, perturbing potential on the ground state electron density, 0 (~ r), of a system. If we apply a perturbation ext (~ r;t) the electron density will be rearranged to screen the perturbation thereby forming an induced charged density ind (~ r;t) which we can relate to an induced potential ind (~ r;t) through the usual relation ind (~ r;t) = Z ind (~ r 0 ;t) j~ r~ r 0 j d~ r 0 : (2.61) This in turn allows us to express the total perturbing potential (~ r;t) = ext (~ r;t) + ind (~ r;t): (2.62) Equivalently to Eq. 2.60 we can now dene the dielectric function"(~ r;~ r 0 ;t;t 0 ) as a function that ensures (~ r;t) = Z d~ r 0 Z dt 0 ext (~ r;t) "(~ r;~ r 0 ;t;t 0 ) : (2.63) Note that the external perturbation does not depend on any material specic properties, therefore all such properties are captured by the dielectric function. If we are dealing with a crystal that has translational invariance and we assume the external potential to 23 2.3. MANY PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY be periodic in time, the dielectric function would only depend on spatial and temporal dierences, "(~ r~ r 0 ;tt 0 ). In this situation it is convenient to take the Fourier transform of Eq. 2.63 to get (q;t) = ext (q;!) "(q;!) : (2.64) Quantum Mechanical expression of the Dielectric function Consider a homogeneous sample of some crystal that can be described with the one-electron Hamiltonian, H 0 , with eigenstates nk and corresponding energies E nk . Now suppose we perturb the system with ext (~ r;t) =A 0 e i(q~ r!t) +c:c: (2.65) where A 0 is the amplitude of the perturbation and c:c: refers to the complex conjugate as usual. Seeing as the perturbation is time dependent it would induce transitions between states that dier in energy by ~!, according to Fermi's golden rule , at a rate R(q;!) = 4 ~ X nk;mk 0 jh nk 0je iq~ r j mk ij 2 (E nk 0E mk ~!)[f(E mk )f(E nk 0)] (2.66) where we add a factor of 2 for spin degeneracy and we make use of the Fermi functions f(E) = 1=(e (E)) + 1) to account for the density of transition states, with being the chemical potential and the inverse temperature. At this point we make the following approximation, known as the random phase ap- proximation (RPA) that will lead us to the Lindhard dielectric function, X nk;mk 0 h nk 0je iq~ r j mk i X nmk h nk+q je iq~ r j mk i: (2.67) Now, the electric eld from the perturbation is given by e ~ E =r(~ r;t) (where e here is the electron charge) i.e. ~ E(~ r;t) =E 0 ~ ee i(q~ r) +c:c: (2.68) where ~ e is the versor of q and E 0 = iqA 0 =e. With such an electric eld acting on an isotropic medium we can reasonably assume an induced current parallel to the eld and proportional to it (this is called linear response), with the same time dependence, therefore the current density, ~ J can be written ~ J(~ r;t) =(q;!)E 0 ~ ee i(q~ r) +c:c: (2.69) Fermi's golden rule states that a perturbation, H 0 , oscillating in time with frequency ! will transition a system from the initial statejii to the set of statesjfi with energy dierence from the initial state of ~! at a rate, i!f = 2=~jhfjH 0 jiij 2 , where is the density of nal states. [84] 24 CHAPTER 2. GENERAL THEORY 2.3. MANY PARTICLE PROPERTIES where (q;!) is the longitudinal conductivity function. We can calculate the power dissi- pated by this current (in a crystal with volumeV), Z V ~ J ~ Ed~ r = Z V [(q;!)E 0 e i(q~ r) +c:c:][E 0 e i(q~ r) +c:c:]d~ r =(q;!)jE 0 j 2 V +c:c: = 2 1 (q;!) q 2 e 2 jA 0 j 2 V (2.70) where 1 indicates the real component of . Eq. 2.70 gives us a classical expression for the power dissipated in the system, while Eq. 2.66 expresses the quantum mechanical rate of energy loss. By simple energy conservation we get 2 1 (q;!) q 2 e 2 jA 0 j 2 V = ~!R(q;!) =) 1 (q;!) = e 2 ~! 2q 2 V R(q;!) jA 0 j 2 (2.71) By denition ~ J = ~ E. Using this along with Eq. 2.60 and the requirement that @ ~ D @t = @ ~ E @t + 4 ~ J; which establishes the relationship " = 1 + 4i=!. This gives us the imaginary part of the dielectric function, " 2 (q;!), as " 2 (q;!) = 4 ! 1 (q;!) = 2~e 2 q 2 V R jA 0 j 2 : (2.72) where we plugged in Eq. 2.71. The real part of the dielectric function, " 1 (q;!), can now be obtained using the Kramers-Kronig relation " 1 (q;!) = 1 + 1 P Z 1 1 " 2 (q;!) ! 0 ! d! 0 : (2.73) Plugging Eq. 2.66 (using Eq. 2.67) into Eq. 2.72 and the result into Eq. 2.73 gives us " 1 (q;!) = 1 + 8 2 e 2 q 2 V X nmk f(E mk )f(E nk+q ) E nk+q E mk ~! jh nk+q je iq~ r j mk ij 2 : (2.74) Finally, using the formal relation lim !0 + 1 xi =P 1 x +i(x) 25 2.3. MANY PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY we can combine Eq. 2.74 and Eq. 2.72 into a single expression: "(q;!) = 1 + 8 2 e 2 q 2 V X nmk f(E mk )f(E nk+q ) E nk+q E mk ~!i jh nk+q je iq~ r j mk ij 2 : (2.75) where plays the role of a small broadening factor, which we will usually denote by 0 + . See Ref. [85] for an in depth discussion of it's physical signicance. Lets return to the total power dissipation, P , expressed in Eq. 2.70 for a moment. From Eq. 2.72 we get 1 (q;!) = ! 4 " 2 (q;!) (2.76) Plugging this back into Eq. 2.70 and recalling that A 0 =A ext =" we get P (q;!) =V !q 2 2e 2 " 2 (q;!) j"(q;!)j 2 jA ext (q;!)j 2 : (2.77) which shows that the amount of energy dissipated from the system due to a driving per- turbation of amplitude A ext is governed by the imaginary part of the dielectric function. We therefore dene the energy-loss function: Im 1 "(q;!) = " 2 (q;!) " 2 1 (q;!) +" 2 2 (q;!) ; (2.78) which is a quantity that can be measured in electron energy-loss spectroscopy (EELS) experiments. The maximum of the energy-loss function denes the plasmonic modes (res- onances of the collective oscillation of electrons). The pairs (!; q) at which these maxima occur dene the plasmonic dispersion. 2.3.2 Phonons Firstly, we'll show the existence of harmonic excitations of the ionic background in a crystal (which are called phonons) and then we'll discuss density functional perturbation theory as a method of calculating material specic phononic properties. In Section 2.2 we discussed methods of nding the ground-state electronic density under the approximation that the lattice ions are xed and therefore their positions were treated as parameters of the electronic Hamiltonian Eq. 2.5. We will now explore what happens when the ions are allowed to move. For a start we'll only be concerned with longitudinal movement since this gives rise to changes in the charge density which aects the Coulomb interaction between electrons and nuclei. We'll also, for simplicity of derivation, assume the nuclei background to be a homogeneous medium with no shear strength (we'll relax this later). 26 CHAPTER 2. GENERAL THEORY 2.3. MANY PARTICLE PROPERTIES We denote by ~ d(~ r) a small displacement of the point ~ r from it's equilibrium position. The change in volume is to rst order in derivatives dV 0 =dV (1 +r ~ d) and the corresponding change in density is 0 = n n 0 =r ~ d; (2.79) where 0 (~ r) = Mn 0 (~ r) is the equilibrium mass density (M is the mass and n 0 (~ r) the density). The medium will support waves (with speed v) that satisfy the wave equation r 1 v 2 @ 2 ~ d @t 2 r 2 ~ d ! = 0: (2.80) Now, in a medium without shear strength and vorticity, ~ d satises the general condition I ~ dd ~ l = 0; for any closed path ~ l, which in turn means that r ~ d = 0: Hence, r 1 v 2 @ 2 ~ d @t 2 r 2 ~ d ! = 0: (2.81) The combination of Eq. 2.80 and Eq. 2.81 tells us the sound waves can be described by the pair of equations 1 v 2 @ 2 ~ d @t 2 r 2 ~ d = 0: (2.82a) r ~ d = 0: (2.82b) Next, we construct a Lagrangian (L 0 =TV ) that results in the above equations of motion, L 0 = 1 2 Z 0 @ 0 X i @d i @t @d i @t B X i;j @d i @r j @d i @r j 1 A d~ r; (2.83) where B = 0 v 2 . We can derive a Hamiltonian from this Lagrangian H 0 = 1 2 Z p 2 0 +B(r ~ d) 2 d~ r; (2.84) 27 2.3. MANY PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY where Eq. 2.82b was used in Eq. 2.83 and the canonical momentum is given by p(~ x;t) = 0 @ ~ d(~ r) @t : Assuming a large box with volumeV and periodic boundary conditions we introduce the normal mode expansion to get p(~ x;t) =(Mn 0 ) 1=2 X k ~! k 2V 1=2 k k (b k e i(k~ r! k t) +b y k e i(k~ r! k t) ) (2.85a) ~ d(~ x;t) = i (Mn 0 ) 1=2 X k ~! k 2V 1=2 k k (b k e i(k~ r! k t) +b y k e i(k~ r! k t) ) (2.85b) where ! k =vk. Plugging these expressions back into Eq. 2.84 we get H 0 = 1 2 X k ~! k (b y k b k +b k b y k ); (2.86) which is exactly the Hamiltonian of a system of non-interacting harmonic oscillators with b k (b y k ) being the annihilation (creation) operator for the pseudo particle (phonon) with momentum k associated with lattice vibration excitations. This shows that phonons obey bosonic statistics, as all harmonic oscillations do. Of course the above is a very simplied approach to the nuclei problem and in order to describe the system properly we have to relax the assumptions we made earlier concerning homogeneousness. In reality the ions are acted upon by the Coulomb interaction between the ions and the electronic density. The large ionic masses do, however, allow us to simplify the problem somewhat by making it possible to solve the problem in two distinct steps. Firstly, the small amplitude of the ionic vibrations means we can consider the ions to be xed in their equilibrium position when calculating the ground state electronic density, as already described. Then in a second step we calculate the change in energy associated with a sequence of stationary ionic congurations displaced slightly from their equilibrium positions. In each step of the sequence we assume the electrons to have the ground-state density associated with the instantaneous ionic conguration. We can do this since the lattice vibrations occur on a much slower scale than the electronic movement. This is called the Born-Oppenheimer approximation [86]. Now, consider a crystal with a composite lattice with nucleic equilibrium position ~ R 0 mi = ~ R m + ~ d i of atomi at ~ d i within unit cellm at ~ R m . For simplicity we'll use the single notation I =fm;ig. We introduce a small displacement ~ u I to each atom such that ~ R I = ~ R 0 I +~ u I : (2.87) 28 CHAPTER 2. GENERAL THEORY 2.3. MANY PARTICLE PROPERTIES This results in a change of the electronic ground state energyE 0 ( ~ R 0 ), which we can expand around~ u I = 0 to get E 0 ( ~ R) =E 0 ( ~ R 0 ) + X I @E 0 ( ~ R) @R I ~ R 0 u I + 1 2 X IJ @ 2 E 0 ( ~ R) @R I @R J ~ R 0 u I u J +::: (2.88) where the indices and label dimensions. The derivatives of an energy with respect to space gives a force in the direction of the derivative, hence, we can write F I = @E 0 ( ~ R) @u I ~ R 0 ; (2.89) but this has to vanish otherwise ~ R 0 would not be the equilibrium conguration. Similarly, we can nd the net force acting on atom J =fn;jg in the direction when atoms I are displaced by u I : F J = X I @ 2 E 0 ( ~ R) @u I @u J ~ R 0 u I =F nj; = X m;i; ij (n;m)u m;i; (2.90) The quantity ij (n;m) can clearly be associated with an harmonic force constant similar to a spring constant between two atoms. In a translationally invariant system it depends only on the distance between atoms, ~ R = ~ R n ~ R m , and as usual it becomes convenient to work in reciprocal space. The Fourier transform gives, ij (q) = 1 N X ~ R ij ( ~ R)e iq ~ R ; (2.91) which is used to dene the dynamical matrix as: D ij (q) = 1 p M i M j ij (q): (2.92) The eigenvectors of the dynamical matrix are called polarization vectors e (q) of the phonon mode , and the eigenvalues dene the phononic dispersion. Note that in d- dimensions we get d-eigenmodes for every (i;j) atom-pair (atoms in unit cell). We see the only necessary ingredient to characterize the phonons of a given crystal is the dynamical matrix (within the harmonic approximation), which can be calculated from the response of the ground-state energy. The ground state energy in turn depends on the ground state density, which we can calculate from DFT. Thus by applying a perturbation to the ground-state density we can calculate the dynamical matrix through what is called density functional perturbation theory (DFPT). Below we brie y discuss DFPT, which is used in this thesis to calculate phononic properties. 29 2.3. MANY PARTICLE PROPERTIES CHAPTER 2. GENERAL THEORY Density Functional Perturbation Theory The DFT scheme discussed in Section 2.2.2 provided us with the ground state density of a non-interacting many electron system described by the Kohn-Sham Hamiltonian, Eq. 2.22. A slight perturbation of the atomic positions changes the eective potential in the Kohn-Sham Hamiltonian. Referencing Eq. 2.30 we see the variation is given by v KS (~ r) =v nucl (~ r) + Z v Hartree (~ r) n(~ r 0 ) + v xc (~ r) n(~ r 0 ) n(~ r 0 )d~ r: (2.93) This leads to a rst-order variation of the wavefunctions i (~ r) = X j6=i h j jv KS j i i i j j (~ r): (2.94) In the same way we can nd j , and using the denition of n(~ r) from Eq. 2.23 we get n(~ r) = X i f i [ i (~ r) i (~ r) + i (~ r) i (~ r)] = X i6=j f i f j i j h j jv KS j i i i (~ r) j (~ r) = 2 X i2v;j2c 1 i j h j jv KS j i i i (~ r) j (~ r); (2.95) where in the last step we used time-reversal symmetry of the system with v labeling the valence (occupied) states and c the conduction (unoccupied) states. We can further avoid summing over conduction states by writing n(~ r) = 2 X i2v i (~ r) i (~ r); (2.96) with i (~ r) = X j2c 1 i j h j jv KS j i i j (~ r): (2.97) The key is now to notice that (H KS i ) i = X j2c j j ih j jv KS j i i = ^ P c v KS j i i = ( ^ P v 1)v KS j i i: (2.98) where ^ P c is the projector onto the conduction states and ^ P v = 1 ^ P c is it's adjoint i.e. the projector onto the valence states. In this way only valence states need to be taken into We use the Hellman-Feynman Theorem [87] that states for a Hamiltonian that depends on a parameter we have dE d =h j dH d j i, where H j i =E j i. 30 CHAPTER 2. GENERAL THEORY 2.3. MANY PARTICLE PROPERTIES account in order to calculate n(~ r). In the specic case where we are dealing with a periodic lattice v nucl can be expressed as a linear combination of atomic potentials centered at the instantaneous position of the ions (again using the notation as in Eq. 2.87) v nucl (~ r) = X mi v i (~ r ~ R I ): (2.99) The rst order variation evaluated at the equilibrium positions is then, q i v nucl (~ r) =e iqr X mi e iq( ~ R 0 I ~ r) r ~ r (~ r ~ R 0 I ): (2.100) The Fourier transform of Eq. 2.96 has the form, q i n(q + G) = 4 V X k hkje i(q+G)~ r q i (k) (2.101) where Eq. 2.97 now becomes, q i (k) = X 2c hk + qj q i v KS jki (k + q) (k) jk + qi: (2.102) And the same as before we can solve for q i (k) and in turn q i n(q + G) by solving the linear equations (see Eq. 2.98): (H k+q KS (k) q i (k) = ( ^ P k+q v 1) q i v nucl jki: (2.103) Eq. 2.93 and Eq. 2.101 form a set of equations that can be solved self-consistently to determine the variation of the electronic density due to a perturbation of the atomic posi- tions for a given q. Once self-consistency is reached the dynamical matrix elements from Eq. 2.91 are given by ij (q) = X G [ q i n(q + G) q j v nucl (q + G) + q i q j v nucl (G)]: (2.104) One very powerful result here is that the phonons are completely dened for independent q-vectors which allows us to calculate the dynamical matrix on a sparse grid of q-vectors and afterwards interpolate the result to a dense grid using a Fourier transform. 31 2.4. PARTICLES INTERACTING CHAPTER 2. GENERAL THEORY 2.4 Particles Interacting The subsequent discussion of particle propagators and Dyson's equation was adapted from the textbook by Fetter and Walecka [81] chapter 3 and Ref. [88]. The discussion of the GW approximation comes mostly from Ref. [89]. The electron-phonon interaction and spec- tral function discussions relied heavily on the book by H. Bruus and K. Flensberg [90] and Ref. [91]. The corresponding chapters in the PhD thesis of Dr. R osner [63] were also closely followed. So far in this thesis we have developed the framework necessary to study the elec- tronic properties of a crystal by constructing an eective one-electron Hamiltonian that allows us to treat the many electrons in the system as independent (non-interacting). We then proceeded to explore the collective excitations that emerge from these many-particle systems and showed that these excitations can in turn be regarded as non-interacting pseudo-particles. In this section we will present the theoretical framework necessary to study the eects brought on by interactions between the electrons in the system (through the Coulomb interaction) as well as the interaction of electrons and the pseudo-particles (plasmons and phonons) discussed earlier. We'll start by discussing particle propagators and the Dyson-equation which provides a general treatment of the interaction of parti- cles (and the consequent screening of interaction). Then we will consider electron-electron interaction, followed by electron-plasmon and nally, electron-phonon interaction. 2.4.1 Particle Propagators and Dyson's Equation Here we present a cursory overview of the Green's function formalism of many-body theory but the curious reader should turn to Ref's [81] or [92] for a proper discussion. Naturally it is impractical to fully account for the interactions between all particles in a system containing billions of particles. Even more so if you consider the fact that the interaction between any two particles aects the interaction of those particles with all the other particles they interact with, which in turn aects all the interactions of those particles, etc. A much better approach is to study the system with the use of propagators or Green's functions. Consider a system of non-interacting fermions given by the eld (~ r) = X k; k (~ r)c k; ; where k is a single particle wave-function and c k; is the fermion operator. We perform 32 CHAPTER 2. GENERAL THEORY 2.4. PARTICLES INTERACTING the canonical transformation to particles and holes by redening the fermion operator as: c k; = ( a k; k>k F particles b y k; k<k F holes (2.105) We can now rewrite the eld as S (~ r) = X k>k F k (~ r)a k; + X k<k F k (~ r)b y k; I (~ r;t) = X k>k F k (~ r)e i! k t a k; + X k<k F k (~ r)e i! k t b y k; (2.106) where we show the Shr odinger picture and the interaction picture. The Hamiltonian of the system is given by H 0 = X k k c y k; c k; = X k>k F k a y k; a k; X k<k F k b y k; b k; + X k>k F k (2.107) where the rst (second) term comes from extra particles (holes) in the system and the third term from the lled Fermi sea. If we x the number of fermions in the system all particles and holes come in pairs meaning that the lled Fermi sea is the ground state, j 0 i. Therefore we have b k; j 0 i =a k; j 0 i = 0 (2.108) The non-interacting fermionic Green's function, or single particle propagator, is dened as: iG 0 (~ rt;~ r 0 t 0 )h 0 j ^ T [ I (~ rt) y I (~ r 0 t 0 )]j 0 i (2.109) where ^ T is the time ordering operator and we assumej 0 i to be normalized. If we plug Eq. 2.106 into Eq. 2.109 using Eq. 2.108 we get with both possible time orderings iG 0 (~ rt;~ r 0 t 0 ) = (2) 3 Z e ik(~ r~ r 0 ) e i k (tt 0 ) [(tt 0 )(kk F )(t 0 t)(k F k)]dk (2.110) Taking the Fourier transform of Eq. 2.110, we get the single particle propagator in mo- mentum space G 0 (k;!) = (kk F ) ! k +i + (k F k) ! k i (2.111) The operator ^ T arranges operators so that earlier time operator appears rst, explicitly for two oper- ators A and B we have ^ T [A(t1)B(T2)] =(t1t2)A(t1)B(t2)(t2t1)B(t2)A(t1), where + is used for bosons and for fermions. 33 2.4. PARTICLES INTERACTING CHAPTER 2. GENERAL THEORY t 1 t 2 t 1 t 2 Figure 2.1: Feynman diagram for screened propagator (left), free propagator (middle) and Hartree-Fock interaction (right) of an electron entering the system at timet 1 and emerging at time t 2 . where we added thei terms in order to ensure convergence of the time integrals. We immediately notice that the excitation spectrum of the system is given as poles of the single particle propagator but it can also be used to calculate the expectation value of any single particle operator on the ground state of the system. The physical interpretation of the propagator is clear from it's denition (Eq. 2.109). We insert a particle (or hole) into the system in the groundstate at (~ r 0 ;t 0 ) and measure the response of the system by calculating the overlap of the resulting state at (~ r;t) with the ground state. Or equivalently we can phrase it as: we insert a particle (or hole) into the system at (~ r 0 ;t 0 ) and nd the probability amplitude of the particle emerging from the system at (~ r;t) [88]. We can similarly consider two-particle propagators where we introduce two particles (or holes) into the system at (~ r 0 1 t 0 1 ;~ r 0 2 t 0 2 ) and consider the probability amplitude of the particles emerging at (~ r 1 t 1 ;~ r 2 t 2 ). With such a propagator in hand we can calculate all two-particle operator expectation values over the ground state. The concept of propagators makes the many-body problem much simpler to solve (ap- proximately) since it aords us a systematic way of accounting for each 'type' of interaction in the system. In his book, Richard D. Mattuck [88], explains very clearly how Feynman diagrams can be used to keep track of the interactions in a system. The interested reader is referred there for more information. Consider the simplest type of interaction an electron propagating through a system could experience, the so-called Hartree-Fock interaction [93], where the electron interacts with the mean-eld caused by the other charges in the system. We can set up a dictionary of diagrams to represent the screened propagator, G, the free propagator, G 0 , and the Hartree-Fock interaction as shown in Fig. 2.1. If we now consider the 'trajectory' of an extra particle added to the system, it would propagate freely for some time until encountering an interaction, which would change it's 34 CHAPTER 2. GENERAL THEORY 2.4. PARTICLES INTERACTING 'trajectory' until it encounters another interaction, and so forth. If the Hartree-Fock inter- action (denoted here HF ) is the only type of interaction in the system, the 'full trajectory' (screened propagator) could be expanded by repeated application of this interaction, = + + + + ... (2.112) We can write this in algebraic form as G =G 0 +G 0 HF G 0 +G 0 HF G 0 HF G 0 +G 0 HF G 0 HF G 0 HF G 0 +::: =G 0 +G 0 HF G (2.113) which allows us to solve for the screened propagator G = G 0 1G 0 HF = 1 (G 0 ) 1 HF : (2.114) If we now use Eq. 2.111 and assume we inserted an extra electron into the system we get G (k;!) = ! k HF +i : (2.115) We therefore see that the interaction in the system shifts the energy spectrum by an amount Re( HF ). For this reason the term is referred to as the self energy, since it is an internal renormalization of the energy. The imaginary contribution of the self energy broadens the spectrum peaks and leads to a nite lifetime for excitations in the system (discussed more in Section 2.4.4). In order to calculate a proper screened propagator we also need to include interactions other than the Hartree-Fock interaction. The diagrammatic technique illustrated above gives us a systematic way to include higher order scattering events which leads to a more accurate self energy and consequently also a more accurate screened propagator. We will now explore the specic interactions we can account for coming from electron- electron scattering due to the Coulomb interaction and due to phonon interactions. 2.4.2 Coulomb Interaction The Hartree-Fock approximation we discussed above is the simplest self-energy one can construct for the electron-electron interaction. Sadly, though, it is not sucient to explain 35 2.4. PARTICLES INTERACTING CHAPTER 2. GENERAL THEORY many eects seen experimentally coming from the Coulomb interaction [94]. We therefore need to include higher order interactions, which we do in form of the random phase ap- proximation (discussed earlier with Eq. 2.67), which gives decent results for metals and most semi-conductors [95]. In terms of Feynman diagrams the RPA-interaction (shown in Eq. 2.116) consist of an innite series of particle-hole excitations screening the Coulomb interaction, Π RPA = + + + ... (2.116) If we now consider the screened Coulomb interaction, we can write a similar Dyson ex- pansion for it as we did in Eq. 2.112 for the Hartree-Fock interaction, where now the polarization is used, W = V + V W Π (2.117) The Feynman diagram in Eq. 2.117 shows the screened Coulomb interaction as W (q;!) =V (q) +V (q)(q;!)W (q;!) W (q;!) = V (q) 1V (q)(q;!) V (q) "(q;!) (2.118) which gives us the expression, "(q;!) = 1V (q)(q;!) (2.119) This expresses the same dielectric function as in Eq. 2.75. Considering that the bare Coulomb interaction is given by V (q) = 4 2 V e 2 q 2 ; where e is the electron charge andV the volume of the unit cell, we are left with an expression for the polarization function within the RPA approximation : RPA (q;!) = 2 X nmk f(E mk )f(E nk+q ) E nk+q E mk ~!i jh nk+q je iq~ r j mk ij 2 : (2.120) The same expression can be derived through the direct use of propagators as done in Ref. [92] 36 CHAPTER 2. GENERAL THEORY 2.4. PARTICLES INTERACTING We can evaluate this once we have the singe particle energies, E nk , and wavefunctions, j nk i. Assuming the RPA-screened Coulomb interaction constitutes the primary scattering process leads to the so-called GW approximation. GW approximation We can write the GW approximation simply as setting the self-energy (~ r;~ r 0 ;) =i~G 0 (~ r;~ r 0 ;)W (~ r;~ r 0 ;); or in energy space as (~ r;~ r 0 ;!) = i~ 2 Z 1 1 G 0 (~ r;~ r 0 ;! +! 0 )W (~ r;~ r 0 ;! 0 )e i! 0 d! 0 (2.121) whereW is the screened Coulomb interaction from Eq. 2.118 andG 0 is the Green's function of the non-interacting system given by Eq. 2.109. By the denition of a Green's function we know that the screened propagator will lead to the quasiparticle equation ^ H 0 (~ r) i (~ r) + Z (~ r;~ r 0 ; i =~) i (~ r 0 )d~ r 0 = i i (~ r) (2.122) where ^ H 0 is the Hamiltonian corresponding to G 0 and i and i are the eigenfunctions and eigenvalues of the interacting system, respectively. This equation closely resembles the Kohn-Sham equations, for which we already have solutions coming from DFT. If we assume that the Kohn-Sham energies and wavefunctions are close matches for the proper non-interacting eigensystem (recall DFT only ensures a match between ground states) we can use rst-order perturbation theory to obtain approximate energies i KS i + KS i ( i =~)v xc KS i ; (2.123) where v xc is the exchange correlation potential (see Eq. 2.28) and the superscript 'KS' is just to stress the eigensystem from the Kohn-Sham Hamiltonian. In order to avoid the necessity to nd the full frequency dependence of the self-energy we use the linear expansion (~ r;~ r 0 ; i =~) (~ r;~ r 0 ; KS i =~) + ( i KS i ) ~ @(~ r;~ r 0 ; KS i =~) @! ; giving i = KS i +Z i KS i ( KS i =~)v xc KS i : (2.124) Here we introduced the quasiparticle renormalization factor: Z i = 1 KS i 1 ~ @( KS i =~) @! KS i 1 = Z j i (~ r)j 2 d~ r (2.125) 37 2.4. PARTICLES INTERACTING CHAPTER 2. GENERAL THEORY which equals the quasiparticle weight. This gives us a starting point to move from DFT results to a more accurate description within the GW approximation. We start with the Kohn-Sham solutions ( KS i and KS i ) to calculate W and use G 0 =G KS 0 to get . Then use Eq. 2.124 and Eq. 2.125 to calculate new energies leading to a new G 0 , etc. Unfortunately this self-consistency is extremely expensive from a computational perspective which leads people to use variants of GW in order to save time. These includeG0W 0 where only the rst iteration is used i.e. the Kohn- Sham system, or G0W where only the screened Coulomb interaction is updated at each iteration while G 0 =G KS 0 is kept xed, or similarly one could use GW 0. The calculation of the screened Coulomb interaction, W , is the most time-consuming step which have led to the use of a plasmon-pole model for the inverse dielectric function [96{98] since it allows for the analytical evaluation of the integral in Eq. 2.121. The plasmon-pole is usually implemented in Matsubara frequencies since this gives a smooth screened Coulomb interaction. 2.4.2.1 Electron-plasmon coupling The fully screened, dynamical Coulomb interaction can be written within the plasmon-pole approximation as [99] ~ W k;k 0(i m ) =W k;k 0(! = 0) + Np X j a j;q 2 ! j;q 2! j;q 2 m +! 2 j;q ! ; (2.126) where i m is the bosonic Matsubara frequency point, W k;k 0(! = 0) is the static screened Coulomb interaction, ! j;q is the j'th plasmon frequency with corresponding a j;q electron- plasmon coupling strength and k 0 = k + q. Throughout this thesis we will focus on the case where only a single plasmon pole is assumed i.e. N p = 1, which allows us to drop the j-index. The plasmon-pole approximation allows us to fully dene the screened Coulomb inter- action with only two free parameters, the plasmon frequency (! q ) and the electron-plasmon coupling strength (a q ). We have already discussed, in Section 2.3.1, how we are able to ex- tract the plasmon dispersion from the calculated dielectric function, "(q;!) (see Eq. 2.119 - 2.120) using the EELS spectrum, Eq. 2.78. Once we have the plasmon frequencies we can extract the electron-plasmon coupling by calculating the screened Coulomb interaction for a single, very large Matsubara frequency point and using, lim m!1 ~ W k;k 0(i m ) =W k;k 0(! = 0) + 2a q ! q )a q = ! q 2 [W k;k 0(! = 0) ~ W k;k 0(i m !1)]: (2.127) 38 CHAPTER 2. GENERAL THEORY 2.4. PARTICLES INTERACTING In a simplied system of free electrons we can derive an exact expression for the electron- plasmon coupling, following the derivation of Caruso and Giustino (2016) [100]. The spec- tral representation of the screened Coulomb interaction is given by, W GG 0(q;!) = v(q + G) Z 2! 0 ! 2 (! 0 ) 2 Im" 1 GG 0 d! 0 (2.128) where v(q + G) is the bare Coulomb interaction. We can decompose the dielectric matrix into it's macroscopic and microscopic parts, " 1 GG 0 (q;!) =" 1 M (q + G;!) GG 0 +" 1 GG 0 (q;!)(1 GG 0): (2.129) We know that the macroscopic dielectric function goes to 0 at the plasmon frequency (by denition) and therefore we can nd the plasmonic contribution to "(q;!) as a rst order Taylor expansion of " M , " p (q + G;!) = @" M @! j !=!p(q) (!! p (q)) +i (2.130) Finally, making use of the identity (a +i) 1 =P (1=a) +i(a) we get, Im 1 GG 0 =(!! p (q)) @ M @! j !=!p 1 In the case of a free electron gas we know the macroscopic dielectric function is given by " M = 1 ! 2 p ! 2 therefore we have Im" 1 GG 0 =(!! p (q)) " 2 ! 2 p ! 3 j !=!p # 1 = ! p 2 (!! p (q)) (2.131) This gives the fully screened Coulomb interaction (Eq. 2.128) as, W GG 0(q;!) = v(q + G) Z 2! 0 ! 2 (! 0 ) 2 ! p 2 (!! p (q))d! 0 = v(q + G)! p 2 2! p ! 2 (! p ) 2 (2.132) which has the same form as Eq. 2.126 with a q = v(q)! p (q) 2 : (2.133) We stress that, this result only holds for the free electron gas. However, it provides a useful benchmark for comparing the electron-plasmon coupling strength obtained from the extraction method, discussed earlier, for the case where a system consisting of essentially free electrons is considered. 39 2.4. PARTICLES INTERACTING CHAPTER 2. GENERAL THEORY |k,σi el |qνi ph |k+q,σi el |k,σi el |k+q,σi el |−qνi ph Figure 2.2: Feynman diagrams showing the scattering processes for electron-phonon inter- actions. 2.4.3 Electron-Phonon Interaction The electron-phonon interaction is discussed thoroughly in Ref. [90], the interested reader should refer to that reference for a detailed discussion as we will only provide key results here. The electron-phonon interaction can be described by, V elph = 1 V X k X q X G g qG c y k+q+G; c k; (b q; +b y q; ) whereg is the electron-phonon coupling strength,c y andc (b y andb) are electron (phonon) creation and annihilation operators. The interaction is divided into normal (G = 0) and umklapp (G6= 0) processes. We assume the umklapp processes are not relevant due to vanishing coupling over such large momenta. This reduces the interaction to V elph = 1 V X k X q g q c y k+q; c k; (b q; +b y q; ): (2.134) The physical interpretation of this is shown in the diagrams in Fig. 2.2. Electrons from any initial statejk;i el can be scattered to a nal statejk + q;i el by either absorbing a phonon from statejqi ph or emitting a phonon into the statejqi ph . 2.4.3.1 Electron-phonon coupling The corresponding electron-phonon coupling matrix elements from within DFPT are given by transitions between the Kohn-Sham states when perturbing the ionic potential d i qnn 0 = k + q;n 0 q i v KS jkni as we discussed in Section 2.3.2. This is easily evaluated once the DFPT scheme has reached convergence since the quantities q i v KS jkni are calculated while solving the Sternheimer 40 CHAPTER 2. GENERAL THEORY 2.4. PARTICLES INTERACTING equation (Eq. 2.103). The electron-phonon coupling elements from Eq. 2.134 are then given by g qnn 0 = X i d i qnn 0 e i (q) p 2M i ! (q) (2.135) where M is the ionic mass, e is the phononic polarization vector and ! (q) is the phonon dispersion. The electron-phonon interaction is seen as the standard mechanism that leads to an attractive electron-electron interaction for very low temperatures. In Appendix A.2 we discuss how this attraction comes about. 2.4.4 Spectral Function It is important that we discuss the spectral function since it is the link that binds all the theory we have discussed so far with experiment, seeing as it can directly be measured through angular resolved photo-emission spectroscopy (ARPES) [101]. The spectral func- tion,A(k;!), represents the probability that a particle in a quantum state k has energy! and is dened by A(k;!) = 1 ImG(k;!): (2.136) If we consider, as an example, the spectral function of a free electron with propagator given by Eq. 2.111, we get A(k;!) = 1 Im 1 ! k +i =(! k ): This shows that in the ideal non-interacting case an excitation with energy ! can only occur when adding an electron to the state k where k =!, as one would expect. When we include interactions the spectral function is broadened due to the redistribution of spectral weight resulting from energy exchanges among particles. As a simple example consider the following propagator that decays in time due to interactions, scattering particles out of their current state, G(k;t) =i(t)e i k t e t= ; (2.137) with being the decay time. The spectral function corresponding to this propagator is given by A(k;!) = 1 Im Z 1 1 G(k;t)e i!t dt = 1 Imi Z 1 0 e i k t e t= e i!t dt = 1 1= (! k ) 2 + (1=) 2 ; and shown in Fig. 2.3. Specically we see the full width at half maximum is equal to 2=, hence the width of the spectral function is related to the characteristic decay time of an excitation or, in other words, the lifetime of an excitation. Below we show that the lifetime of an excitation is determined by the self-energy. 41 2.4. PARTICLES INTERACTING CHAPTER 2. GENERAL THEORY −1/τ ǫ k +1/τ A(k,ω) Figure 2.3: Sketch of spectral function, showing broadening due to interaction. We have already discussed the correction to the free particle propagator due to inter- actions that appears in the form of a self-energy. In general the interacting propagator is given by (see Eq. 2.115), G (k;!) = ! k +i ; (2.138) which has the spectral function, A(k;!) = 1 X n Im nn (k;!) [! k Re nn (k;!)] 2 + [Im nn (k;!)] 2 ; (2.139) where we assumed the self-energy is diagonal in the band basis. We can then also dene the total spectral function as [91] A(!) = X k A(k;!); (2.140) which can also be measured experimentally. The spectral function is a very powerful tool to determine which interactions in a system play the most important role since we can discriminately apply the self-energy coming from individual interactions when calculating the spectral function. These results can then be compared to experiment to determine which interactions are important. 42 Chapter 3 Superconductivity Theory The following discussion of superconductivity came mostly from Ref. [90] chapters 4 and 16, Ref. [57] chapter XVIII and the SC-DFT articles Ref. [102] and Ref [99]. The PhD thesis of Miguel Marques [103] was also heavily referenced. It should be noted that the derivation of McMillan's equation from SC-DFT and the subsequent inclusion of plasmonic eects to the equation are original to this dissertation. In order to discuss the theory of superconductivity we rst have to say a few general things concerning phase transitions. Consider a system governed by a Hamiltonian, H, with some symmetry, re ected by the symmetry operator, S. Naturally, the Hamiltonian commutes with S i.e. [H;S] = 0 and there exists a set of eigenstates that simultaneously diagonalizes H and S. To proceed more concretely, lets assume the symmetry at hand is translational invariance so S = T ( ~ R) = exp i ~ R P , where P is the total momentum operator. The orthogonal basis set of denite total momentumjPi thus forms a common set of eigenvectors for H and T . Now consider the density operator, (Q) = X k c y k c dagger k+Q : In a charge density wave phase the Fourier transform of the density operator would have a nite value, but hc y k c dagger k+Q i = 1 Z X P e P hPjc y k c dagger k+Q jPi = 0; since c y k c dagger k+Q jPi!jP Qi which is orthogonal tojPi. Hence, there is a contradiction seeing as charge density waves are possible in nature contrary to what is predicted above. This points to a breakdown 43 CHAPTER 3. SUPERCONDUCTIVITY T >T c T <T c ⇒ Figure 3.1: Sketch of energy diagram of a phase transition showing above the transition temperature the well-dened minimum (left) while below the transition temperature the phase space is split into separate regions (right). Below the transition temperature the overall symmetry is still maintained but not within the individual pockets where the system has a lower symmetry. Table 3.1: Examples of phase transitions accompanied by spontaneous symmetry breaking and the corresponding order parameters. The table was copied from Ref. [90]. Phenomena Order parameter - mathematical Order parameter - physically crystal P k hc y k c k+Q i density wave ferromagnet P k hc y k" c k" c y k# c k# i magnetization Bose-Einstein condensate ha y k=0 i population of k = 0 state superconductor hc k" c k# i pair condensate of ergodicity for the system in a charge density wave phase implying the phase space of the system is divided into physically separated regions and therefore the sum of states in the thermodynamic average should be restricted. This can be visualized with the energy diagram of a phase transition shown in Fig. 3.1. The critical temperature for a phase transition is dened as the temperature where the system develops a nite value for the expectation value of some macroscopic quantity (called the order parameter) that has a lower symmetry than the original Hamiltonian. This is called spontaneous symmetry breaking. In Table 3.1 we show several examples of phase transitions due to broken symmetry and the corresponding order parameter. We will now demonstrate specically how the breaking of global gauge symmetry leads to the remarkable superconducting phase. In such a phase a metal displays zero resistivity, meaning that a current would ow in the metal for years without measurable decay. Furthermore, in the superconducting phase, metals have perfect diamagnetism (Meissner eect) meaning that a magnetic eld (below a critical value) cannot penetrate into the interior of the metal. As we mentioned above the symmetry of interest for the onset of superconductivity 44 CHAPTER 3. SUPERCONDUCTIVITY is the global U(1) gauge symmetry, which is to say, the Hamiltonian is invariant under a phase change of the electrons i.e. c n !c n e i' =) H!H The superconducting order parameter has the formhc n c n 0i, and as discussed above, this order parameter would have a nite value if the U(1) gauge symmetry is broken. Now, obviously, a global phase change cannot have any eect but if we apply a position depen- dent phase shift, '(~ r), things become more interesting. We dene the following unitary transformation to study the phase change U = exp i Z (~ r)'(~ r)d~ r ; (3.1) where we used the density operator, (~ r) = y (~ r) (~ r). To see that this unitary indeed shifts the phase of a quantum eld operator, notice, ~ (~ r) =U (~ r)U 1 = exp i Z (~ r)'(~ r)d~ r (~ r) exp i Z (~ r)'(~ r)d~ r = (~ r)e i'(~ r) ~ y (~ r) =U y (~ r)U 1 = y (~ r)e i'(~ r) coming from the boundary condition ~ '=0 = and the relation @ @'(~ r 0 ) ~ (~ r) =iU[(~ r 0 ); (~ r)]U 1 =i ~ (~ r)(~ r~ r 0 ): If we calculate the transformed partition function, we get ~ Z = Tr 0 [e ~ H ] = Tr 0 [Ue H U 1 ] (3.2) where we use Tr 0 to stress that the trace is restricted, similar to the restriction discussed earlier, due to the spontaneous symmetry breaking. In a mean eld Hamiltonian the transformation only aects the kinetic energy term since all other interactions only involve the electronic density operator,(~ r), which is obviously invariant under the transformation. The kinetic energy term is given by ~ T = 1 2m e Z y (~ r)e i'(~ r) ~ i r +e ~ A 2 e i'(~ r) (~ r)d~ r = 1 2m e Z y (~ r) ~ i r +e ~ A~r'(~ r) 2 (~ r)d~ r =T~ Z r'(~ r) ~ J(~ r)d~ r + ~ 2 2m e Z (~ r)(r'(~ r)) 2 d~ r (3.3) The order parameter includes two annihilation operators, however, we assume the system is connected to an electron reservoir such that the particle number would remain constant 45 3.1. BCS-THEORY CHAPTER 3. SUPERCONDUCTIVITY where ~ A is the vector potential and ~ J the resulting current density. Let '(~ r) be such that it varies extremely slowly so that a dierence is only observed over macroscopic distances. In this case the phase is a macroscopic quantity governed by thermodynamic laws and therefore aects the free energy of the system. In fact, from Eq. 3.3 it is evident that'(~ r) does not aect the energy but ratherr'(~ r) does. By minimizing the free energy we can nd the equilibrium properties of the system, @F @r' =~h ~ J(~ r)i + ~ 2 m e h(~ r)ir'(~ r) = 0; (3.4) showing that in equilibrium the system supports a current h ~ J(~ r)i = ~ 0 m e r': A current can only be carried in equilibrium if it is dissipationless, indicating zero resistivity in the system i.e. a superconducting phase. There is still an energy cost to carrying a current and the system will revert back to the normal state if that energy cost equals the superconducting gap. 3.1 BCS-theory The rst proper explanation of superconductivity was due to Bardeen Cooper and Schrieer (1957) [15]. In the following we will review their ndings in establishing the idea of Cooper pairs and consequently the condensation of these pairs leading to a superconducting state. A crucial milestone in the development of BCS theory was Cooper's demonstration that the Fermi sea would become unstable if an arbitrarily small, attractive interaction operates between electrons in the system [104]. 3.1.1 Cooper instability Consider a metal in the ground state i.e. with electrons in all states with k< k F , where k F is the Fermi-vector. Now assume we add two electrons to the system that interact with each other through an attractive two-body interaction, V (~ r 1 ;~ r 2 ). The electrons will reside in some statesjk 1 i andjk 2 i from which we can form the two-body wavefunction through a Slater determinant. The two-body state will have total momentum ~K = ~k 1 + ~k 2 . Since the interaction,V (~ r 1 ;~ r 2 ), necessarily conserves momentum the system can only tran- sition to other states with total momentum K. Now, since we are interested in the lowest energy behavior of the system we focus on the case when the center of mass of the two electrons are still, i.e. K = 0 which implies k 1 =k 2 . Furthermore, it is well known that the combination of two fermions lead to a singlet and triplet state of which the singlet state has spin 0. For our purposes it is sucient to focus on this case, we will therefore be 46 CHAPTER 3. SUPERCONDUCTIVITY 3.1. BCS-THEORY interested in the consequences of an attractive interaction between electron pairs occupying time-reversed statesjk"i andjk#i. The total interaction can be represented as the pair scattering vertex ( ~ k; ~ p) : k↑ −k↓ p↑ −p↓ Λ = k↑ −k↓ p↑ −p↓ + k↑ −k↓ V p↑ −p↓ + k↑ −k↓ V V p↑ −p↓ + ... where we see can be written in a Dyson equation-like form, Λ = + Λ with the equivalent algebraic formula ( ~ k; ~ p) =V ( ~ k ~ p) + 1 V X ~ q [V ( ~ k ~ q)]G 0 " (~ q)G 0 # (~ q)( ~ k; ~ p): (3.5) In order to proceed we have to make some assumptions regarding the properties of the attractive interaction V . We know that in any real material electrons separated by large q (corresponding to small separation in real-space) will strongly repel each other due to the Coulomb interaction. Thus any attractive interaction mediated through some other mechanism will only exist for small momenta separations. Anticipating phonons to play the mediating role we assume the attraction happens only for energies below the characteristic phonon frequency, namely, the Debye frequency, ! D . The assumed form of V is therefore: V (q;iq n ) = ( V 0 forjiq n j< ~! D 0 forjiq n j> ~! D (3.6) where V 0 as a positive constant. This allows us to simplify Eq. 3.5, giving (k; p) =V 0 + 1 V X jiqnj<~! D X q V 0 G 0 " (q;iq n )G 0 # (q;iq n )(k; p): (3.7) Now, notice that k never appears in the summand, hence we can drop the index. Fur- thermore, the momentum p only appears in which means we can still nd a consistent solution when we assume (k; p) = is a constant. We can then solve Eq. 3.7 to get = V 0 1 V 0 V P jiqnj<~! D P q G 0 " (q;iq n )G 0 # (q;iq n ) : (3.8) we use the notation ~ k = (k;ikn) with ikn being a fermionic Matsubara frequency 47 3.1. BCS-THEORY CHAPTER 3. SUPERCONDUCTIVITY We can now interpret the aect of an attractive interaction between electrons on the system. At high temperatures (small ) we have V 0 and since we assume V 0 is small the aect is negligible. However, if the temperature is decreased ( increased) the denominator in Eq. 3.8 becomes smaller and smaller leading to an increasing value for . As the denominator approaches zero, the scattering amplitude diverges, indicating a resonance in the system, which in this case is the formation of a bound state between the time-reversed electron pair. This bounded pair of electrons is called a Cooper pair. At this critical temperature where diverges the whole Fermi surface becomes unstable and Cooper pairs are formed between all time-reversed pairs of electrons within the ~! D energy window. This phenomenon is called the Cooper instability and it marks the transition from a normal to a superconducting state. The critical temperature, T c , at which the transition occurs can be found by setting the denominator of Eq. 3.8 equal to zero i.e. 1 = V 0 V c X jiqnj<~! D X q G 0 " (q;iq n )G 0 # (q;iq n ): (3.9) For simplicity we assume a high density of states such that we can make the replacement P q ! R d" where is the unit cell volume, meaningV =N withN being the number of unit cells in the supercell. We can also write the normalized density of states per unit cell as n(") = 1=N. Finally, recall the form of fermionic Matsubara frequencies, q n = 2 (n + 1 2 ); and the free electron propagator G 0 (q;iq n ) = 1 (iq n " q ) : Substituting all this into Eq. 3.9 gives 1 = Vn(" F ) c X jiqnj<~! D Z 1 1 1 iq n " 1 iq n " d" = V 0 n(" F ) c X jiqnj<~! D Z 1 1 1 q 2 n +" 2 d" = V 0 n(" F ) c X jiqnj<~! D jq n j (change sum-variables: q n to n) = V 0 n(" F ) 2 0 @ 1 2 c~! D X n=0 1 n + 1 2 2 1 A V 0 n(" F ) 2 ln 2 c ~! D ; (3.10) 48 CHAPTER 3. SUPERCONDUCTIVITY 3.1. BCS-THEORY from which we can solve for T c according to, 2 ~! D k B T c = exp 2 V 0 n(" F ) ) k B T c = 2 ~! D e 2=(V 0 n(" F )) (3.11) We can learn three important things concerning the transition temperature in conven- tional superconductivity from the T c expression in Eq. 3.11, namely, 1. T c increases if the phonon frequency ! D increases 2. T c increases strongly if the strength of the attractive interaction increases 3. a high density of states at the Fermi-level leads to a high T c . Also, note that the transition temperature is a non-analytic function of V 0 which means that it is not possible to describe the Cooper instability by any order of perturbation theory as is usually the case - it is not possible to cross a phase boundary through perturbation theory. 3.1.2 BCS Gap Equation After the discovery of the Cooper instability it was pretty clear that superconductivity is the result of a Bose-Einstein condensate formed by Cooper-pairs. In order to further under- standing of the superconducting state BCS suggested the following eective Hamiltonian that they anticipated would capture the important properties of the state, H BCS = X k k c y k c y k + X kk 0 V kk 0c y k" c y k# c y k 0 " c y k 0 # ; whereV kk 0 is a small coupling of the form in Eq. 3.6 and k =" k " F , is the electron energy relative to the Fermi-level. We assume the Coulomb repulsion is already included in the k term. We now use a mean-eld approximation , assuming the pair operator b k =c k# c k" to have a nite expectation value and small deviance (dened by c k# c k" hc k# c k" i), so that the Hamiltonian becomes H BCS = X k k c y k c y k X k y k c y k" c y k# X k k c y k# c y k" + X kk 0 V kk 0hc y k" c y k# ihc y k 0 # c y k 0 " i (3.12) where k = X k 0 V kk 0hc k 0 # c k 0 " i (3.13) In the mean-eld approximation, if we have an interaction given by the product of two operators, H =AB we simplify the situation by assuming A couples to the mean eld of B and vice versa such that H MF =AhBi +hAiBhAihBi. 49 3.1. BCS-THEORY CHAPTER 3. SUPERCONDUCTIVITY The Hamiltonian is easily solved with the use of the Bogoliubov transformation [105]. To this end we write the Hamiltonian in matrix form H BCS = X k c y k" c y k# k k k k c y k" c y k# ! + X k k + X kk 0 V kk 0hc y k" c y k# ihc y k 0 # c y k 0 " i; (3.14) but the last two terms are constants and can therefore be ignored leaving us with H BCS = X k a y k H y k a y k : (3.15) Let U k be the transformation matrix that diagonalizes the Hamiltonian, and let U 1 k a k = u k v k v k u k 1 c y k" c y k# ! = b k = y k" y k# ! : (3.16) Combining Eq. 3.16 with the form of a k and H k from Eq. 3.14 and using U y k H y k U y k = E k 0 0 ~ E k ; we nd after some algebra the follwing solutions for the coecients u and v, ju k j 2 = 1 2 1 + k E k jv k j 2 = 1 2 1 k E k (3.17a) and the following solutions for the energies of the superconducting state E k = ~ E k = q 2 k +j k j 2 : (3.17b) The solution of the BCS Hamiltonian shows that in the superconducting phase a gap ( k ) is created in the electron energy spectrum, which is appropriately called the super- conducting gap. By combining the results from above we get k = X k 0 V kk 0hc k 0 # c k 0 " i = X k 0 V kk 0 u k 0v k 0h y k 0 # y k 0 # iv k 0u k 0h y k 0 " y k 0 " i = X k 0 V kk 0u k 0v k 0[1 2f(E k 0)] (3.18) 50 CHAPTER 3. SUPERCONDUCTIVITY 3.1. BCS-THEORY where f(E) is the Fermi-Dirac distribution value at energy E. Plugging in Eq. 3.17a we get an expression that has to be solved self-consistently k = X k 0 V kk 0 1 2 1 2 k 0 2 k 0 + 2 k 0 1=2 [1 2f(E k 0)] = 1 2 X k 0 V kk 0 k 0 q 2 k 0 + 2 k 0 [1 2f(E k 0)]: Finally, we use the identity 1 2f(E) tanh E 2 ; to get the BCS equation for the superconducting gap: k = 1 2 X k 0 V kk 0 k 0 q 2 k 0 + 2 k 0 tanh q 2 k 0 + 2 k 0 2 : (3.19) If we now insert the form ofV kk 0 from Eq. 3.6 and we assume the density of states over the window [! D ;! D ] is constant and equal to n(E F ). And if we further assume the gap to be constant in the same window and equal to 0 , Eq. 3.19 simplies to 1 = 1 2 V 0 n(E F ) Z ~! D ~! D 1 p " 2 + 2 0 tanh p " 2 + 2 0 2 d": (3.20) Note that the gap has it's maximum value at T = 0 and decreases for increasing tempera- ture. The temperature at which the metal transitions back to the normal state is the same temperature at which the gap becomes 0. We can therefore use Eq. 3.20 to get an estimate for the transition temperature, since at T c 1 =V 0 n(E F ) Z ~! D 0 1 " tanh " 2k B T c d" which we rewrite as 1 V 0 n(E F ) = Z ~! D =(k B Tc) 0 1 x tanh x 2 dx: For ~! D >>k B T c the integral is equal to ln 2 exp( ) ~! D k B Tc where is the Euler constant, which gives k B T c = 2 exp( ) ~! D e 1=(V 0 n(E F )) 1:13~! D e 1=(V 0 n(E F )) : (3.21) We see that Eq. 3.11 and Eq. 3.21 give very similar values forT c but more importantly they show the same functional dependence on ! D , V 0 and n(E F ). These characteristic dependencies have been experimentally conrmed in countless simple superconductors. 51 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY 3.2 DFT for Superconductors The BCS-theory presented so far played an enormously important role in our basic under- standing of the phenomenon of superconductivity. But in order to calculate the properties of realistic superconducting systems we require a more in-depth theory. The go-to theory in this case is Eliashberg theory (see Ref. [106] for a detailed description). In short Eliash- berg theory provides a proper Green's function treatment of the theory of superconductors. This unfortunately has the draw back of being computational intensive, even to the point of being infeasible to use for complicated systems. Specically, the Coulomb interaction is generally treated in an approximated, averaged manner (using ) which neglects the contribution of variational screening on the superconducting state. In Section 2.2.2 we discussed how DFT can be used to simplify the complicated problem of many-body interactions in crystals. In essence the DFT technique turns the full Green's function theory of interacting systems into a computationally tractable problem. Oliveira, Gross and Kohn in 1988 [107] developed a similar approach to solving the superconducting system and aptly named it density-functional theory for superconductors (SC-DFT). In the following we will provide an overview of the approach along with some original ndings by the author regarding the theory. We start with the following Hamiltonian to describe the superconducting system: H = X Z y (~ r) ~ 2 2m r 2 +v ext (~ r) y (~ r)d~ r +U +W (3.22) where is the chemical potential, U is Coulomb repulsion and W is some attractive electron-electron interaction (as discussed in the previous section). The following term is added to the Hamiltonian, as well, ^ = ZZ (~ r;~ r 0 ) " (~ r) # (~ r 0 )d~ rd~ r 0 +h:c: (3.23) as a way to mathematically break the gauge invariance that is spontaneously broken in nature at the onset of superconductivity (see Table 3.1). The term can be physically thought of as a pairing eld to a nearby superconductor that induces superconductivity or simply a mathematical artifact that can be taken to zero at the end of the derivation. Mermin (1965) [108] showed that a version of the Hohenberg-Kohn theorem (Section 2.2.2) is also valid for a system at thermodynamic equilibrium. The electron density and so-called 'anomalous density' is given by, n(~ r) = X h y (~ r) y (~ r)i (3.24a) (~ r;~ r 0 ) =h " (~ r) # (~ r 0 )i: (3.24b) 52 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS Analogous to how we obtained Eq. 2.25, the grand-canonical thermodynamic potential of the Hamiltonian above, can be written as [n;] =F [n;] + Z n(~ r)[v ext (~ r)]d~ r ZZ [ (~ r;~ r 0 )(~ r;~ r 0 ) +c:c:]d~ rd~ r 0 ; (3.25) with F [n;] =T [n;] +U[n;] +W [n;] 1 S[n;]: The kinetic energy, T [n;], and entropy, S[n;], functionals are of the interacting system while U[n;] and W [n;] are density functionals of the thermal averages of the respec- tive interactions. Similarly as we did in the derivation of DFT we now construct a non- interacting Hamiltonian (like the Kohn-Sham Hamiltonian of Eq. 2.22) such that the non- interacting densities n 0 (~ r) and 0 (~ r) match the interacting densities. The non-interacting Hamiltonian can be written H 0 = X Z y (~ r) ~ 2 2m r 2 +v 0 (~ r) y (~ r)d~ r ZZ 0 (~ r;~ r 0 ) " (~ r) # (~ r 0 )d~ rd~ r 0 +h:c: ; (3.26) where v 0 (~ r) =v ext (~ r) + Z n(~ r 0 ) j~ r~ r 0 j d~ r 0 +V xc [n;](~ r) (3.27a) 0 (~ r;~ r 0 ) = (~ r;~ r 0 ) j~ r~ r 0 j + xc [n;](~ r;~ r 0 ): (3.27b) Notice the introduction of the so-called anomalous Hartree potential (rst term in Eq. 3.27b). As before, the connection between the interacting and non-interacting systems is in the exchange-correlation potentials that are found through the variational procedure (see Eq. 2.28) as, V xc [n;](~ r) = F xc [n;] n(~ r) (3.28a) xc [n;](~ r;~ r 0 ) = F xc [n;] (~ r) ; (3.28b) with (similar to Eq. 2.27) F [n;] =T 0 [n;] 1 S 0 [n;] +E ee H [n;] +F xc [n;]; where E ee H [n;] = 1 2 ZZ n(~ r)n(~ r 0 ) j~ r~ r 0 j d~ rd~ r 0 + ZZ j(~ r;~ r 0 )j 2 j~ r~ r 0 j d~ rd~ r 0 (3.29) 53 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY The Hamiltonian in Eq. 3.26 is bilinear in y ; y so according to the Bogoliubov-de Gennes ndings [105], the Hamiltonian can be diagonalized by writing (~ r) = X i ( i u i (~ r) + y i v y i (~ r)) where y and y are the creation and annihilation operators of some new eld and u i (~ r) and v i (~ r) are particle and hole amplitudes that satisfy the Bogoliubov-de Gennes (BdG) equations ~ 2 2m r 2 +v 0 (~ r) u i (~ r) + Z 0 (~ r;~ r 0 )v i (~ r 0 )d~ r 0 = ~ E i u i (~ r) (3.30a) ~ 2 2m r 2 +v 0 (~ r) v i (~ r) + Z 0 (~ r;~ r 0 )u i (~ r 0 )d~ r 0 = ~ E i v i (~ r): (3.30b) A full numerical solution of these equations would require a very high resolution since typical values of 0 are orders of magnitude smaller than typical electronic energies. How- ever, this large energy scale separations allows us to simplify the equations through the decoupling approximation [109]. Within this approximation we use the solutions of the normal-state (non-superconducting) Kohn-Sham equation (also given in Eq. 2.29) ~ 2 2m r 2 +v 0 (~ r) i (~ r) = i i (~ r); (3.31) to expand the particle and hole amplitudes. We assume we only need the rst term in the expansion, thereby giving u i (~ r)u i i (~ r) v i (~ r)v i i (~ r): (3.32) Finally, the matrix elements of the superconducting gap, in the Kohn-Sham basis is dened as i = ZZ i (~ r) 0 (~ r;~ r 0 ) i (~ r 0 )d~ r 0 d~ r: (3.33) Plugging Eq. 3.31 - 3.33 into Eq. 3.30 gives i i i i u i v i = ~ E i u i v i (3.34) where we used the new variable i i . Diagonalizing the matrix we get ~ E i =E i with eigenvalues y E i = q 2 i +j i j 2 : The normal-state solutions,i, are found through the usual DFT method we discussed in Section 2.2.2 y The same result as we got in Eq. 3.17b 54 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS We can also solve for the coecients and doing so gives u i = 1 p 2 sgn( ~ E i )e i' s 1 + i ~ E i v i = 1 p 2 s 1 i ~ E i ; (3.35) where e i' = i =j i j is a phase. Also within the approximation, using Eq. 3.32 and Eq. 3.35 along with Eq. 3.24 we can nd the normal and anomalous densities n(~ r) = X i 1 i E i tanh 2 E i j i (~ r)j 2 (3.36a) (~ r;~ r 0 ) = 1 2 X i i E i tanh 2 E i i (~ r) i (~ r 0 ): (3.36b) Recall that the index i are labels for the quantum state and particularly, in a crystal, the momentum k is a good quantum number, so we can rewrite all the equations above with k as the labeling index. If we plug Eq. 3.36b into Eq. 3.27b and the resulting equation for 0 into Eq. 3.33 we get the matrix elements of the superconducting gap ( i ) as a functional of the chemical potential and i itself. This can be solved self-consistently as follows: 1. using standard DFT code we obtain the electronic energies and eigenfunctions 2. using standard DFPT we can nd phonon eigenmodes which can be incorporated in W 3. the gap-equation is iterated until self-consistency is achieved. The superconducting transition temperature is dened as the temperature where 0 ! 0, since at this temperature the anomalous density vanishes. Hence we can linearize the gap equation in , giving i = 1 2 X j F Hxci;j [] tanh[( c =2) j ] j j ; (3.37) where we used c to stress that this equation is only valid at the transition temperature. The anomalous Hartree exchange-correlation kernel is also introduced and reads F Hxci;j [] = Hxci j =0 = 2 (E ee H +F xc ) i j =0 : (3.38) 55 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY We now split the kernel,F Hxci;j [], into two parts (mostly intuitively), a diagonal term, Z i , and the rest,K i;j i =Z i i 1 2 X j K i;j [] tanh[( c =2) j ] j j : (3.39) This transforms the gap equation into a similar form as the BCS equation (Eq. 3.19) where K ij plays the role of V kk 0 i.e. the Cooper-pair binding interaction. AndZ i resembles the mass-renormalization term in Eliashberg theory (see Ref. [106]). The big convenience of SC-DFT lies in this simple form of the gap equation. By studying the kernel,K, we can learn which mechanism are responsible for superconductivity in particular systems. We can also easily include more-and-more eects inK and see how the transition temperature changes. The fact that! is not explicitly included in the gap equation makes the SC-DFT gap equation much less numerically intensive to solve than the Eliashberg equations, but retardation eects are still included in the theory through the exchange-correlation terms. Thus SC-DFT is not a mean-eld theory such as BCS but in principle should be able to account for much of the same eects as Eliashberg theory. We explore this point more in Section 5.5. The diculty now is that without knowledge of the interaction kernel we can't pro- ceed. Seeing as the theory is built on the DFT strategy we here have the same prob- lem discussed in Section 2.2.2: in order to nd the true Kohn-Sham system we need the exchange-correlation functionals, which we don't have. And same as before the solution is to nd appropriate approximations for the functional that give results well matched with experimental values. 3.2.1 Approximation for the exchange-correlation kernel The technique to approximate the exchange-correlation functional is laid out in Ref. [102] and discussed in detail in Ref. [103]. Here we will present in broad strokes the strategy and main results. The starting point is to frame the problem in terms of perturbation theory . To this end we dene the unperturbed Hamiltonian as the Kohn-Sham Hamiltonian (Eq. 3.26) with a phononic part, H ph , added to account for lattice vibrations given by H ph = X ;q ;q b y ;q b y ;q + 1 2 ; where ;q is the phonon frequency and b y ;q (b y ;q ) the phononic creation (annihilation) operator. The unperturbed Hamiltonian is therefore H s =H el +H ph +V Hxc + 0 ; specically the Kohn-Sham perturbation theory of G orling and Levy [110] 56 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS to which we then apply the perturbation H 1 =U +H elph V Hxc xc : As in standard perturbation theory we expand the dierence in the grand-canonical po- tential, , between the perturbed and unperturbed systems as a power series. The grand- canonical potential is dened as = 1 ln(Z) where Z is the partition function given by Z = Tr e H : It is therefore easy to see that s = 1 ln Z Z s In this case the perturbed grand-canonical potential is given by Eq. 3.25 (with added term for phononic parts) while the unperturbed grand-canonical potential is given by a similar term for the Kohn-Sham Hamiltonian. Here we refer to Ref. [103] for the calculation details which eventually gives s = 1 where refers to the irreducible self-energy i.e. the sum of all closed Feynman diagrams that are topologically distinct and can't be separated into two parts by cutting a single propagator line. The propagators that appear in the Feynman diagrams of the self-energy can be represented as G 0 σσ ′(~ rt,~ r ′ t ′ ) = (~ r ′ t ′ σ ′ ) (~ rtσ) sgn(σ)F 0 σσ ′(~ rt,~ r ′ t ′ ) = (~ r ′ t ′ σ ′ ) (~ rtσ) sgn(σ)F 0† σσ ′(~ rt,~ r ′ t ′ ) = (~ r ′ t ′ σ ′ ) (~ rtσ) 57 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY with the algebraic form for each given by G 0 0(~ rt;~ r 0 t 0 ) = ; 0 X i u i (~ r)u i (~ r 0 ) i! n E i + v i (~ r)v i (~ r 0 ) i! n +E i F 0 0(~ rt;~ r 0 t 0 ) = ; 0sgn( 0 ) X i v i (~ r)u i (~ r 0 ) i! n +E i u i (~ r)v i (~ r 0 ) i! n E i F 0y 0 (~ rt;~ r 0 t 0 ) = ; 0sgn( 0 ) X i u i (~ r)v i (~ r 0 ) i! n +E i v i (~ r)u i (~ r 0 ) i! n E i ; (3.40) and the bosonic propagator given by D 0 ;q ( n ) = 2 ;q 2 n + 2 ;q : (3.41) Note that as usual in nite temperature Green's function formalism we make use of imagi- nary (Matsubara) frequencies and denote by! n the odd (fermionic) Matsubara frequencies ! n = (2n + 1) ; and by n the even (bosonic) Matsubara frequencies n = 2n : Considering interactions up to rst order in the electron-boson coupling and including phononic and plasmonic interactions there are only four diagrams contributing to F xc . We will only draw the two diagrams coming from the electron-phonon interaction although the same diagrams will be taken for the electron-electron (Coulomb) interaction. We take F ph xc = + where we can label the two terms (diagrams) as F (a) xc and F (b) xc , respectively. The analytic form of these terms are given by F (a) xc = 1 2 X ;q X ij jg ij ;q j 2 1 + i j E i E j I(E i ;E j ; ;q ) + 1 i j E i E j I(E i ;E j ; ;q ) (3.42a) 58 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS F (b) xc = 1 2 X ;q X ij jg ij ;q j 2 i j E i E j [I(E i ;E j ; ;q )I(E i ;E j ; ;q )]; (3.42b) where the function I is dened as I(E 1 ;E 2 ; ;q ) = 1 2 X ! 1 ! 2 1 i! 1 E 1 1 i! 2 E 2 2 (! 1 ! 2 ) 2 + ; (3.43) andjg ij ;q j 2 are the electron-phonon coupling elements. The frequency sums in Eq. 3.43 can be done analytically, giving the simplication, I(E 1 ;E 2 ; ;q ) =f (E 1 )f (E 2 )n ( ) " e E 1 e (E 2 + ) E 1 E 2 e E 2 e (E 1 + ) E 1 E 2 + # ; (3.44) where f (E) is the Fermi-Dirac function and n (E) the Bose-Einstein function. We now have F ph xc = F (a) xc +F (b) xc = F ph xc [; i ], showing the exchange correlation function as a functional of and i , but we know these quantities are also functionals of the densities n and . By multiple applications of the chain rule with Eq. 3.38 we are able to get approximations for the exchange-correlation kernel,F Hxci;j (that enters the gap equation, Eq. 3.37). As we did before, we split the exchange-correlation kernel into a diagonal part, Z i , and a non-diagonal part,K i;j . This works out conveniently since the kernel derived from F (a) xc is diagonal and therefore contributes only toZ ph i , giving Z ph i = 1 tanh[(=2) i ] X j X ;q jg ij ;q j 2 [J( i ; j ; ;q ) +J( i ; j ; ;q )] (3.45) where J( 1 ; 2 ; ) = ~ J( 1 ; 2 ; ) ~ J( 1 ; 2 ; ) with ~ J( 1 ; 2 ; ) = f ( 1 ) +n ( ) 1 2 f ( 2 )f ( 1 ) 1 2 f ( 1 )f ( 2 + ) : The part of the kernel coming from F (b) xc contributes toK ph i;j and is given by K ph i;j = 2 tanh[(=2) i ] tanh[(=2) j ] X ;q jg ij ;q j 2 [I( i ; j ; ;q )I( i ; j ; ;q )]: (3.46) This gives us the ingredients necessary to include electron-phonon interactions in our calcu- lation of the superconducting T c . We similarly include the Coulomb interaction by taking into account, we should note that this is not the form ofZ used in this thesis but instead we use a numerically more stable version developed by Akashi and Arita in Ref. [111] 59 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY F el xc = + where instead of using the bare Coulomb interaction we use the more accurate screened Coulomb interaction (discussed in detail in Sec. 2.4.2). This gives, similarly as before, a contribution to the diagonal ofF Hxci;j given by Z el i = 1 tanh[(=2) i ] X j a pl ij [J( i ; j ;! pl ij ) +J( i ; j ;! pl ij )] (3.47) where! pl is the plasmon frequency anda pl the electron-plasmon coupling. The o-diagonal contribution is also mostly similar to before but due to the dierence in the phononic propagator (Eq. 3.41) and the screened Coulomb propagator (given in the plasmon-pole approximation by Eq. 2.126) there are two extra terms inK el ij compared toK ph ij (as shown in Ref [99]). For convenience we can split the kernel up into a static part coming from direct electron-electron interactions through the static Coulomb interaction and a dynamic part coming from electron-plasmon interactions. This then gives, K el i;j =K elel i;j +K elpl i;j =W ij (! = 0) + 2a pl ij 1 ! pl ij + I( i ; j ;! pl ij )I( i ; j ;! pl ij ) tanh[(=2) i ] tanh[(=2) j ] ! (3.48) where W (! = 0) is simply the static, screened Coulomb interaction. 3.2.2 Energy Averaged Kernels The above expressions for the exchange-correlation kernels fully allow us to study the eect of varying electron-phonon, electron-plasmon and Coulombic screening on the supercon- ducting transition temperature. However, since the superconducting gap, i , is several orders of magnitude smaller than the typical electronic energies it is numerically very hard to reach self-consistency in Eq. 3.39, thereby requiring us to use very dense k-grids. In order to circumvent this problem we recast the gap-equation and kernels in energy space and use the transformation as a natural interpolation scheme going from coarse momen- tum grids to very dense energy grids. Specically, we can do this if we consider only s-type superconductors where the pair-potential has s-wave symmetry. The procedure for this is shown below usingK elpl as an example. Specically, the labels (i;j) used above are meant to label quantum numbers; in the setting of crystals we usually use the momentum as the primary quantum number, so we 60 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS have, K elpl k;k 0 = 2a pl k;k 0 1 ! pl k;k 0 + I( k ; k 0;! pl k;k 0 )I( k ; k 0;! pl k;k 0 ) tanh[(=2) k ] tanh[(=2) k 0] ! (3.49) We now use the Dirac-delta function to change variables (according to f(t) = R f(x)(t x)dx) which gives, K elpl (; 0 ) = 1 N()N( 0 ) X k;k 0 ( k )( k 0)K elpl k;k 0 ; (3.50) and when plugging in Eq. 3.49 we get K elpl (; 0 ) = 1 N()N( 0 ) X k;k 0 ( k )( k 0)a k;k 0 2 ! k;k 0 + 2 I(;;! k;k 0)I(;;! k;k 0) tanh[(=2)] tanh[(=2)] : (3.51) Now perform the same change of variables with ! k;k 0!! to yield, K elpl (; 0 ) = 1 N(0) Z 2 F pl (; 0 ;!) 2 ! + 2 I(; 0 ;!)I(; 0 ;!) tanh[(=2)] tanh[(=2) 0 ] d! (3.52) where we introduce the function 2 F pl , dened as 2 F pl (; 0 ;!) N(0) N()N( 0 ) X k;k 0 a k;k 0( k )( 0 k 0)(!! k;k 0): (3.53) Commonly in the calculation of 2 F one assumes the most important contributions come from the Fermi surface and so you set = 0 = 0. Here it is crucially important to note that this approximation can only be made in the case where 2 F ( = 0 ;!! 0)! 0 which happens if a k;k 0 +q ! 0 as q! 0. This is not the case for pure 2d plasmons which is discussed in depth in Section 3.2.4. However, let's rst consider the simple case where the standard assumption holds (as is the case for plasmons in 3d). If we also assume a k;k 0 and ! k;k 0 only depend on q = k k 0 , we can simplify Eq. 3.53 to, 2 F pl (!) = 1 N(0) X k;q a q ( k )( kq )(!! q ); (3.54) in accordance with the standard denition of the Eliashberg function. The energy averaged dynamic, electronic kernel can now be separated into two simple terms as, K elpl (; 0 ) = 1 N(0) Z 2 2 F pl (!) ! d! + Z 2 2 F pl (!) I(; 0 ;!)I(; 0 ;!) tanh[(=2)] tanh[(=2) 0 ] d! = pl N(0) +K pl (; 0 ) (3.55) 61 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY where we dened a new functionK pl and made use of the standard, averaged bosonic coupling parameter pl = Z 2 2 F pl (!) ! d!: (3.56) Doing the same averaging for the phononic kernel from Eq. 3.46, gives K ph (; 0 ) = 2 N(0) Z 2 F ph (!) I(; 0 ;!)I(; 0 ;!) tanh[(=2)] tanh[(=2) 0 ] d!; (3.57) where the phononic Eliashberg function now enters. In this form the designation ofK pl should make sense as it exactly matches the form ofK ph with only the proper Eliashberg function substituted. Considering also the static Coulomb kernel,K elel , we get K elel (; 0 ) = 1 N()N( 0 ) X k;k 0 ( k )( k 0)K elel k;k 0 = 1 N()N( 0 ) X k;k 0 ( k )( k 0)W k;k 0(! = 0): If we again assume only electrons close to the Fermi-level contribute and recall the denition of the static Coulomb parameter [49] traditionally used in Eliashberg theory, = 1 N(0) X k;k 0 ( k )( k 0)W k;k 0(! = 0) we get K elel (; 0 ) = N(0) (3.58) which is constant (does not depend on or 0 ). Here it is interesting to note that the Coulomb interaction in fact contributes two constants to the full interaction kernel, one coming fromK elel and the other from the pl =N(0) term ofK elpl . Lastly, we also show the energy-space versions of the diagonal parts of the exchange correlation kernel, for the sake of completeness Z ph () = Z 2 F ph (!) Z N( 0 ) N(0) J(; 0 ;!) +J(; 0 ;!) tanh[(=2)] d 0 d! (3.59) and Z pl () = Z 2 F pl (!) Z N( 0 ) N(0) J(; 0 ;!) +J(; 0 ;!) tanh[(=2)] d 0 d! (3.60) where again the only dierence is in substituting the proper Eliashberg function. The gap equation can now be written in the much cleaner form () =Z()() 1 2 Z 1 1 N( 0 )K(; 0 ) tanh[(=2) 0 ] 0 ()d 0 : (3.61) 62 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS 3.2.3 Derivation of McMillan's equation in SC-DFT It is simple to show that the celebrated McMillan equation, T c = ! ln 1:20 exp 1:04(1 +) (1 + 0:62) (3.62) commonly used to approximate the Eliashberg solution for the critical temperature, T c , can be derived within SC-DFT. Here we will follow much the same steps as McMillan himself took in deriving the form of his equation in his publication, Ref. [112]. For simplicity we assume the kernel consists of two interactions, 1) the static, repulsive Coulomb interaction with interaction strength given by the parameter, and 2) an attractive phonon mediated interaction with strength . The Coulomb interaction is present for all electron pairs whereas the phonon mediated interaction only happens for electrons within a energy window given by the characteristic phonon frequency, we'll call it ! ln (to match with standard Eliashberg notation). Hence we assume the kernel can be decomposed into two parts, K(; 0 ) =K stat +K dyn (; 0 ) with the following forms for each part: K stat = N 0 (3.63a) K dyn (; 0 ) = 8 > < > : =N 0 ; ifjj;j 0 j! ln =N 0 ! ln ( 0 +! ln ) 1 tanh(=2 0 ) ; jj! ln ;j 0 j>! ln 0; otherwise (3.63b) We note here that the inclusion of the coupling between one electron inside thejj! ln window with another electron outside the window is very important in recovering the exact form of the McMillan equation. Without this the following derivation will result in the standard BCS-like equation for T c . In order to justify the form of this high-energy tail, we show in Fig. 3.2 a comparison between the approximated kernel above with one calculated from the phonon kernel in Eq. 3.57. The full numerically calculated kernel is for a Einstein-model of phonons (see Section 4.1.2) with phonon frequency 10meV and uniform electron-phonon couplingg 2 = 50meV. In order to nd a solution for the superconducting gap, we plug a trial function for () into Eq. 3.61, perform the integral to compute () and require the two to be consistent. In accordance with McMillan's strategy, we choose the trial function, () = ( 0 ; jj! ln 1 ; jj>! ln (3.64) 63 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY 10 −6 10 −5 10 −4 10 −3 10 −2 10 −1 10 0 10 1 ξ ′ −10 −8 −6 −4 −2 0 K(ξ = 0,ξ ′ ) Approximate Kernel Function Full Kernel, ω E = 10meV Approx Kernel Figure 3.2: Full numerical kernel and kernel approximated by the form given in Eq. 3.63b. 64 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS By plugging Eq. 3.63 and Eq. 3.64 into Eq. 3.61 we nd the necessary conditions on 0 and 1 . There are three contributions to 0 , (0) 0 =Z(0) 0 (1) 0 = Z ! ln 0 N 0 K(! ln ; 0 ) tanh(=2 0 ) 0 0 d 0 (2) 0 = Z 1 ! ln N 0 K(! ln ; 0 ) tanh(=2 0 ) 0 1 d 0 and two contributions to 1 , (0) 1 = Z ! ln 0 N 0 K(>! ln ; 0 ) tanh(=2 0 ) 0 0 d 0 (1) 1 = Z 1 ! ln N 0 K(>! ln ; 0 ) tanh(=2 0 ) 0 1 d 0 : Lets rst consider the 1 terms since they are simpler, (0) 1 = Z ! ln 0 N 0 ( N 0 ) tanh(=2 0 ) 0 0 d 0 = 0 Z ! ln 0 tanh(=2 0 ) 0 d 0 (1) 1 = Z 1 ! ln N 0 ( N 0 ) tanh(=2 0 ) 0 0 d 0 = 0 Z 1 ! ln tanh(=2 0 ) 0 d 0 : We dene the following for the sake of simplicity, p = Z ! ln 0 tanh(=2 0 ) 0 d 0 q = Z 1 ! ln tanh(=2 0 ) 0 d 0 ; which gives 1 =p 0 q 1 : (3.65) Now, we return to the 0 terms. The (1) 0 is straight forwardly (1) 0 =p 0 +p 0 ; while the nal term is given by, (2) 0 = Z 1 ! ln N 0 K(! ln ; 0 ) tanh(=2 0 ) 0 1 d 0 = Z 1 ! ln N 0 N 0 N 0 ! ln ( 0 +! ln ) 1 tanh(=2 0 ) tanh(=2 0 ) 0 1 d 0 =q 1 + 1 ! ln Z 1 ! ln ( 0 +! ln ) 1 d 0 0 65 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY McMillan replaces the integral above by introducing an average boson frequencyh!i where R 1 ! ln ( 0 +! ln ) 1d 0 0 =h!i=! ln . However, the integral can be easily evaluated and we nd, ! ln Z 1 ! ln ( 0 +! ln ) 1 d 0 0 = ln(2): Putting all the 0 terms together we get 0 =Z(0) 0 p 0 +p 0 q 1 + ln(2) 1 : (3.66) Combining Eq. 3.66 and 3.65 we get a system of equations given by 0 1 = Z(0)p +p q + ln(2) p q 0 1 (3.67) We are interested in the critical temperature at which the Eq. 3.67 has a non-trivial solution or equivalently, the temperature at which the above matrix (M) has an eigenvalue equal to 1. This is found where det(1 M) = 0: After some simple algebra we nd p = (1 +Z(0)) ( 1 +q (1 + ln(2))) 1 : (3.68) Now we return to the expressions forp andq. Integrals of this sort were solved in the BCS paper [15] and the given solutions are p = ln ! ln 2k B T c + ln 4e q = ln E B ! ln (3.69) whereE B indicates the highest band-energy in the system, a natural cuto for the integra- tion. Interestingly, we note that with the proper substitution of q we recover the Coulomb potential of Morel and Anderson [113], commonly used in Eliashberg theory, 1 +q = 1 + ln E B ! ln = Finally, we get an expression for T c by solving Eq. 3.68 using the proper value of p, k B T c = 1:13! ln exp 1 + (1 + ln(2)) : (3.70) This equation has exactly the same form as the celebrated McMillan equation, proving that the SC-DFT theory contains all the parts necessary to reproduce any eect seen from it has been demonstrated by Luders et. al. [102] that Z(0) 66 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS the McMillan equation (Eq. 3.62). We are able to do more though. Knowing the form of the electron-plasmon interaction kernel given in Eq. 3.55, we see that the inclusion of plasmons in a superconducting system contributes in 2 ways: 1. the static term inK elpl can be seen to increase the static Coulomb repulsion i.e. mapping ! tot = + pl 2. the diagonal contributionZ pl and non-diagonal contributionK pl both share the exact same form as the phononic kernel contributions, thus the plasmons increase the overall coupling strength i.e. mapping ! tot = ph + pl In this simple way we can obtain an equation for the critical temperature of a supercon- ducting system with both phonon and plasmon interactions included: k B T c = 1:13! tot ln exp 1 + tot tot tot (1 + ln(2) tot ) (3.71) where ! tot ln = exp 2 tot Z 2 F ph (!) + 2 F pl (!) ! log(!)d! : We can explore the consequences of the above equation by assuming Einstein modes for the phonons and plasmons i.e. the bosons exist at specic frequencies ! ph and ! pl . This would simplify ! tot ln to ! tot ln = exp h ( ph log(! ph ) + pl log(! pl ))= tot i : The equation constructed above should only be considered for the sake of qualitative argu- ments. A quantitative analysis would require solving the gap equation for an actual system considering the proper forms of the interaction kernels. However, we can consider whether the inclusion of plasmons lead to an increased or decreased critical temperature compared to the situation where plasmons are neglected from the system. In order to do this we x the frequencies of the two modes, vary the coupling strengths and consider the ratio, R, of the resulting T c to the T c obtained when pl = 0 i.e. R = ! tot ln ! ph ln exp 1 + ph ph ph 1 + tot tot tot < 1: (3.72) The results for this are shown in Fig. 3.3. Interestingly we notice that the inclusion of plasmons to a system already containing phonons can in some instances increase T c while in others decrease it. This interplay has not, to my knowledge, been seen by any other authors although there have been publications showing either an increased T c [114] when 67 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY Figure 3.3: Plot showing the T c ratios calculated, using Eq. 3.71, between a system including plasmons and phonons and a system including only phonons. We take ! pl = 200meV and! ph = 2meV. Note the regions where T c is increased (ratio > 1) by plasmons and those where T c is decreased (ratio < 1). including plasmons or a decreasedT c , but so far not both possibilities. The increase inT c is somewhat expected since another bosonic mode allows stronger coupling between Cooper pairs. The plasmonic reduction of the full critical temperature is less obvious but can be understood as follows: Within the chosen parameters the ratio ! tot ! ph > 1 can not be responsible for R being smaller than 1. Therefore, the argument of the exponential function must me smaller than 0, 1 + ph ph ph 1 + tot tot tot < 0: (3.73) In the limit ph pl we can approximate tot ph and can reformulate the statement to 1 + ph ph ph < 1 + ph ph tot ; (3.74) which is true for tot > ph : (3.75) 68 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS This, in turn is equivalent to ph + pl 1 + ln(E F =! tot ) > ph 1 + ln(E F =! ph ) (3.76) , ph + pl ph | {z } &1 > 1 + ln(E F =! tot ) 1 + ln(E F =! ph ) | {z } <1 since ! tot >! ph : (3.77) tot is thus slightly enhanced in comparison to ph due to pl and the increased eective frequency ! tot . The plasmonic reduction for enhanced phononic coupling is therefore a result of the plasmonic enhancement of the pseudo Coulomb potential tot , whereas the plasmonic contribution to the total eective coupling becomes negligible and tot ph . On the other side, if the phononic coupling is weak, the plasmonic contribution to the total eective coupling becomes relevant. Here, the eective coupling can overcome the eec- tive Coulomb potential leading to a net enhancement of T c in comparison to T ph+ c . It is thus a very subtle interplay between all interactions which controls the resulting transition temperature. Finally, we say a few words regarding the numerics of nding T c for a given system. In principle the technique is to self-consistently solve the gap equation (Eq. 3.61) and nd the lowest temperature at which the only solution is the trivial solution, () = 0. This is computationally a very expensive exercise seeing as each iteration of the gap equation requires summing over a large number of energy points. A dierent approach is to recast the gap equation as an eigenvalue problem [115] K = (3.78) using the generalized kernel matrix K ; 0 = ( Z() 1 2 N()K(;) tanh[(=2)] = 0 1 2 N( 0 )K(; 0 ) tanh[(=2) 0 ] 0 6= 0 : Using this we can ndT c as the temperature at which the leading eigenvalue is = 1. This method also has an advantage over the iterative method in that we can make an educated guess for which temperature to consider next by using Newton's method. 3.2.4 Plasmonic Eliashberg function in 2d Now we turn to the more complicated case where the simplifying assumption made previ- ously concerning 2 F (!) (from Eq. 3.53) fails. That is, the case where it is not permittable to only consider 2 F (!) = 2 F ( = 0; 0 = 0;!) in the evaluation of the exchange- correlation kernel. This is precisely the situation for the electron-plasmon interaction in a 69 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY 2d system, since the coupling is given by (see Eq. 2.133) a q = U(q)! p (q) 2 ; where ! p (q) = (v 2 F N 0 e 2 ) 1=2 p q and U(q) = e 2 =q, which gives a q !1 as q! 0. In such a situation 2 F ( = 0; 0 = 0;!! 0) diverges and it is necessary to consider the full 2 F (; 0 ;!) function in the exchange-correlation kernel. Therefore, we will discuss, here, the small q behavior of the full Eliashberg function to show that this is well behaved. As given in Eq. 3.53, the full Eliashberg function is dened as 2 F pl (; 0 ;!) = N(0) N()N( 0 ) X k;kq a k;kq ( k )( 0 kq )(!! k;kq ): If we consider the case wherejqj is small we generally have kq k v F q ? , where v F is the Fermi velocity and q ? is the component of q perpendicular to k. Plugging this and the expressions for a q and ! p (q), form above, back into the denition of the Eliashberg function we get, 2 F pl (; 0 ;!) = N(0) N()N( 0 ) X k;kq (v 2 F N 0 ) 1=2 e 3 2 p q ( k )( 0 k v F q ? )(! (v 2 F N 0 e 2 ) 1=2 p q) = N(0) N( 0 ) X q (v 2 F N 0 ) 1=2 e 3 2 p q ( 0 v F q ? )(! (v 2 F N 0 e 2 ) 1=2 p q) = N(0) N( 0 ) Z Z ( 0 v F q ? ) ! (v 2 F N 0 e 2 ) 1=2 (q 2 ? +q 2 k ) 1=4 (v 2 F N 0 ) 1=2 e 3 2 q (q 2 ? +q 2 k ) 1=2 dq k dq ? = N(0) N( 0 ) Z 0 @ ! (v 2 F N 0 e 2 ) 1=2 " 0 v F 2 +q 2 k # 1=4 1 A N 1=2 0 e 3 2 0 v F 2 +q 2 k 1=4 dq k (3.79) where the integration over q ? and the rst function gives q ? = ( 0 )=v F . In order to evaluate the q k integral, recall for a function g(x) we have (g(x)) = X i (xx i ) jg 0 (x i )j where x i is the i'th root of the function g(x) and g 0 (x) is it's rst derivative. Take g(q k ) =! (v 2 F N 0 e 2 ) 1=2 " 0 v F 2 +q 2 k # 1=4 ; this expression, for the general 2d plasmon dispersion, comes from Ref. [31], taking the small q limit 70 CHAPTER 3. SUPERCONDUCTIVITY 3.2. DFT FOR SUPERCONDUCTORS then we have g 0 (q k ) = (v 2 F N 0 e 2 ) 1=2 q k 2 0 v F 2 +q 2 k 3=4 ; with the roots of g(q k ) given by q i = s ! 4 (v 2 F N 0 e 2 ) 2 0 v F 2 : Plugging this back into Eq. 3.79 and realizing that ! < p j 0 j(v F N 0 e 2 ) results in an unphysical imaginary component of 2 F , yields 2 F pl (; 0 ;!) = 8 > > < > > : N(0) N( 0 ) e 2 v F 1 q 1 ! min ! 4 ; !>! min 0; otherwise (3.80) where ! 2 min j 0 j(v F N 0 e 2 ): 0 1 2 3 4 5 ω α 2 F(ξ = 0,ξ ′ ,ω) Energy dependent plasmonic Eliashberg function ξ ′ = 0.0 ξ ′ = 1.0 ξ ′ = 2.0 Figure 3.4: Full energy dependent plasmonic Eliashberg function for a 2d system, using = 0 and varying 0 . This function is shown in Fig. 3.4 for = 0 and dierent values of 0 (assuming a constant DOS and setting all other constants to 1). The function is gapped up to ! min at which point it diverges logarithmically (from the right). For larger ! the Eliashberg 71 3.2. DFT FOR SUPERCONDUCTORS CHAPTER 3. SUPERCONDUCTIVITY function approaches a constant. The full energy dependent Eliashberg function still shows a divergence at! =! min , however, this logarithmic divergence disappears when integrating over ! in the evaluation of the kernel (Eq. 3.52). Specically, we can dene an energy dependent average coupling, (; 0 ), similarly as before, plas (; 0 ) = Z 2 F plas (; 0 ;!) ! d! = 1 2 log 0 @ q 1! 4 min + 1 ! 2 min 1 A ; which is nite for 6= 0 i.e. ! min > 0. On the diagonal of the kernel the divergence is cancelled by the correspondingly divergent attractive interaction. Resultingly the gap function is well behaved and a critical temperature can be calculated in the same way as discussed before. 72 Chapter 4 Model System - 2d Square Lattice In this chapter we will explore SC-DFT and it's extension to include plasmonic eects using a model system, in the form of the 2d square lattice. In particular, we show how plasmonic excitations can be precisely tuned by changing the dielectric environment or doping level of a layered metal, and under which circumstances these excitations strongly couple to the electrons. This allows us to easily isolate the eects specically coming from the electron- plasmon coupling since we can precisely control the electronic properties of the system. Furthermore, we use a simple model for the phonon contribution to the superconducting system in order to not cloud the results with phenomenon stemming from electron-phonon coupling. In the rst section of this chapter we will discuss the model system and it's electronic properties and in the second we'll outline the model phonons. The third section will be devoted to discussing plasmon properties and how we can manipulate the electron- plasmon coupling and plasmon dispersion. Finally, we'll discuss how the superconducting critical temperature of the system is aected by all the aforementioned. We demonstrate that the critical temperature can be enhanced as well as reduced due the electron-plasmon coupling, as we demonstrated analytically in Section 3.2.3. This chapter is entirely based on our publication - Ref. [116]. 4.1 Model Description We use the eective heterostructure dened by a two-dimensional metal sandwiched be- tween substrates, with tuneable dielectric constants, to study plasmonic properties. The metal is dened by a square lattice (see Fig. 4.1) since close to half-lling the square lattice electrons behave very similar to the free-electron gas, and therefore the system is well suited to isolate plasmonic eects. Furthermore, in order to obtain a realistic description of the Coulomb interactions, we imagine the metallic band to be part of a multi-band structure, e.g. formed by s, p, and d orbitals. While we will entirely concentrate on the low-energy subspace around the Fermi level, it is important to realize that the neglected high-energy 73 4.1. MODEL DESCRIPTION CHAPTER 4. SQUARE LATTICE t U(r 1 ,r 2 ) Figure 4.1: Square lattice with lattice constant a 0 is shown. We indicated by green line the nearest neighbor hopping parameter, t. The Coulomb interaction, U(r 1 ;r 2 ), is present between all sites and is shown with a blue line. The reciprocal lattice (also a square lattice) is shown in red. parts of the band structure have a screening in uence on the Coulomb interaction in the low-energy subspace [117]. If the metallic band is the only band crossing the Fermi level, we can readily approximate the polarization function of the neglected or remaining part of the band structure as a static function rest q . The total polarization function then reads total q (!) = q (!) + rest q , and the full dynamic Coulomb interaction is given by (recall Eq. 2.118) W q (!) = v q 1v q q (!) + rest q = U q 1U q q (!) = U q " q (!) ; (4.1) where v q is the bare interaction, and U q the background-screened interaction, dened by U q = v q 1v q rest q = v q " rest q : (4.2) The background dielectric function" rest (q) now renders all screening eects resulting from the neglected part of the band structure and also those resulting from the dielectric envi- ronment. This function can be derived from classical electrostatics. In the case of a layered system with thicknessd and embedded between two semi-innite dielectric substrates with dielectric constants " env it reads [118,119] " rest (q) =" 1 1 ~ " 2 e 2qd 1 + 2~ "e qd + ~ " 2 e 2qd ; (4.3) with ~ " = " 1 " env " 1 +" env : (4.4) 74 CHAPTER 4. SQUARE LATTICE 4.1. MODEL DESCRIPTION This function smoothly interpolates between the long-wavelength limit " rest (q! 0) =" env and the short-wavelength limit " rest (q!1) = " 1 of the layered system. Thus, U q is mostly aected by" env for small momenta, while the internal" 1 controls short wavelengths. If we additionally dene the bare interaction as v q = 2e 2 =[A(q + q 2 )], whereA is the unit cell size and is an eective form-factor rendering eects from the non-zero height of the layer, the realistic background-screened Coulomb interaction is fully dened by the parameters A, , " 1 , and d. Now, we choose the parameters of the model such that the system's energy scale (band- with and Coulomb interaction strength) are comparable to the corresponding values found in transition metal dichalcogenides [120,121]. To this end we choose the nearest neighbor electron hoppingt = 0:125 eV with lattice constanta 0 = 3 A. We take the eective thickness d = 5 A and " 1 = 25 (used in Eq. 4.3). This set of parameters results in a model band width of D = 1 eV and a background-screened on-site Coulomb interaction of U r=0 = 1= R IBZ dqU q 0:8 eV (at" env = 1), with being the area of the rst Brillouin zone (IBZ) and U q as dened in Eq. 4.2. This is similar to the situation in metallic TMDC's [122], although the involved Coulomb interaction is relatively reduced in order to get a moderately correlated model system [123,124]. This simple model allows us to study the two important eects. Firstly, we can investi- gate how the plasmon dispersion is changed by changing the dielectric environment around the metal. We achieve this by simply changing the value of" env entering Eq. 4.4. Studying this in uence is important since screening in 2d is highly aected by the environment, as shown in Fig. 4.2. For particles separated by large distances the electric eld lines of the Coulomb interaction between them, lie entirely outside of the 2d layer, therefore the mate- rial itself screens this interaction very lightly. It is due to this that 2d materials experience a much stronger Coulomb interaction than 3d materials. It also means that we are able to control the screening between charges in layered materials to a great extent. This ability is in big part the motivation for the work done in this thesis. ! Figure 4.2: Sketch showing how the majority of electric eld lines connecting charges in a 2d material extend into the space surrounding the layer. Adding a dielectric material above or below the 2d sheet thus greatly aects the Coulomb interaction between charges in the 2d layer. 75 4.1. MODEL DESCRIPTION CHAPTER 4. SQUARE LATTICE Secondly, if we treat doping as a rigid shift of the system's Fermi-level, we can also study how the plasmon dispersion is changed when we change the number of electrons in the system. This is useful since the superconducting transition temperature is greatly aected by the electronic density of states at the Fermi level (see Eq. 3.11). 4.1.1 Electronic Dispersion and DOS The tight-binding method was discussed in Section 2.2.1. The Hamiltonian is given in real-space, re ecting the hopping from one atomic site to another. We nd the electronic dispersion (in k-space) by diagonalizing the Fourier transform of the tight-binding Hamil- tonian. The density of states is calculated in the usual way N(E) = 1 Z (EE k )dk = 1 jkj X (EE k ); where is the area of the rst Brillouin zone andjkj is meant to be the cardinality of the set of k-points considered. Γ X M Γ −0.6 −0.4 −0.2 0.0 0.2 0.4 0.6 ξ k (eV) Band Structure DOS(ξ) Figure 4.3: Square lattice electronic band structure and density of states shown for the half-lling case (red dotted line shows Fermi-level). The inset shows the path taken through the BZ. The band structure of the model is shown in Fig. 4.3. In plotting the band structure we follow a path crossing points of high symmetry, !X!M! (shown in inset of Fig. 76 CHAPTER 4. SQUARE LATTICE 4.1. MODEL DESCRIPTION 4.3). At half-lling the model shows a Van-Hove singularity characterized by a saddle point in the band structure and a diverging density of states. Due to this divergence we avoid the exact half-lling scenario and instead we will study variable doping levels, ranging from 0.6 to 0.9 lling. 4.1.2 Phononic Properties In order to keep the eects stemming from electron-phonon coupling as simple as possible we assume an Einstein model for the system phonons [125]. In this model we only include optical phonons and assume the phonon dispersion is constant, meaning kk 0 =! E 8k; k 0 . The phonon density of states is thus given by N ph (E) =d(E! E ); where d is the number of phonons in the system. Furthermore, we also assume constant coupling between electrons and phonons which means we set jg kk 0 q j 2 =g 2 ; with g 2 being a parameter of our model. The phononic Eliashberg function is given by 2 F ph (!) = 1 N(0) X k;k 0 jg kk 0 q j 2 ( k )( k 0)(! kk 0); (4.5) where N(0) is the electronic density of states at the Fermi-level. If we plug the Einstein model into Eq. 4.5 we get 2 F ph (!) = 1 N(0) X k;k 0 g 2 ( k )( k 0)(!! E ) =g 2 N(0)(!! E ): (4.6) The electron-phonon coupling parameter ph can also be found, ph = 2 Z 2 F ph (!) ! d! = 2g 2 N(0) ! E : (4.7) We are thus able to tune the value of ph directly by changing theg 2 parameter in our model system. Here it is important to note that the superconducting critical temperature for the phonon only case, depends on ph and the density of states at the Fermi-level. Thus, even our simple model to include phonons in the system will have non-trivial behavior when we change doping levels. 77 4.2. PLASMONIC PROPERTIES CHAPTER 4. SQUARE LATTICE Figure 4.4: Plots showing the polarization function, q (!), for q-points on a path from to K for various doping levels and environmental screening scenarios. The extracted plasmon dispersions are shown as dashed, red lines. N 0 is the density of states at the Fermi-level. Notice the clear appearance of the particle-hole continuum: negative (blue) area of the real part of the polarization function. 4.2 Plasmonic Properties As we discussed in Section 2.3.1 the plasmonic properties of the system are completely captured by " q (!), given in Eq. 2.119 and copied here for ease of reading, " q (!) = 1U q q (!); whereU q is the background-screened Coulomb interaction given in Eq. 4.2 and q (!) is the polarization function given in Eq. 2.120. Once we have the eigenvalues and eigenvectors of the tight-binding Hamiltonian, we have everything we need to calculate the polarization function and consequently the dynamical, dielectric screening function. In Fig. 4.4 we show the polarization function calculated for doping levels of 0:6, 0:75 and 0:9 lling. We show results for various environmental screening situations in which we change " env (used in Eq. 4.4) from 1 to 5 to 25. The plots are focused on smalljqj values where the plasmon dispersion is properly separated from the particle-hole continuum. The particle- hole continuum arises due to electrons on (or just below) the Fermi-level being scattered to states immediately above the Fermi-level i.e. particle-hole pair creation and annihilation 78 CHAPTER 4. SQUARE LATTICE 4.2. PLASMONIC PROPERTIES events. If we assume a isotropic Fermi-surface, a perfect circle, the particle-hole excitations are scattering processes that connect two points on the Fermi-surface. It is easy to see that such excitations can occur with 0 energy separation for any momentum transfer, q, such that 0 q 2k F where k F is the radius of the Fermi-surface. For !6= 0 the range of possible momenta transfers rigidly shifts by an amount, q 0 , dictated by the eective mass, m , of electrons in the system such that ! = q 2 0 2m . This is assuming a low doping level so that the energy dispersion is still quadratic. In our case the doping-levels considered are high enough that the electrons disperse linearly (see Fig. 4.3) which is why the particle-hole continuum seen in Fig. 4.4 increases linearly with !. It is important to note the distinction between single-particle excitations that lead to the particle-hole continuum and collective excitations, namely plasmons. A crucial dierence is that the single particle excitations are not directly aected by electron-electron interactions and therefore are not aected by environmental screening. This we clearly see in Fig. 4.4, where the particle-hole continuum is unchanged when varying " env . The same is not true for the plasmon dispersion as we discuss next. 4.2.1 Plasmon Dispersion As we have discussed already, plasmonic excitations are dened as collective oscillations of the Fermi-sea that comes about due to a diverging Coulomb interaction, leading to absolute electron-electron correlation. This divergence of the Coulomb interaction is the result of vanishing dielectric screening i.e. " q (!) = 0; which can be seen from Eq. 4.1 would cause an innitely strong interaction. The (q;!) pairs for which this occurs correspond to maxima of the EELS spectrum, given by Eq. 2.78. In Fig. 4.5 we show the EELS spectra for dierent environmental screenings and dierent doping-levels. We connect the maxima in the spectrum with red dotted lines showing the plasmon dispersion, ! pl (q). In all cases we nd gapless modes with p q-like dispersions for small momenta, which atten at intermediate q before smoothly hybridizing with the particle-hole continuum, which is characteristic for plasmons in 2d metals [21, 25, 126]. Most importantly, Fig. 4.5 illustrates the sensitivity of the plasmon dispersion to the layer's environment and to the doping level. Increasing the environmental screening from " env = 1 to 25 decreases the dispersion for q. 0:25 A 1 only, which directly results from the decreased Coulomb interaction U q in this range. A similar eect has been predicted for heterostructures consisting of graphene and hexagonal boron nitride [127]. Increasing the electron-doping level from n 0:6 to n 0:9 induces a reduction of the plasmon frequencies for the entire q-range. In order to understand why the plasmon dispersion is aected by the doping we show here the approximate solution for " q (!) = 0 given by the Coulomb interaction plays an indirect role since it aects the spectrum of the Kohn-Sham Hamil- tonian 79 4.2. PLASMONIC PROPERTIES CHAPTER 4. SQUARE LATTICE Figure 4.5: Plots showing the electron energy loss spectrum (EELS) for dierent" env values and dierent doping levels. The red lines trace the maxima of the EELS, hence showing the plasmon dispersion while the blue, dotted lines show that upper limit of the particle- hole continuum. Notice the strong dependence on doping and substrate screening of the plasmon dispersion while the particle hole continuum is unaected by these factors. Ref. [31]. They nd ! pl (q) =qv F s 1 + (N 0 U q ) 2 0:25 +N 0 U q ; (4.8) thus, decreasing the density of states at the Fermi level (N 0 ) directly in uences the plas- monic dispersion. We would like to stress the important result seen here, namely, since the plasmons originate due to the Coulomb interaction we can manipulate the plasmon dispersion by altering the environmental screening surrounding the metallic layer. This provides us with a tuning knob in the study of plasmonics and plasmonic superconductivity since any eect due to electron-plasmon coupling would also be aected by environmental screening. 4.2.2 Plasmonic Eliashberg Function In order to nd the plasmonic Eliashberg function we rst have to nd the electron-plasmon coupling strength (discussed in Section 2.4.2.1). Once we have the plasmon dispersion we 80 CHAPTER 4. SQUARE LATTICE 4.2. PLASMONIC PROPERTIES −4 −2 0 2 4 ε env = 1 ε env = 5 ε env = 25 n≈ 0.60 N 0 ≈ 1.5eV −1 −4 −2 0 2 4 log(a q ) n≈ 0.75 N 0 ≈ 1.0eV −1 0.00 0.15 0.30 −4 −2 0 2 4 0.00 0.15 0.30 q ( ˚ A −1 ) 0.00 0.15 0.30 n≈ 0.90 N 0 ≈ 0.8eV −1 Figure 4.6: Plots showing the log of extracted electron-plasmon coupling, a q , for various doping levels and environmental screening scenarios, in red, dashed line. The expected, analytical coupling for a free electron gas with the same plasmon dispersion, is shown with solid, blue lines. can nd the electron-plasmon coupling, a q , using Eq. 2.127 (copied here for reading ease) a q = ! pl q 2 [W q (0)W q (1)]: In Fig. 4.6 we show the electron-plasmon coupling for the same environmental screening and doping cases used in Fig. 4.4 and 4.5. We also show the expected coupling associated with a free electron gas with the same plasmon dispersion as our model, given by Eq. 2.133 (again, copied here for reading ease) a q = U(q)! pl q 2 : For q-points where the plasmon dispersion is signicantly removed from the particle-hole continuum the extracted coupling matches our expectation from the free electron gas very well. As soon as the plasmon frequencies get closer to the continuum, the tted and analytic values fora q start to dier and the tted electron-plasmon coupling becomes smaller than 81 4.2. PLASMONIC PROPERTIES CHAPTER 4. SQUARE LATTICE the analytical results. As soon as we have the coupling strengths we can calculate the plasmonic Eliashberg function (given in Eq. 3.53) 2 F pl (; 0 ;!) N(0) N()N( 0 ) X k;k 0 a k;k 0( k )( 0 k 0)(!! k;k 0); which describes the average eciency of plasmons with frequency ! to scatter electrons from state to state 0 . As discussed in Sec. 3.2.4, for small energy dierences between and 0 the plasmonic 2 F function diverges at a specic ! value. This means converging the numerically cal- culated 2 F function requires very high resolution momentum space calculations. In fact, such high resolution that it is not feasible. Fortunately, we have access to the analytic form (Eq. 3.80), which is valid for smallj 0 j and small !. On the other hand, high energy transitions are accurately described within the numerical calculations. We can thus dene a cross over point below which we have condence in the analytic form and above which we can trust the numerical solution of 2 F (; 0 ;!). This crossing point is determined by the smallest non-zero plasmon frequency coming from the q-grid. We call this! q min and it sets the crossing point on the !-axis. Thej 0 j value at which the analytic solution is no longer trustable is determined by the validity of approximating kq k v F q ? , for which it is harder to determine an exact number. Trial and error shows that a cuto value ofj 0 j = 0:02eV is a good mark. With these boundaries in mind we establish a region around them in which we linearly mix the analytic and numerical solutions based on how close we are to the respective sides of the boundary. This ensures a smooth 'stitching' of the two solutions and maximizes accuracy. For the mixing region we use 1 2 ! q min <! < 3 2 ! q min on the !-axis and 0<j 0 j< 0:02 on the other. The above described stitching routine allows us to deal with the numerical instabilities arising from having a nite minimum value of q in our grids. However, handling the calculation of 2 F pl (; 0 ;!) for, on the order of, 10 6 pairs of (; 0 ) remains a techinical diculty which we are working on. In the interim as well as for the sake of simpler comparison to the analytic results from Sec. 3.2.3 we now make the following non-physical simplication: we assume the electron-plasmon coupling is zero forjqj values smaller than the norm of the minimum q in our numerical grid. Explicitly, we make the simplication a q = 0 if jqj< 0:02 A where 0:02 A is the smallest detail we can resolve on our numerical grids. IMPORTANT: This simplication is not physical and reduces what follows to a purely theoretical exer- cise. The benet is that the results can easily be understood in the regular language of superconductivity theory and the hope is that the results will provide a stepwise insight into the eect plasmons have on superconductivity. This simplication is not made in our journal publication concerning the same topic. Within this severe simplication, we are 82 CHAPTER 4. SQUARE LATTICE 4.2. PLASMONIC PROPERTIES Figure 4.7: Plasmonic Eliashberg functions, 2 F pl (w), for dierent j 0 j situations stagered. The blue line marks ! q min and grey line indicate the boundaries of the region in which the analytic and numerical solutions are mixed. able to neglect the and 0 dependencies in plasmonic Eliashberg function but instead use 2 F pl (!) = 2 F pl ( = 0; 0 = 0;!). In Fig. 4.8 we show the 2 F pl (!) function corresponding to dierent environmental screenings " env and doping-levels n. In the case of varying " env (left panel) we nd the peaks in 2 F pl (!) shifting to lower energies and decreasing in amplitude with increasing environmental screening (decreasing 1=" env ). This shift in energy results from the gen- eral decrease of the plasmon frequencies with increasing screening as already observed in Fig. 4.5. The strongly decreasing peak height originates from the electron-plasmon cou- pling which is approximately given by a pl q / ! pl q U q . The coupling is thus reduced due to the decreasing frequencies! pl q and, most importantly, due toU q , which is strongly reduced by " env . These trends are also observed in the eective electron-plasmon coupling pl and in the eective plasmon energy! pl log , which both decrease with increasing" env (we'll return to this later). 83 4.2. PLASMONIC PROPERTIES CHAPTER 4. SQUARE LATTICE 0 5 10 15 α 2 F pl (ω) n≈ 0.8 0 1 1/ε env ε env = 1 1.8 0.7 N 0 0.0 0.1 0.2 0.3 ω (eV) 0.0 2.5 5.0 7.5 F pl (ω) 0.0 0.1 0.2 0.3 ω (eV) Figure 4.8: Plasmonic Eliashberg functions for dierent screening environments (top left) and dierent doping levels (top right), with bottom plots showing the plasmonic density of states. Notice the strong increase in plasmon numbers for higher electron doping (bottom right). This strong increase in plasmonic DOS explains the roughly constant height of 2 F pl seen in the top right panel. For the case of varying doping level n and thus changing N 0 (top right panel in Fig. 4.8) we primarily nd a shift of 2 F pl (!) to higher energies when increasing N 0 , with- out any major changes in its spectral weight. The shift is again explained by the raised plasmon frequencies. At the same time the electron-plasmon coupling a q is increased due to its ! pl q dependence. Thus, a q actually increases the peak height of 2 F pl (!) for higher N 0 . However, the pure plasmon density of state F pl (!) = P q (! pl q !) shifts to higher frequencies and looses weight for small frequencies (see bottom right panel Fig. 4.8) which compensates for the enhancement from a q for larger N 0 to 2 F pl (!) yielding just slightly elevated pl for comparably large N 0 . 4.2.3 Self-energy Contribution We check the in uence of the electron-plasmon interaction on the spectral function by using the anomalous electronic self energy in the Nambu space [128], given within the 84 CHAPTER 4. SQUARE LATTICE 4.2. PLASMONIC PROPERTIES G 0 W 0 approximation as k (i! n ) = 1 X qm G 0 q (i! m )W kq (i! n i! m ): (4.9) When we use a single plasmon pole t to the full dynamic Coulomb interaction (see Eq. 2.126) we can split the above into a static and dynamic part, where stat k = X q N( q ) " W q (0) + 2a pl q ! pl q # ; and dyn k (i! n ) = 1 X k 0 m G 0 k (i! m )a pl q 2! pl q (i! n i! m ) 2 (! pl q ) 2 : (4.10) We are only really interested in the dynamic part, and we can evaluate the Matsubara sum analytically [92], yielding dyn k (!) = Z BZ dqa pl q (4.11) " n q +f k+q ! k+q +! pl q +i + + n q + 1f k+q ! k+q ! pl q +i + # : This corresponds to the self-energy given in Ref. [100]. From the electron-plasmon self- energy we calculate the interacting spectral functions (see Eq. 2.139), A(k;!) = 1 jIm dyn k (!)j [! k Re dyn k (!)] 2 + [Im dyn k (!)] 2 ; and interacting DOS A(!) = X k A(k;!); as shown for dierent scenarios in Fig. 4.9. These calculations were performed on 80 80 k-meshes and equivalent q-meshes using Gaussian functions instead of Dirac delta functions (for the energies) with 15meV broadening and with a value of 0:3eV for + . We observe a reduction of the band-width in combination with arising plasmonic satel- lites above and below the non-interacting band edges. At the same time the former van- Hove singularity shifts to higher energy and, most important, the density of states at the Fermi level is enhanced. The latter is equivalent to an eective electron mass enhance- ment, which is why we decided to take the plasmonic mass enhancement factorZ pl (!) into 85 4.2. PLASMONIC PROPERTIES CHAPTER 4. SQUARE LATTICE Figure 4.9: Interacting Spectral Functions. The doping (left to right) and substrate screen- ing (top to bottom) dependencies of the electron-plasmon renormalization of the spectral function, A(k;!), is shown along with the renormalized electronic DOS. The white lines show the non-interacting spectral function and red lines the non-interacting DOS. In the DOS plots the blue lines show the interacting DOS. 86 CHAPTER 4. SQUARE LATTICE 4.3. SUPERCONDUCTIVITY account . All of these eects are controlled by the dielectric environment and the doping level. If we compare this to the spectral ngerprints of plasmons in a 3d (semiconducting) silicon system (see Ref. [129]) we nd a dierent scenario. In 3d, due to the non-dispersive plasmonic modes at energies around 5 to 10 eV, there is just a replica of the original non- interacting band structure shifted to lower energies (in the amount of the non-dispersive plasmonic energy). This might also be seen as a plasmonic satellite. But, this satellite band structure is strongly reduced in its spectral weight since in 3d the coupling scales as ja pl q j 2 / 1=q 2 , while in 2d we haveja pl q j 2 / 1=q. In 2d we thus observe two major eects: (a) low-energy plasmonic modes which result in low-energy changes to the band structure around the Fermi level and (b) an enhanced coupling of the low energy modes, which enhance the low-energy spectral ngerprints. It is well known that the electron-phonon coupling renormalizes the electronic spectrum within the Debye window around the Fermi level. Thereby, it enhances the eective electron mass, which needs to be considered in form of the phononicZ ph () in the superconducting state. In the case of plasmonic excitations and their coupling to the electronic states a corresponding shift of spectral weight is also expected [123, 130], but has so far been neglected within SC-DFT treatments. Neglecting this plasmonic contribution to Z() is indeed a reasonable approximation for 3d dispersion-less plasmons with energies on the order of the electronic band width. However, for the case of 2d gapless plasmons there is a signicant eect on the electronic spectrum around the Fermi energy, and we need to take Z pl () into account. Otherwise, we would overestimate the plasmonic contributions to the critical temperature. 4.3 Superconductivity Based on the discussion in Section 3.2.3 we expect a decent match between the critical temperature we calculate using our SC-DFT implementation and the value one would get from McMillan's equation (Eq. 3.62). Thus as a rst step we compare exactly these things in Fig. 4.10, where we x the electron-phonon coupling and assume no Coulomb interac- tion in the system. The correspondence in trends gives us condence in the accuracy of our SC-DFT implementation as well as the correctness of the derivation in Section 3.2.3. For further benchmarks of our SC-DFT implementation see Appendix C. Knowing that Eq. 3.71 gives us a good handle to understand how the electron-phonon and electron-electron interactions in uence T c , we start by examining the behavior of the averaged coupling parameters as a function of environmental screening and doping. This is shown in Fig. 4.11 where the color of each data-point matches the color of the corresponding 2 F (!) curve from Fig. 4.8. The largest eective plasmonic couplings are found for small in Ref. [99] the authors dealt with 3d materials and therefore were able to safely neglectZ pl 87 4.3. SUPERCONDUCTIVITY CHAPTER 4. SQUARE LATTICE 0.04 0.08 0.12 0.16 ω E (eV) 10 20 30 40 T c (K) McMillan SCDFT Figure 4.10: We show the McMillanT c values for dierent! E 's and theT c calculated from our SC-DFT implementation (using a xed electron-phonon coupling and = 0). Notice the maximum value around ! E = 0:04. This is due to the interplay of / 1=! E while T c /! E exp[1=]. N 0 . This results from the constant shift of 2 F pl (!) to small frequencies with decreasingN 0 (which in turn results from a strongly enhanced F pl (!) for smallN 0 and small!) together with the denition in Eq. (3.56) which favors small frequencies. For intermediate N 0 , the eective plasmon frequency is comparably large, whilea pl q is not big enough to compensate for the reduced plasmonic density of states, yielding a reduced eective coupling pl here. After discussing the plasmonic characteristics, we now turn to their in uence on super- conducting properties. Since the electron-plasmon coupling matrix element a pl q is rather isotropic, we neglect the superconducting gap anisotropy over the Fermi surface and assume s-wave symmetry of the latter [131]. Under this assumption, we can solve the linearized SC-DFT gap equation in energy space [52,102] given in Eq. 3.61. We numerically solve the SC-DFT gap equation, including fully energy dependent ex- pressions forZ() andK(; 0 ) (as given in Section 3.2.2) subject to the in uence of dierent doping levels and varying dielectric environments. In Fig. 4.12 we show resultingT c for electron-doping levels increasing fromn& 0:5 to 1:0 88 CHAPTER 4. SQUARE LATTICE 4.3. SUPERCONDUCTIVITY 0 2 4 6 8 λ plasmon phonon total 0.0 0.1 0.2 0.3 ω log 0.0 0.5 1.0 1/ε env 0.0 0.5 1.0 1.5 2.0 μ ∗ 1.0 1.5 N 0 Figure 4.11: Plots showing the various coupling strengths ( pl , ph and ) and averaged bosonic frequency (! Log ) for dierent screening (left) and doping (right) scenarios. For the left panels we used n 0:8 (N 0 0:9 eV 1 ) and " env = 1 for the right panels. Dashed (solid) gray lines correspond to the phononic (total) parameters using a constant frequency ! ph = 45 meV and an electron-phonon coupling matrix element g 2 = 50 meV 2 . (corresponding toN 0 1:5 eV 1 to 0:7 eV 1 ) and for varying dielectric environments. The phononic properties are set to! ph = 40 meV (similar to the optical modes in TMDC's [132]) and g 2 = 45 meV 2 yielding a realistic ph 1:5 3:5. Both, the doping level and the dielectric environment, aect the critical temperature. We nd enhanced T c for large N 0 and " env as well as for small N 0 and " env . These regimes are labeled as \phononic" and \plasmonic", respectively. In the phononic regime, pl is relatively small (see inset of Fig. 4.12), while ph = 2N 0 g 2 =! ph constantly increases with N 0 . The increasing trend of the critical temperature thus follows the eective phononic coupling in this regime. In the plasmonic regime, pl is large and increases with decreasing environmental screening. Since " env has no eect on the phononic properties in our model, the increasing T c follows pl . For large N 0 and small " env , and thus rather large electron-plasmon and electron-phonon couplings, we also nd enhanced critical temperatures. Although the total coupling seems to be the strongest here, we do not nd the highest critical temperatures in this regime. In order to understand this interplay between the two coupling channels in more detail, we now x the doping level to n 0:8 (N 0 0:9 eV 1 ) and study the eects of simultane- ously varying electron-phonon coupling g 2 and environmental screening " env . The former now solely controls ph and the latter pl . In Fig. 4.13 we show results for ! ph = 60 meV and g 2 22 76 meV 2 , yielding ph 0:65 2:2. Besides T c , we also show the total 89 4.3. SUPERCONDUCTIVITY CHAPTER 4. SQUARE LATTICE Figure 4.12: T c as a function of doping (rendered by N 0 ) and environmental screening, including the eects of phonons, plasmons, and static Coulomb interaction. Solid black lines correspond to constant values ofT c . \Phononic" and \plasmonic" regimes are marked by blue and red dotted lines, respectively. The inset depicts the eective electron-plasmon coupling pl . 90 CHAPTER 4. SQUARE LATTICE 4.3. SUPERCONDUCTIVITY Figure 4.13: T c as a function of phononic coupling and environmental screening. The upper panel shows the full critical temperature including the eects of phonons, plasmons and static Coulomb interaction. The lower panels show the total eective coupling tot , total eective frequency ! tot log , and total eective Coulomb repulsion tot . 91 4.3. SUPERCONDUCTIVITY CHAPTER 4. SQUARE LATTICE eective coupling tot , frequency ! tot and static repulsion tot dened by tot = ph + pl (4.12) ! tot = exp log(! ph ) ph + log(! pl ) pl ph + pl (4.13) tot = tot 1 + ln(E F =! tot ) with tot = + pl : (4.14) It is important to note that pl needs to be taken into account in both, tot and tot . This is a direct result of the two terms of the dynamic Coulomb contribution dened by Eq. (3.55). The rst (constant) term contributes to tot , which reduces T c , and the second (dynamic) term to tot , which enhancesT c . In Fig. 4.11 we show how these total eective parameters depend on the environmental screening and the doping level. Most importantly, we see here that the total eective coupling is simply enhanced, while the total eective frequency and static repulsion are strongly reduced compared to the pure plasmonic quantities. Based on the interplay of these three parameters, we can now qualitatively understand the behavior of T c shown in Fig. 4.13. In the phononic regime neither the total eective frequency nor the static repulsion change drastically with increasing ph (see lower panels in Fig. 4.13). Only the total eective coupling increases with ph , which is responsible for the increasingT c trend, here. Similarly, in the plasmonic regimeT c increases towards small environmental screenings due to the enhancement of tot and ! tot log , which is driven by the increasing trend of the plasmonic pl and ! pl log , respectively. Here, however, tot and tot also increase with pl , which reduces the increasing trend inT c . In the remaining regime the total eective coupling tot is enhanced by both coupling channels. Interestingly, T c seems to be mostly controlled by ph here (lines of constantT c are vertical), whereas" env and thus pl seem to have a negligible eect. However, by studying the eective total parameters, we realize that T c is simultaneously enhanced and reduced by the counteracting trends in tot and tot with increasing pl . It is thus an interplay between both coupling channels which is responsible for T c here. Finally, we draw attention to the ratio map in Fig. 4.14, which depicts R =T c =T ph+ c , whereT ph+ c is the critical temperature including only phononic and static Coulomb eects (no plasmonic contributions). This quantity reveals those regimes where the plasmons in- crease T c (R> 1) in comparison to the situation without any plasmonic in uence (R = 1) and those regimes where the plasmons reduce it (R< 1). We nd a narrow regime of plas- monic enhancement for most dielectric screenings and rather small phononic couplings. The amount of the plasmonic enhancement within this stripe is controlled by the environmental screening. In the freestanding situation (" env = 1), we nd full critical temperatures which are enhanced by factors up to 5. By increasing the environmental screening, we decrease this enhancement. If we increase the phononic coupling we also nd a decreasing plasmonic enhancement, which is in line with previous ab initio data by Akashi and Arita [32]. If we further increase 92 CHAPTER 4. SQUARE LATTICE 4.4. CONCLUSION Figure 4.14: Plasmonic enhancement-reduction map. We show the ratio of the full transi- tion temperatureT c and the phononic oneT ph+ c (excluding plasmonic eects). Red parts represent the regime of plasmonic enhancement and blue the plasmonic reduction. Notice the similarity between the full numerical ratios here and the analytically expected ones from Fig. 3.3 the phononic coupling we arrive in a situation with plasmonic reduction, which holds for all environmental screenings. Similar charge- uctuation induced reductions of the critical temperature has also been reported in Ref. [133]. 4.4 Conclusion We studied the interplay of phonons and plasmons in the creation of a superconducting state in a model system. We applied the drastic simplication of assuming vanishing electron- plasmon coupling below a critical momentum value which allowed us to treat phonons and plasmons with a similar dictionary. We have shown how the plasmon frequencies and the electron-plasmon coupling of a layered metal can be tuned by electron doping and envi- ronmental screening. While the latter always reduces both, the plasmon frequencies and the couplings, the former has more subtle consequences. Most importantly, we found that the coupling is enhanced while the frequency is reduced for lowered density of states at the Fermi level. These strongly screening and doping dependent plasmonic properties imprint strong changes to the full superconducting transition temperature, allowing for external 93 4.4. CONCLUSION CHAPTER 4. SQUARE LATTICE control of the latter. Thereby, we demonstrated that the electron-plasmon coupling does not only contribute to the total attractive coupling, but also to the total repulsive static Coulomb term. As a consequence, there is both plasmonic enhancement and reduction depending on the doping level and dielectric environment. In order to nd a sweet spot for plasmon-enhanced superconductivity in an experiment, the major task is to enhance the eective electron-plasmon coupling pl (; 0 ) while trying to keep the eective plasmon frequency ! pl low. The latter is important since pl neces- sarily also adds to the eective Coulomb repulsion , which is, however, decreased by a decreasing ! pl . These properties are indeed not contradicting as evident from Eq. 3.56: an enhanced plasmonic spectral weight around small frequencies enhances the eective coupling. Based on these considerations we can dene two rules: (I) Use two-dimensional metals with gapless plasmonic modes. (II) Use metals with a reduced density of states at the Fermi level, i.e. with a small eective mass. Both rules guarantee a decreased! pl and an enhanced pl . In order to further enhance the eective plasmonic coupling we can dene two additional rules: (III) Reduce the environmental screening. (IV) Avoid inter-band polarization eects to reduce Landau damping, i.e. try to nd a single free-standing metallic band embedded in an electronic band gap. The former actually increases ! pl which, however, can be compensated by a strongly in- creased pl at small ph (see top left corner of the ratio map in Fig. 4.14). The latter is discussed in Ref. [134]. This also leads to the last rule: (V) Use a material with small intrinsic electron-phonon coupling. Promising candidates which fulll at least some of these guidelines include slightly electron or hole doped semiconducting TMDC monolayers [135, 136] (metallic TMDC's might show too strong electron-phonon interactions [137], too strong inter-band polar- izations [138], and too high density of states at the Fermi level), monolayers of recently proposed 1T-AlCl 2 [134], and in general single-layersp-electron systems, such as hexagonal boron nitride or functionalized graphene [139,140]. In order to disentangle plasmonic and phononic eects from each other, it would be best if the eects of both, the environmental screening and the doping level, could be experi- mentally studied. If there are strong changes in T c by changing the dielectric environment or if T c increases by decreasing the density of states at the Fermi level, our results show that it is very likely that T c is signicantly controlled by the coupling between electrons and plasmons. 94 CHAPTER 4. SQUARE LATTICE 4.4. CONCLUSION These ndings clearly show that static and dynamic Coulomb interaction eects need to be accurately considered in order to explain superconducting properties from a theoretical point of view. At the same time they point towards exciting new directions in the eld of on-demand material-property design using layered systems. Here, a sophisticated choice of materials can increase critical temperatures by combining advantageous properties from dierent materials. 95 Chapter 5 Material Realistic Calculations for MoS 2 In this chapter we apply the theory discussed and developed so far, to material realistic cal- culations for MoS 2 . We have shown in chapter 4 that plasmonic excitation can be strongly tuned in a layered system and we have also shown that the electron-plasmon coupling can have a signicant impact on the superconducting properties of the system. An analogous but even richer phenomenology can be expected in the structurally related monolayer tran- sition metal dichalcogenides (TMDC's) MX 2 , where M stands for a transition metal and X for a chalcogen atom. These materials host rich plasmonic physics including an interplay of plasmons with charge density waves [141{143] and the rst plasmon based applications have already been pro- posed [144{146]. In the rst section of this chapter we will present the results of ab initio calculations of the electronic and phononic properties of monolayer molybdenum disulphide (MoS 2 ). In Section 5.4 we discuss some of the plasmonic properties of MoS 2 coming from our publication Ref. [126]. The overall interest in plasmonics has been on the rise with ex- perts such as Feliciano Giustino calling for "a systematic investigation of electron-plasmon couplings in a wide class of materials" [100]. The transition metal dichalcogenides cer- tainly constitute an interesting class of materials which is why we have focused on setting up realistic ab initio calculations of the screened Coulomb interaction in these materials, along with the realistic treatment of the static Coulomb interaction both with [119] and without [147] substrates. In Section 5.5 we consider the ab initio treatment of the static Coulomb parameter, , that enters superconductivity calculations, as discussed in our publication Ref. [49]. We compare the relative contributions to coming from inter- and intra-valley scattering to answer whether a static Coulomb driven unconventional pairing mechanism can be present in the superconducting phase of MoS 2 . 96 CHAPTER 5. MOS 2 CALCULATIONS MoS 2 is a very well studied material and a single layer eld eect transistor has even been fabricated using the material [39]. Two stable congurations are known for the MoS 2 crystal lattice although only the 2H phase [148] can be used to produce monolayers of the material. Seeing as we are interested in pseudo 2d materials in this thesis, we will only focus on the 2H phase. The crystal structure is shown in Fig. 5.1. It is hexagonally structured and a single layer of the crystal consists of a S-Mo-S trilayer. In bulk MoS 2 the sulfur layers only interact weakly through Van der Waals forces thereby allowing the exfoliation of single layers of MoS 2 [6]. a Figure 5.1: Left: Ball and stick representation of the 2H phase of MoS 2 , showing the plane of molybdenum atoms (black) sandwiched between layers of sulfur atoms (yellow). Right: Top view of the crystal structure, notice the hexagonal structure. Also shown is the FBZ of the reciprocal lattice in blue outline, as well as the lattice parameter, a. In Fig. 5.1 we also show the FBZ of the reciprocal lattice, consisting of a diamond containing one Mo atom and two S atoms. The lattice parameter has been experimentally determined to be a = 3:125 A while the interlayer distance is given by c = 12:066 A [149]. One interesting, open question concerning MoS 2 (and other TMDC's) is the fact that experiments have found the superconducting transition temperature to be lower in the monolayer ( 2K) than in the bulk ( 10K) [150]. This trend is opposite to what is seen in materials such as FeSe on SrTIO 3 , in which the highest T c value is found for the monolayer [151,152]. One possible reason for this behavior could be due to the reduction of Coulomb screening for thinner samples or due to the dierent screening scenarios coming from substrates. This is one question we attempt to answer in this work. Electron doping in these compounds is conventionally achieved by intercalation of dopant atoms, which however changes the band structure of the host compound and intro- duces signicant impurity scattering. Recently, an alternative experimental approach has been established which relies on the application of a strong electric eld, whose magnitude can be controlled and greatly magnied by combining solid gates with liquid ions. This way, eective 2d electron carrier densities of up to 10 13 cm 2 with high mobility can be 97 5.1. AB INITIO CALCULATIONS CHAPTER 5. MOS 2 CALCULATIONS achieved [153]. In Section 5.4, we model such eld-eect doping by an adjustment of the chemical potential. In reality, since the DFT approach specically revolves around the ground state electron density, the calculation results naturally changes if the electron car- rier density changes. This has been found to strongly aect the phonon properties of the system [154] and therefore also the superconducting transition temperature. For this rea- son it is important to recalculate the electronic- and phononic dispersions for each doping level considered, using a relaxed lattice parameter for the relevant electron density. 5.1 Ab initio Calculations We use DFT and DFPT calculations as implemented in Quantum ESPRESSO [155,156] in order to obtain the ab initio electronic and phononic properties of MoS 2 . The parameters used in the calculations can be found in Appendix B.1. 5.1.1 DFT - electronic dispersion The electronic dispersion (band-structure) for the high symmetry path ! M! K! ! through the FBZ is shown in Fig. 5.2. Notice the existence of a direct band gap of 1:87 eV at the K-point (compared to the experimental optical gap of 1:80 eV [46]), which makes monolayer MoS 2 a promising candidate for the active material in low-power electronics [157]. In the bulk form MoS 2 has an indirect band gap [158] and is therefore less ideal for electronic components. For completeness we also show the electronic DOS in Fig. 5.2. The very low energy bands (around15eV) come mainly from the 3s orbitals of sulfur. The bands grouped just below the band-gap arise due to the hybridization of 4d orbitals coming from molybdenum and 3p orbitals coming from sulfur, while the lowest conduction bands are mostly made up of the remaining Mo d orbitals' characters. This separation of energy bands coming from dierent classes of atomic orbitals is what allows us to construct a minimal model using Wannier functions as we discuss later. 5.1.2 DFPT - phonon dispersion The phonon dispersion and phonon DOS is shown in Fig. 5.3 along the same ! M! K! ! path as before. The phonon frequencies are nowhere imaginary which indicates a stable crystal structure. There are 3 acoustic phonon branches stemming from in plane and out of plane modes. The in plane modes have linear dispersions around while the out of plane mode goes as q 2 . the optical gap is smaller than the band gap due to the excitonic binding energy 98 CHAPTER 5. MOS 2 CALCULATIONS 5.1. AB INITIO CALCULATIONS Γ M K Σ Γ −15 −10 −5 0 5 E−E F (eV) DOS(E) Figure 5.2: Electronic dispersion of monolayer MoS 2 from DFT calculations along high symmetry points through the FBZ. The calculation was performed for 10% electron doping. The bands in blue, surrounding the Fermi-level (set such that E F = 0), are the ones used in our minimal model discussed in Section 5.2. Notice the direct band gap at the K-point. 99 5.2. MODEL - WANNIER FUNCTIONS CHAPTER 5. MOS 2 CALCULATIONS Γ M K Σ Γ 0 10 20 30 40 50 60 Ω (meV) DOS ph Figure 5.3: Phononic dispersion of monolayer MoS 2 from DFPT calculations along high symmetry points through the FBZ. The calculation was performed for 10% electron doping. 5.2 Model - Wannier Functions The true strength of our calculations lie in the fact that we can handle extremely dense k-grids (consisting of up to 10 5 k-points) which is made possible by the use of eective models to describe real materials. This is done by using Wannier constructions that extract the tight-binding hopping parameters coming from the projection of the full electronic system onto only a couple of atomic orbitals. Specically, we use the Mo d z 2, d xy and d x 2 y 2 orbitals to construct a tight-binding model that very accurately describes the highest valence band and two lowest conduction bands of MoS 2 , as is shown in Fig. 5.4. We also show the DOS in the side panel with the DOS coming from the full ab initio calculations overlayed to show that our model fully captures the low-energy properties of the electronic system. The use of a minimal basis is computationally extremely valuable as most many- body calculations scales with the number of eigenstates in the system. In the case of MoS 2 this allows us to decrease the number of involved states from 23 to only 3. We see in Fig. 5.4 a peak in the highest valence band at the K-point. Thus if the material is hole doped it would open a Fermi surface around the K-points in the Brillouin 100 CHAPTER 5. MOS 2 CALCULATIONS 5.3. SUBSTRATE SCREENING Γ M K Σ Γ −3 −2 −1 0 1 2 E−E F (eV) DOS(E) Figure 5.4: Electronic dispersion of 10% electron doped MoS 2 from the minimal model obtained with Wannier constructions, along high symmetry points through the FBZ. The Wannier construction is done using the Mo d z 2, d xy and d x 2 y 2 orbitals. The density of states in blue comes from the 3 bands shown while the grey line is the full density of states (as shown in Fig. 5.2). zone. Similarly, with low electron doping the valley of the conduction band at the K-point would become occupied thereby also creating a circular Fermi surface at the K-point. If more electrons are added to the system the conduction band minimum at the -point would also become occupied. A graphical depiction of these possible Fermi surfaces are shown in Fig. 5.5. Notice, however, that the orbital characters of the dierent valleys and peaks vary and therefore the response of the system to spin-orbit coupling diers under dierent doping scenarios. We discuss this point in greater detail in Section 5.4. 5.3 Substrate Screening The available screening channels can be divided into internal and external channels, which both contribute to the strength of the Coulomb interactions. Here, internal processes refer 101 5.3. SUBSTRATE SCREENING CHAPTER 5. MOS 2 CALCULATIONS Figure 5.5: Sketch of the Fermi surfaces in hole (left) and electron doped (right) monolayer MoS 2 . The dierent orbital characters are indicated by red (d z 2) and blue (d xy andd x 2 y 2) lled surfaces. Points of high symmetry are indicated by dierent markers. to the screening due to transitions between electronic states within the MoS 2 layer. The external screening arises due to the polarizability of adjacent substrates or adsorbates with dielectric constant " sub . Parts of the calculations on and parametrization of the Coulomb interaction were pre- viously described in Ref. [147]. Here, we follow a similar procedure and make use of the Wannier function continuum electrostatics approach (WFCE) [119] to include the screening eects of substrates, as described in the following. The bare interaction matrix U (q) in the orbital basis , 2fd z 2; d xy ; d x 2 y 2g of the Mo atoms is obtained for the freestanding undoped material via RPA calculations using the Spex and Fleur software codes [159, 160]. To parameterize the Coulomb interaction, we use the sorted eigenbasis of the bare interaction to diagonalize it, U diag (q) = 0 B @ U diag 1 (q) 0 0 0 U diag 2 0 0 0 U diag 3 1 C A: (5.1) Here, the diagonal matrix elements are given by U diag i =he i jUje i i (5.2) using the eigenvectors of ^ U(q! 0) e 1 = 0 @ 1= p 3 1= p 3 1= p 3 1 A ;e 2 = 0 @ p 2=3 1= p 6 1= p 6 1 A ;e 3 = 0 @ 0 1= p 2 1= p 2 1 A : (5.3) In the sorted systemU diag 1 (q) is the leading eigenvalue of the bare interaction and the other two eigenvalues are approximately constant. For the leading eigenvalue, we obtain a t of the form U diag 1 (q) = 3e 2 2" 0 A 1 q(1 + q) (5.4) 102 CHAPTER 5. MOS 2 CALCULATIONS 5.3. SUBSTRATE SCREENING Table 5.1: Parameters describing the Coulomb interaction in MoS 2 . Parameter Value ( A) 2.091 U diag 2 (eV) 0.810 U diag 3 (eV) 0.367 " 1 9.253 d ( A) 9.136 " 2 3.077 " 3 2.509 with the area of the 2d hexagonal unit cellA = p 3 2 a 2 and the lattice parametera = 3.18 A. The standard least-squares method is used to determine the parameter . The factor 3 in Eq. (5.4) arises from the fact that we use three orbitals to describe the system and treat the Coulomb interaction in the eigenbasis of the bare interaction. describes how the eective height aects short wavelengths, which means that it is a structure factor and becomes important at large wavevectors q. The value of is given in Tab. 5.1. The screened matrix elements in the eigenbasis of the bare interaction are then obtained for the undoped system via V diag i (q) = h " diag i (q) i 1 U diag i (q) (5.5) where " diag i (q) accounts for the material specic interband polarizability and the polariz- ability of the substrate. Its diagonal representation is given by " diag (q) = 0 @ " 1 (q) 0 0 0 " 2 0 0 0 " 3 1 A (5.6) where the constants " 2 and " 3 (see Tab. 5.1) describe microscopic screening eects which are similar to the bulk. The macroscopic eects are described by the leading eigenvalue via " 1 (q) =" 1 1 1 2 e 2qd 1 + ( 1 + 2 )e qd + 1 2 e 2qd (5.7) with i = " 1 " sub;i " 1 +" sub;i : (5.8) The involved parameters are derived from ts to the ab initio calculations (see Tab. 5.1) for the freestanding layer. The surrounding substrates have dielectric constants " sub;1 above 103 5.3. SUBSTRATE SCREENING CHAPTER 5. MOS 2 CALCULATIONS 0 0.2 0.4 0.6 0.8 1 1 5 10 15 q (1/ ˚ A) Macroscopic Screening ε (q) ε ub =∞ ε b = 1 (b) Figure 5.6: Macroscopic element ( 1 ) of the dielectric matrix in eigenbasis of the Coulomb interaction. From bottom to top, the substrate dielectric constants used for " 1 (q) are " sub = 1, 5, 10, 50,1. and" sub;2 below the monolayer which can be varied using Eq. (5.8). In the case of vacuum surrounding the monolayer (" sub;1 =" sub;2 = 1) Eq. (5.7) simplies to " 1 (q) =" 1 " 1 + 1 (" 1 1)e qd " 1 + 1 + (" 1 1)e qd : (5.9) The macroscopic screening for various substrates are shown in Fig. 5.6. In the long- wavelength limit," 1 (q) is fully determined by the dielectric background: " 1 (q! 0) =" sub . In the opposite limit (q & 1 A 1 ), screening in the undoped case is solely due to the microscopic inter-band polarizability of the MoS 2 layer itself, " 1 (q) 9:3 " 1 , but unaected by the dielectric environment [119,127]. Once we have obtained the diagonal dielectric matrix " diag (q) we can calculate the screened Coulomb interaction in the eigenbasis using Eq. (5.5) together with Eqs. (5.1) and (5.4) and the parameters in Tab. 5.1. Afterwards, we can transform to the orbital basis using the eigenvectors in Eq. (5.3). This analytic description allows to evaluate the bare and screened Coulomb matrix elements at arbitrary momenta q in the rst Brillouin zone and for arbitrary dielectric environments. To get the real space values we need to do a simple Fourier transform, resulting in the on-site bare and screened Coulomb matrix elements given in Tab. 5.2. In addition to the external and inter-band screening, we also have to include metallic intra-band screening by the conduction electrons in the case of electron doping. We then arrive at the fully screened The function "1(q) describes macroscopic screening eects while the functions "2;3, describe only mi- croscopic screening 104 CHAPTER 5. MOS 2 CALCULATIONS 5.4. PLASMONIC PROPERTIES Coulomb interaction ^ W (q;!) = ^ V (q) 1 ^ V (q) ^ (q;!) 1 , where ^ (q;!) is the intra-band polarizability and where ^ W (q;!), ^ V (q) and ^ (q;!) are matrices in the Wannier function basis. The polarizability is obtained using RPA for the lowest conduction band [126]. Table 5.2: Bare on-site U as well as background screened on-site V and fully screened on-site Coulomb matrix elements W for the three important orbitals in real space. Values for W are in the range of low electron doping x 0:04 (x being the number of additional electrons per unit cell) in the fth (W low , K is occupied) and for high electron doping x 0:13 in the last column (W high , K and are occupied). bare undoped doped orbitals U (eV) V (eV) W low (eV) W high (eV) d z 2 d z 2 9.11 1.55 0.82 0.68 d z 2 d xy 8.30 1.29 0.58 0.44 d z 2 d x 2 y 2 8.30 1.29 0.58 0.44 d xy d xy 8.89 1.49 0.80 0.64 d xy d x 2 y 2 8.52 1.35 0.65 0.51 d x 2 y 2 d x 2 y 2 8.89 1.49 0.80 0.64 5.4 Plasmonic Properties There is substantial spin-orbit coupling (SOC) in TMDC's [161], with a primary eect on the low-energy physics by introducing a splitting of the valleys in the lowest conduction band and of theK valleys in the highest valence band. Although all of these Fermi surface characteristics can be experimentally sampled by means of eld eect electron or hole doping [153], the resulting impact on the plasmonic dispersions have not been studied before. Here, we present an extensive study of the plasmon dispersion at arbitrary momenta along paths throughout the whole Brillouin zone for dierent doping levels. Specically we are interested in inter-valley plasmons which have not been studied in TMDC's so far. In order to highlight the multi orbital character of the Fermi surface (see Fig. 5.5) and the presence of spin-orbit coupling we consider hole and electron doped cases. In the hole doped example we show how spin-orbit coupling aects the inter-valley plasmons while the electron doped case is used to study the in uence of the multi pocket structure of the Fermi surface. Thereby we gain a comprehensive and realistic picture of the most impor- tant contributions to the low energy plasmon modes in monolayer TMDC's. The collective plasmon modes are described by the polarization and dielectric functions, which we evaluate in several steps, starting with a G 0 W 0 calculation to determine the 105 5.4. PLASMONIC PROPERTIES CHAPTER 5. MOS 2 CALCULATIONS electronic band structure for the undoped system. We then obtain an eective 3-band model by projecting to a Wannier basis spanned by the Mo d x 2 y 2, d xy and d z 2 orbitals, which has been found to accurately describe the highest valence band and the two lowest conduction bands with tight-binding hopping matrix elements t , where and are the orbital indices (as was discussed in Section 5.2). The same projection is used to obtain the static part of the Coulomb interaction in the Wannier basis, which is screened by all bands including those, which are not included in the minimal 3-band-model . This procedure leads to an eective material-specic model with screened Coulomb U and hoppingt matrix elements in the orbital basis, describing the undoped system in its ground state. We nd that this treatment is essential to derive material realistic plasmonic dispersions upon doping, see Appendix C.3 for a comparison to simple kp models. In contrast to simpliedkp models [164,165], which utilize bare Coulomb matrix elements at this stage, our interaction matrix elements are strongly reduced due to screening eects from the electronic bands which are neglected in the kp models. As a result of the 2d layer geometry, these dielectric properties cannot be modeled by a simple dielectric constant but have to be described as a q-dependent dielectric function [147,166]. In order to obtain the dynamic response in the doped system, we determine the dynamic susceptibility within the 3-orbital basis by evaluating the polarization in the random phase approximation (RPA), which is given for a single spin channel by (q;!) = X 1 2 k M 1 2 f 2 (k + q)f 1 (k) ! +i +E 2 (k + q)E 1 (k) ; (5.10) where q and k are wave vectors from the rst Brillouin zone, i band indices, f i (k) Fermi functions for the energies E i (k) andi a small broadening parameter. The overlap matrix elements are given by M 1 2 = c 1 (k)c 1 (k) c 2 (k + q)c 2 (k + q), where c i (k) is the expansion coecient of the eigenfunction corresponding to E i (k) in the orbital basis. Here, we already reduced the polarization tensor of 4 th order to a matrix to describe density-density correlations only. Hence, we neglect orbital exchange (Fock-like) matrix elements as well as elements with three or even four dierent orbital contributions. A detailed analysis of the full background screened Coulomb tensor U shows, that these elements are in general one order of magnitude smaller or even vanish due to symmetries, which convinces us to stay with density-density like elements. Using the full density-density polarization (q;!) = " (q;!) + # (q;!) the dielectric function is obtained via the following matrix equation "(q;!) = 1U(q)(q;!); (5.11) In all ab initio calculations we used an inter-layer separation of 35 A. The G0W0 calculations are performed with the Vienna ab initio simulation package (VASP) [162, 163], while the Coulomb matrix elements are obtained from the SPEX code [159] with FLAPW input from the FLEUR code [160] as described in [147]. For the involved Wannier projections we use the Wannier90 package [62]. 106 CHAPTER 5. MOS 2 CALCULATIONS 5.4. PLASMONIC PROPERTIES where the background screened Coulomb interaction enters via U(q). By including an eective spin-orbit coupling [167] the spin degeneracy is removed but time reversal sym- metry is preserved. Then, the spin resolved band structure still obeys E " (k) = E # (k) and the total polarization including the spin summation can be written as (q;!) = " (q;!) + " (q;!). The dielectric function describes the screened Coulomb matrixV (q;!) =" 1 (q;!)U(q) and implicitly denes the plasmonic dispersions by " m (q;!) = 0, where " m is the macro- scopic part of the dielectric function [119]. The most promising experimental method to map these plasmon modes is electron energy loss spectroscopy (EELS) from Eq. 2.78, measuring the imaginary part of the inverse dielectric function EELS(q;!) = Im 1 " m (q;!) ; (5.12) which is sensitive to both collective and single-particle excitations (visible as maxima in the EELS spectra) [168]. The combination of our material realistic description of the undoped system and the very accurate band structure for the RPA evaluation yields indeed quite accurate plasmon dispersions compared to full ab initio results, as we show for NbS 2 in Appendix C.3. Hole Doped We x the chemical potential such that there are holes in the valence band in the K and K 0 valleys only. The resulting Fermi surfaces consists of circle-like areas around the K points (see Fig. 5.5), which have mainly d xy =d x 2 y 2 character and depend on spin-orbit coupling. Hence, we expect low energy plasmon modes for q (intra-valley) and q K (inter-valley), which are possibly in uenced by SOC. In Fig. 5.7 (a) we show an intensity plot of the real part of the polarization function for scattering withind xy orbitals along the complete path ! K without SOC . Next to some band-like structures (red) we clearly see the particle-hole continuum (blue). In comparison to the corresponding EELS data in Fig. 5.7 (c), we see that for higher momentum transfers (away from ) the EELS maxima closely follow the band-like characteristics of the polar- ization function. For small momenta around we nd a clearly separated band in the EELS spectra, which can not be seen in the real part of the polarization. This separated band arises from the well known p q-dispersive intra-valley plasmon mode for smalljqj in 2d [164], while we nd a linear-dispersive mode around K stemming from an inter-valley plasmon [29]. These activation laws are consistent with the generalized expression for the See Appendix B.4 for calculation details. 107 5.4. PLASMONIC PROPERTIES CHAPTER 5. MOS 2 CALCULATIONS Figure 5.7: Real and imaginary parts of the polarization functions (d xy =d xy channel) and EELS spectra for hole doped MoS 2 without (top row) and with (bottom row) spin orbit coupling. The insets in (a) and (d) illustrate the Fermi surface pockets around K andK 0 . plasmon dispersion relation dened by the dielectric function via [31], !(q) = ~v F q s 1 + [N 0 U(q)] 2 (1=4) +N 0 U(q) ; (5.13) where v F is the Fermi velocity, N 0 the density of states at the Fermi level and U(q) the macroscopic background screened Coulomb interaction of the undoped system. In the long-wavelength limit (q ! 0), the Coulomb potential remains unscreened, i.e. in leading order U/ 1=q, resulting in a square-root renormalization of the otherwise linear dispersion. However, in the opposite short-range limit, i.e. at the zone boundary, the screened Coulomb potential approaches a constant, and therefore the resulting dispersion of the dielectric function is linear in q, same as the polarization function itself. Thus for momenta away from it is sucient to study the polarization function to understand how the resulting plasmon dispersion will behave. Of special interest are damping eects, which are known to attenuate plasmon modes which merge with the particle-hole continuum. Here the square-root mode around be- haves in a distinctly dierent manner compared to the linear modes originating at K. At suciently small momentum transfers q < q c the square-root modes are more sepa- rated from the nearby particle-hole continua [Fig. 5.7 (c) and (f)], and therefore better protected from decomposition via hybridization and Landau damping [expressed as non- 108 CHAPTER 5. MOS 2 CALCULATIONS 5.4. PLASMONIC PROPERTIES vanishing imaginary parts of the polarization as shown in Fig. 5.7 (b) and (e)] compared to the linear modes originating at nite momenta. In contrast, the linear plasmon modes are much closer to their neighboring continua [Fig. 5.7 (c)], which leads to attenuation eects, re ected in reduced oscillator strength and broadening of the peaks. There is a signicant dierence in the oscillator strengths of these modes, which can be several orders of magnitude apart as can be seen in Fig. 5.7 (c) and (f). Hence, in order to clearly detect these linear plasmon modes in experiments, it may prove practical to use a logarithmic scale to shield the dominant square-root mode around q = , as shown in Fig. 5.7 (c) and (f). When we account for spin-orbit coupling the relative depth of the K and K 0 pockets shifts. In this case momentum transfer of q = K no longer connects points on the Fermi surface belonging to dierent hole pockets, which results in two clearly visible characteris- tics in the polarization of Fig. 5.7 (d): (1) At q = K the scattering process is possible only for a nite energy dierence, which opens a nite energy gap of 250 meV. (2) The Fermi surfaces at K and K 0 are now of dierent sizes but can still be connected with slightly smaller and larger q, resulting in gap-less linear modes originating slightly shifted from K as seen in Fig. 5.7 (d). We conclude that the plasmonic features in hole doped MoS 2 are qualitatively similar to graphene as long as SOC is not taken into account and the K valley is occupied solely. Upon inclusion of SOC the linear plasmon mode around K is shifted leading to a gapped excitation spectra at this point. Electron Doped The lowest conduction band is characterized by two prominent minima around K and . Without SOC these minima are separated by only 90 meV. Hence, in contrast to the hole doped case, small variations in the electron doping can change the Fermi surface drastically. In order to study these changes, we will neglect the SOC for the beginning and choose two doping levels, resulting in Fermi surfaces comparable to the hole doped case (i.e. K valley occupation only) and a surface with additional pockets at , labeled by low- and high-doping respectively. Since the K valley is described by d z 2 orbitals and the valley predominately by d xy and d x 2 y 2 states, we focus on corresponding diagonal orbital channels in in the following. O-diagonal elements betweend z 2 andd xy =d x 2 y 2 orbitals are negligible here (o-diagonal terms betweend xy =d x 2 y 2 states are similar). The corresponding polarization functions are shown along the path KM through the whole Brillouin zone in Fig. 5.8. Analogous to the hole doped case, we observe around q = the expected resonances arising from intra-valley scattering. This is naturally present in both high and low electron doping cases. By inspection of Fig. 5.5 we can understand the structure of the polarization for larger q. The momenta q = ; K and M connect dierent valleys. Therefore, in the high 109 5.4. PLASMONIC PROPERTIES CHAPTER 5. MOS 2 CALCULATIONS Figure 5.8: (Color online) (a) The polarization function for d z 2=d z 2 scattering at low elec- tron doping concentration (only K valleys are partially occupied) without SOC. (b) The polarization function for d xy =d xy scattering at elevated electron doping concentration (K and valleys are partially occupied) without SOC. In (c) the same situation as in (b) is shown, but with the eect of spin-orbit interaction. The insets show illustrations of the Fermi surfaces. electron doping case (with partially occupied) we expect additional inter-valley plasmon branches close to these momenta. At q K we observe plasmon bands in both, high and low doping cases since this momentum transfer allows inter-valley scattering between K and pockets. As we are calculating orbital resolved polarization functions, the observed low energy excitation in Fig. 5.8 (a) is due to K$ K 0 and thus d z 2 scattering, whereas in Fig. 5.8 (b) it is due to $ and correspondingly d xy =d x 2 y 2 scattering (K$K 0 scattering is obviously still present, but can only be seen in the d z 2 polarization). Finally, momentum transfers q = M and can connect dierent valleys, and there- fore we nd a gap-less linear inter-valley plasmon mode originating at this point only in the high doping case. In the low doping case, we observe a gapped ( 0:1eV ) excitation at q = M, originating from a K$ excitations. While the SOC has a negligible eect on the d z 2 valley at K it splits the d xy =d x 2 y 2 valleys at resulting in minima at comparable energies. The corresponding Fermi surface for a single spin component is indicated in the inset of Fig. 5.8 (c). The six points decompose into two distinct sets, and 0 . Fermi pockets within each of these subsets are mutually connected by 2=3 rotations and remain equivalent after inclusion of SOC, while the degeneracy of and 0 is lifted by SOC. As a consequence, the phase space for $ 0 is lost and the gap-less excitations at q and q K must vanish, but $ scattering processes are still possible. Consequently, we see in the corresponding polarization for the d xy channel with SOC in Fig. 5.8 (c) gap-less modes only at and M (since M connects $ valleys). Since the Fermi surface around K is not changed drastically upon SOC, the corresponding polarization for the d z 2 channel is very similar to the one obtained without SOC. 110 CHAPTER 5. MOS 2 CALCULATIONS 5.5. SUPERCONDUCTIVITY We can also explore how environmental screening aects the MoS 2 plasmon dispersion (similarly to the study in Sec. 4.2). The EELS spectra for low-momentum plasmons under a range of doping and environmental screening situations are shown in Fig. 5.9. We see here the same suppression of the plasmon dispersion, with increased environmental screening, seen earlier for the model system. This eect appears to be present as a general trend of low-dimensional materials. Figure 5.9: Low momentum plasmon dispersion of MoS 2 under various doping and envi- ronmental screening scenarios. Spin orbit coupling is also included. 5.5 Superconductivity In this section, we apply our quantitative theory to study how Coulomb interactions af- fect the superconducting transition temperature in MoS 2 as a representative example of TMDC's. We focus on the aect substrate screening and doping has on the eective static Coulomb interaction, , that enters most superconductivity theories. 111 5.5. SUPERCONDUCTIVITY CHAPTER 5. MOS 2 CALCULATIONS Realistic static Coulomb parameter Using ab initio calculations we calculate the eective static Coulomb coupling constant, (usually treated as a free parameter), where we account for extrinsic and intrinsic screening within the random phase approximation (RPA). On this basis, we show that a purely Coulomb driven superconducting phase with an order parameter that has op- posite signs in dierent valleys [169] is not favored in MoS 2 . We nd this to be true independently of the dielectric environment of the substrate. For the scenario of phonon mediated SC [133, 135, 136], we show that the phonon mediated electron-electron attrac- tion generally overcomes the Coulomb repulsion when a Lifshitz transition takes place and additional Fermi pockets become available. The intrinsic screening of TMDC's at their su- perconducting transition is shown to be typically so large that it renders external substrate screening rather unimportant for the SC transition despite the atomic scale proximity of the substrate. In the following we specically examine the in uence of the static screened Coulomb interactions W 0 (q) = W (q;! = 0), on superconducting pairing, which depends on the electron doping level and on the dielectric environment. Using the fully screened interaction matrix ^ W 0 (q), we can compute the full Coulomb coupling constant = 1 N(E F ) X kk 0 W kk 0( k E F )( k 0E F ); (5.14) which is the Fermi surface average of the screened Coulomb interaction, including scattering processes with initial statesfk; k 0 g and nal statesfk 0 ; kg (i.e. q = k k 0 ). In Eq. (5.14) ^ W (q)!W kk 0 =h(k 0 ; k)j ^ W (kk 0 )j(k; k 0 )i has been transformed from the orbital basis to the band basis, and we only consider the lowest conduction band for ^ W (q) since it is the only band that crosses the Fermi level for the electron doping concentrations considered here. The resulting eective Coulomb coupling constants are shown in Fig. 5.10(a) for dierent dielectric environments of the MoS 2 layer and in dependence of the electron doping concentration. In the low-doping regime,x. 0:07, where only two Fermi pockets around K and K 0 are present, the coupling is renormalized by up to 30% via external screening. In contrast, at higher doping concentrations is clearly much less sensitive to its dielectric environment, and variations of due to external screening are limited to 10%. If we account for the multi-valley structure of the Fermi surface (see Fig. 5.5), is no longer a simple scalar but becomes a matrix in the electronic valleys. To further investigate the eect of the dielectric environment, we discuss this matrix structure of . For low doping concentrations, where only the two valleys around the K points are occupied, we obtain the following structure low = intra inter inter intra ; (5.15) 112 CHAPTER 5. MOS 2 CALCULATIONS 5.5. SUPERCONDUCTIVITY 0 0.02 0.04 0.06 0.08 0.1 0.12 0.1 0.15 0.2 0.25 0.3 electron doping concentration x Averaged Coulomb potential μ ε = 1 ε = 10 ε = 20 ε = 50 ε =∞ (a) 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0 0.05 0.1 0.15 electron doping concentration x Entries of μ matrix ε = 1 ε = 10 ε = 20 ε = 50 ε =∞ (b) filled: μintra empty: μinter Figure 5.10: Dependence of Coulomb coupling constants on the electron doping concen- tration x, subject to various " sub of the dielectric environment. (Left) The full coupling constant is shown. (Right) We present the values of the 2 2-matrix (cf. Eq. 5.15) for low doping concentrations. where the statesfk; k 0 g in Eq. (5.14) are in the same valley for intra while they are in dierent valleys for inter . The sum of all matrix elements yields the total coupling constant . A comparison of external screening eects on intra- and inter-valley Coulomb scattering [Fig. 5.10(b)] shows that essentially only the intra-valley scattering is aected by the dielectric environment. These observations can be explained intuitively. External screening is most eective when the separation ( 1=q) of the interacting charges inside the monolayer is larger than the distance 1 2 d to their image charges in the environment but smaller than the internal Thomas-Fermi screening length 1=q TF , i.e. forq TF <q< 2 d . As a consequence, the in uence of the substrate weakens as soon asq TF & 2 d . Using the eective thickness of d 9:1 A, a Thomas-Fermi wave vector q TF = 2e 2 N(E F )=(A" 1 (q TF )), and a background dielectric constant on the order of " 1 (q)" 1 = 9:3 for q > 2 d , we nd that the substrate in uence is minor as soon as the density of states at the Fermi level exceeds N(E F ) 0:19 /eV per unit cell. In MoS 2 , we have N(E F ) 0:4 eV 1 and N(E F ) 2 eV 1 for low (x < 0:07) and high electron doping concentrations (x > 0:07), respectively, which means that substrate in uence is weak, especially in the regime of high doping concentrations. However, for suciently low doping concentrations, the scattering inside the same K- or K 0 -valley can be controlled via the substrates. In the framework of Eliashberg theory [170], the Allen-Dynes formula [171] yields an estimate of the critical temperature, T c = ~! log 1:2k B exp 1:04(1 +) (1 0:62 ) ; (5.16) which accounts for the competition of the phonon driven attractive interaction (entering via the eective coupling strength and the typical phonon frequency ! log ) with the 113 5.5. SUPERCONDUCTIVITY CHAPTER 5. MOS 2 CALCULATIONS repulsive Coulomb interaction expressed by the Morel-Anderson parameter [113]. The phononic parameters for MoS 2 have been calculated in Refs. [135,136]. The coecient that describes the Coulomb repulsion is obtained using the formula given by Morel and Anderson [113] for the retarded Coulomb potential, = 1 + ln[ E F ! log ] : (5.17) In Fig. 5.11, we plot the dependence of on the electron doping level for freestanding MoS 2 and MoS 2 embedded in a perfect metallic environment. For free-standing MoS 2 we observe a decrease of from > 0:25 to . 0:15 for x. 0:07, which is caused by the corresponding decrease in and the decrease in the phonon frequency! log (see Ref. [136]). At larger electron doping concentrations, is basically constant with 0:13. For MoS 2 embedded in a metallic environment, shows essentially the same trend with the only dierence in comparison to the free-standing case being a slight reduction of , particularly at low doping, which could shift the onset of the superconducting phase to a lower doping concentration than the critical concentration for the freestanding layer. A signicantT c is only reached when the exponent in Eq. (5.16) is close to1 or larger, especially when the electrons mainly couple to acoustic phonons and thus ! log is rather small, e.g. T c & ~! log 1:2k B e 2 . To achieve this, > 3 has to be realized for the range of 0:1< < 0:3 found here. From the comparison of and =3 in Fig. 5.11, we see that a signicant T c (as occurring for > 3 ) can only be observed once x& 0:07, i.e. when 0 0.02 0.04 0.06 0.08 0.1 0.12 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 electron doping concentration x Coulomb and el.-phon. coefficient μ ∗ , ε sub = 1 μ ∗ , ε sub =∞ λ/3 Figure 5.11: Retarded Coulomb potential used in Eliashberg theory and eective electron-phonon coupling strength in dependence of the electron doping concentrationx. for the freestanding monolayer (" sub = 1) and in the presence of a metallic environment (" sub =1). is scaled by a factor of 1=3. 114 CHAPTER 5. MOS 2 CALCULATIONS 5.6. CONCLUSIONS both valleys in the conduction band are occupied by electrons. We thus conclude that the frequent use of a constant for the Coulomb pseudopotential, e.g. = 0:13 [112,119,135], is not sucient in the case of electron doped MoS 2 to describe the in uence of the Coulomb interaction directly at the transition to the superconducting phase. However, the drop in the critical temperature of TMDC's [150,172{175] when going from the bulk or multilayer-system to a monolayer cannot be caused by enhanced Coulomb interactions, because the values of the electron-phonon coupling are much larger than the 0:13, which we nd in the region of optimal doping independently of the dielectric environment of the MoS 2 monolayer. 5.6 Conclusions We found that the low energy dynamical screening in MoS 2 is controlled by both inter- and intra-valley scattering processes. These give rise to plasmons with a square root dis- persion at small q and linear dispersion for higher momentum transfers which connect separate valleys on the Fermi surface. In general, inter-valley plasmon modes are observ- able, although their oscillator strengths are strongly reduced in comparison to zone center modes. Due to the multi-orbital character of the wave functions and spin-orbit coupling, which leads to spin-valley coupling in monolayer TMDC's, not all inter-valley scattering processes are allowed. As a consequence of spin-valley coupling some inter-valley plasmon modes are shifted and gapped out, while the 2=3 rotation symmetry protects certain low energy modes at M. We speculate this selective gapping out of collective modes could have consequences for the realization of many-body instabilities towards superconducting or charge density wave phases in monolayer TMDC's. The microscopic description of the static Coulomb interactions in the electron doped monolayer MoS 2 developed here reveals a clear decrease of the retarded Coulomb poten- tial with increasing doping, which renders the frequent use of a constant, doping and material independent questionable. Comparing the values for the electron-phonon in- teraction in Ref. [136] with the retarded Coulomb potential and the valley decomposed electron-electron interaction coupling constants presented here, we conclude that the su- perconductivity in MoS 2 is electron-phonon driven and has its onset at electron doping levels for which both valleys in the conduction band are occupied. The eects of substrates turn out to be relatively small, at least around optimal doping, and we nd that the exper- imentally observed reduction of the critical temperature upon approaching the monolayer limit [150] is not caused by enhanced Coulomb interactions, i.e. lack of screening as the dimensionality of the material is reduced. This conclusion should be generally applica- ble also to other superconducting 2d materials such as NbSe 2 [175] and particularly the electron doped TMDC's like WS 2 or MoSe 2 [176] because of their similar electronic and phononic structure. 115 Chapter 6 Conclusions and Outlook A quantitative description of the Coulomb interactions was developed for two-dimensional superconducting materials, enabling us to compare intrinsic with external screening eects, such as those due to substrates (see Sec. 5.3). Using the example of a doped monolayer of MoS 2 embedded in a tunable dielectric environment, we demonstrated that the in uence of external screening on the static Coulomb interaction is limited to a length scale, bounded from below by the eective thickness of the quasi two-dimensional material and from above by its intrinsic screening length. We also calculated the Morel-Anderson pseudo potential, traditionally, used to account for the static Coulomb repulsion in superconducting systems, using ab initio techniques. We applied this to the case of MoS 2 and found that unconven- tional, static, Coulomb driven superconductivity cannot be induced in MoS 2 by tuning the substrate properties alone (see Sec. 5.5). Using the same quantitative description we were able to study the rich phenomenology of plasmonic excitations in the dichalcogenides as a function of doping. The many-body polarization, the dielectric response function and electron energy loss spectra were cal- culated using an ab initio based model involving material-realistic Coulomb interactions, band structure and spin-orbit coupling. Focusing on the representative case of MoS 2 , a plethora of plasmon bands were observed, originating from scattering processes within and between the conduction or valence band valleys. We discussed the resulting square-root and linear collective modes, arising from long-range versus short-range screening of the Coulomb potential. We showed that the multi-orbital nature of the bands and spin-orbit coupling strongly aects inter-valley scattering processes by gapping certain two-particle modes at large momentum transfer. Knowledge of the dynamically screened Coulomb interaction was then applied to study the eect of electron-plasmon interactions on superconductivity on such low dimensional materials. We applied SC-DFT to study this in uence and was able to show that the McMillan equation can be derived from the framework of SC-DFT (see Sec. 3.2.3) when including 'high' energy interactions in Cooper-pair formation. The McMillan equation is 116 CHAPTER 6. CONCLUSIONS AND OUTLOOK therefore captured in a G 0 W 0 -theory. We extended the energy space formulation of SC- DFT to include G 0 W corrections and presented a McMillan-like equation for the super- conducting critical temperature for a combined phononic and plasmonic superconducting state. Furthermore, we derived an analytic approximation to the plasmonic Eliashberg function for a 2d system (Sec. 3.2.4), thereby showing that the full energy dependent Eliashberg function is needed to capture the relevant physics of plasmonic induced super- conductivity in 2d. We then utilized model calculations to address general consequences of including plasmonic eects in superconductivity and formulated experimental guidelines to nd plasmon induced enhanced critical temperatures (see Ch. 4). Seeing as this is a theoretical thesis it is very fortunate that we were also able to corrob- orate some of these conclusions with experimental results. Specically, a recent experiment by the group of Dr. Wei Wu explored the plasmonic excitation frequency of gold nano- particle dimers placed close to each other and coated in dielectric material. They invented a technology to fabricate gap plasmonic structures with sub-nanometer resolution, high reli- ability and high throughput, deterministically using collapsible nano-ngers. That enabled us to experimentally investigate the eects of gap size and tunneling barrier height on the plasmonic properties of the system. These systems are of interest since such nanostructures have many unique properties, for example, they can incorporate light energy into a small volume at nanometer scale [177{179]. That makes them a promising platform for various applications including optical communication [180], disease diagnosis [181], and chemical sensing [182]. The technology from Dr. Wu's group, is based on collapsible nano-ngers fabricated using nanoimprint lithography (NIL), reactive-ion etch (RIE) and atomic-layer deposition (ALD). The SEM (scanning electron microscope) images in Fig. 6.1 shows the gap plasmonic nanostructures they fabricated. For detailed information please see our publication Ref. [183]. Once the nano-particles were fabricated and brought close to one another, the system was excited by an external light source while scanning through frequencies. A known result is that the system responds strongly at the plasmon frequency, resulting in a strongly amplied eld at the midpoint between the nanoparticles [184, 185], as indicated in Fig. 6.2. We noticed that the resonant frequency shifted depending on the amount of dielec- tric material separating the nano-particles as well as the type of dielectric material used. This phenomenon turned out to be easily understood in terms of environmental screening brought on by the coating materials. Leading to the expression max =A D n + 0 ; (6.1) where max is the incident wavelength that maximally excites the system,D is the distance between the nano-particles,n the index of refraction of the dielectric material used, 0 the USC Department of Electrical Engineering 117 CHAPTER 6. CONCLUSIONS AND OUTLOOK Figure 6.1: (a) SEM image of nano-ngers before collapse. (b) SEM images of nano-ngers after collapse. The inset is magnied image of the same nano-ngers from a dierent viewing angle. The diameter and height of each nger are 50 nm and 350 nm, respectively. Image taken from Ref. [183] Figure 6.2: Schematic of controllable hot spots created using nano-ngers. A pair of metal- lic nanoparticles with ALD-coated dielectric layer contact each other when exible nano- ngers beneath them collapse. Gap size is dened by twice the dielectric layer thickness and hot spot is right at the gap center. Image taken from Ref. [183] 118 CHAPTER 6. CONCLUSIONS AND OUTLOOK 6.1. OUTLOOK plasmon frequency of a standalone nanoparticle and A is a geometric parameter that is constant for systems of this type. A comparison between the theoretical expression and experimental results is shown in Fig. 6.3, and it clearly explains the observed trends very well. 0.0 0.5 1.0 1.5 2.0 2.5 3.0 D/n(nm) 550 600 650 700 λ max (nm) Red-shift of plasmonic frequency WO 3 TiO 2 SiO 2 Figure 6.3: The measured red shift for dierent values ofD=n from 3 dierent nanoparticle composites compared to the theoretical expression given in Eq. 6.1. The dots represent experimental data while the solid lines are theoretical predictions. Image taken from Ref. [183] This experiment provided a valuable validation of one of the major conclusion of this thesis that environmental screening strongly aects plasmonic properties, and more specif- ically that increased environmental screening leads to a suppression of the plasmon fre- quency, as was shown in Fig. 4.5. I hope that a future group will perform an experiment measuring the critical temperature of a plasmon induced superconductor in order to test the trends concluded from our SC-DFT studies. 6.1 Outlook First and foremost our next step is to apply the SC-DFT framework, here developed, to study the in uence of electron-plasmon coupling on real material systems, such as mono- layer MoS 2 . Doing this requires us to relax the assumption of a gapped electron-plasmon coupling and instead make use of the full energy dependent plasmonic Eliashberg function. We have already solved the technical diculties for this (as discussed earlier) but bench- marking the implementation is outstanding at this time, which is why no prelimary results are shown. Once the full energy dependent SC-DFT implementaion is benchmarked, we 119 6.1. OUTLOOK CHAPTER 6. CONCLUSIONS AND OUTLOOK can combine the ab initio calculation based phononic Eliashberg function and the realistic plasmonic Eliashberg function in the SC-DFT kernel to do a detailed study of the interplay between those modes in the superconducting phase. The importance of understanding and properly treating the Coulomb interaction in physical systems was a key motivation behind most of the work done in this thesis. As is generally the case in scientic work the progress made here has led to many other questions concerning the Coulomb interaction. Continuing to study these questions forms the other major outlook to future work. A particular interesting direction we are now, and will continue to explore, is the cal- culation of plasmonic properties of systems in real-space and how environmental screening aects them. Real-space calculations are interesting since it allows us to study systems that lack translation invariance and therefore we can study nite systems, or systems with patterned dielectric substrates as shown in Fig. 6.4. Figure 6.4: Cartoon of a 2d system between dielectrics with a junction of dierent dielectric materials. The general machinery for these types of calculations very naturally extends from the work in this thesis. Specically we use a tight-binding model for the electronic Hamiltonian, with s-like orbitals centered at each atomic site (although it is trivial to extend to other orbitals). Once the single-particle electronic states are known the RPA can be used to calculate the polarization , from which the dielectric function and plasmon modes are calculated, as before. We calculate the plasmon eigenmodes according to the method described in Ref. [187]. The dielectric matrix in real space is given by, hrj ^ "(!) r 0 = r r 0 Z hrj ^ V C r 00 r 00 ^ (!) r 0 dr 0 where ^ V C is the Coulomb interaction and ^ the polarization. Wang et al. [188] describes how to visualize the plasmon modes in real-space, by a spectral decomposition of the the real-space form of the RPA polarization is described in Ref. [186] 120 CHAPTER 6. CONCLUSIONS AND OUTLOOK 6.1. OUTLOOK dielectric matrix: ^ "(!) = X n n (!)j n (!)i: The eigenfunction with corresponding eigenvalue that yields the largest value of Im[1= n (!)], is taken as the dominant plasmon mode at the given!. Call this eigenvectorj 0 i, then the real space plasmon mode can be visualized with Re[hrj 0 i] as demonstrated in Fig. 6.5. As already discussed, the standard experimental way to measure plasmonic excitations is through EELS measurements. We are able to recover the EELS spectrum from the real-space calculations through the Fourier transform of the dielectric matrix, hqj ^ "(!)jqi = 1 (2) d Z d d r Z d d r 0 hrj ^ "(!) r 0 e iq(rr 0 ) ; since the EELS spectrum is simply given by Im[1=hqj ^ "(!)jqi]. In Fig. 6.6 the EELS spectra (and real and imaginary polarizations) are shown for the same 20 3 site system as before, with dierent values of environmental screening used in the evaluation of the Coulomb interaction. The Fourier transformed results show very similar trends to those seen forq-space calculations (as expected) and we see a strong suppression of the plasmon dispersion with increased environmental screening, same as before. The goal is to use the real-space calculations to accurately describe systems where the dielectric substrate changes at some boundary, or perhaps at several boundaries (see Fig. 6.4). Such heterostructures have shown promise for the construction of type-II heterojunc- tions [189]. The possibility of engineering boundary modes has interesting applications to 0 10 20 30 40 r x −4 −2 0 2 4 6 8 r y Charge distribution −0.15 −0.10 −0.05 0.00 0.05 0.10 0.15 Figure 6.5: Real space visualization of plasmon mode for a 20 3 site system with tight- binding hopping parameter, t = 0:5, at ! = 4eV. 121 6.1. OUTLOOK CHAPTER 6. CONCLUSIONS AND OUTLOOK Figure 6.6: EELS spectra (and real and imaginary polarizations) for a 20 3 site system with tight-binding hopping parameter, t = 0:5. The environmental screening of " env = 1 was used for the top case and" env = 5 for the bottom. The yellow line traces the calculated plasmon frequency in all cases. Notice the strong suppression of the plasmon dispersion for the higher " env case. plasmonic circuits [190]. The possibility of topological order in the plasmonic state will also be investigated within this framework where the localization of plasmons at a dielectric substrate interface occurs. The other big advantage of the real-space formulation is that we can handle clusters of atoms that do not exhibit any regular periodicity. These include limited systems or quasicrystals. Quasicrystals are structures that are ordered but not periodic and have been shown to exist widely in many metallic alloys [191]. Such systems are under-investigated since most solid state computational tools are designed for periodic systems. There are clearly still a lot to learn concerning the Coulomb interaction and the screen- ing of it, that has deep reaching consequences in the observed properties of matter. Under- standing these eects better not only increases our understanding of the world around us but also increases our ability to engineer material properties to better t our ever expanding technological needs. 122 Appendix A Theory Overview Sections A.1 Hartree-Fock Equations The following description of the Hartree-Fock approximation comes from Ref. [90]. The Hartree-Fock approximation is a mean eld approach to the interaction of electrons in a system. Consider a system described by the Hamiltonian H = X i i c y i c y i + 1 2 X ijkl V ijkl c y i c y j c y k c y l ; where c y j (c y k ) is the electron creation (annihilation) operator and V the interaction matrix element. The mean eld approach dictates that we only take into account terms where n ij =hc y i c y j i is non-zero. We thus write the operators in second term of the Hamiltonian in terms of a deviation from an average c y i c y j c y k c y l =c y i (c y j c y k hc y j c y k i)c y l +c y i c y l hc y j c y k i: As long as j6=l we can commute c y l with the bracket in the rst term, and we can again write c y i c y l as it's deviation from average, c y i c y l = (c y i c y l hc y i c y l i) +hc y i c y l i: Plugging this into the above expression we get c y i c y j c y k c y l = (c y i c y l hc y i c y l i)(c y j c y k hc y j c y k i) +c y i c y l hc y j c y k i +c y j c y k hc y i c y l ihc y i c y l ihc y j c y k i c y i c y l n jk +c y j c y k n il n il n jk (A.1) where we neglect the rst term since it goes as the square of the deviation from average and is therefore assumed small. The Hartree approximation for the interaction term is therefore given by: V Hartree = 1 2 X ijkl V ijkl n jk c y i c y l + 1 2 X ijkl V ijkl n il c y j c y k 1 2 X ijkl V ijkl n il n jk (A.2) 123 A.2. ATTRACTIVE ELECTRON-ELECTRON INTERACTION FROM PHONON SCATTERING APPENDIX A. THEORY OVERVIEW SECTIONS The Fock term is still missing and comes from the fact that we also get a contribution to the Hamiltonian (within the mean eld approximation) ifn ik orn jl is non-zero. Doing the same with these terms as we did before yields the Fock term: V Fock = 1 2 X ijkl V ijkl n ik c y j c y l 1 2 X ijkl V ijkl n jl c y i c y k + 1 2 X ijkl V ijkl n ik n jl : (A.3) The Hartree-Fock equations then come from diagonalizing the mean eld Hamiltonian given by H HF =H 0 +V Hartree +V Fock ; whereH 0 = P i i c y i c y i . The combined interaction termV Hartree +V Fock is what is meant by v Hartree in Eq. 2.30. Furthermore, if we are considering a crystal system with translational invariance, the expectation valuehc y k c y k 0 i is diagonal and we get the eigenvalues of the Hartree-Fock Hamiltonian as E HF k = k +V (0)N X k 0 V (k k 0 )n k 0: Usually, the second term is canceled by an equal and opposite term coming from the lattice or external potential. A.2 Attractiveelectron-electroninteractionfromphononscat- tering The following was taken from Ref. [90]. If we assume a 'jellium-model' for the phonons (similar to how we derived the plasmon dispersion in Section 2.3.1) we nd the dispersion- less phonon modes q = = s Ze 2 N 0 MV ; where ion =N=V is the ion density and Z is the atomic number for an ion with mass M. Furthermore, as derived in most introductory solid state courses, one nds the electron- phonon coupling constant for these jellium phonons are given by, 1 V jg q j 2 = ~e 2 2 0 q 2 = 1 2 W (q) : We thus get the 'bare' phonon-mediated interaction 1 V jg q j 2 D 0 (q;iq n ) =W (q) (iq n ) 2 2 : 124 Ω ω −2 0 2 V eff (q)/W(q) Figure A.1: Eective electron-electron interaction with Coulomb repulsion and phonon- mediated interactions included. Thus to a very low degree we can write the eective electron-electron interaction coming from the Coulomb repulsion and phonon-mediated interaction as, V eff (q;iq n ) =W (q) +W (q) 2 (iq n ) 2 2 =W (q) (iq n ) 2 (iq n ) 2 2 ; which we can analytically continue to the real frequency axis to get: V eff (q;!) =W (q) ! 2 ! 2 2 +i (A.4) In Fig. A.1 we show the resulting eective electron-electron interaction from Eq. A.4, notice the negative (attractive) interaction for!< . We thus see that the electron-phonon interaction could lead to an eective attraction between electrons. Appendix B Computational Parameters B.1 MoS 2 Ab initio Calculations The complete input le for the Quantum ESPRESSO 'self-consistent eld' calculation (scf.in) is shown below. This le is given as input during the execution of the "pw.x" 125 B.1. AB INITIO CALC. APPENDIX B. COMPUTATIONAL PARAMETERS program. 1 &c o n t r o l 2 c a l c u l a t i o n =' scf ' , 3 p r e f i x ='MoS2' , 4 pseudo dir = ' . . / . . / pp / ' , 5 outdir = './ ' , 6 t p r n fo r = . true . , 7 t s t r e s s = . true . , 8 e t o t c o n v t h r = 1.0 d5 9 f o r c c o n v t h r = 1.0 d4 10 / 11 &system 12 ibrav = 4 , 13 celldm (1) = 5 . 9 , 14 celldm (3) = 4.237288136 , 15 nat= 3 , 16 ntyp = 2 , 17 ecutwfc = 40.0 , 18 smearing = 'mp' , 19 occupations = ' smearing ' 20 degauss = 0.02 , 21 t ot c ha rg e = 0.100 , 22 nbnd = 23 , 23 la2F = . true . , 24 / 25 &e l e c t r o n s 26 d i a g o n a l i z a t i o n = ' david ' 27 mixing mode = ' plain ' 28 mixing beta = 0.7 29 conv thr = 1.0 d12 30 / 31 ATOMIC SPECIES 32 Mo 95.960000 Mo. pwmt fhi .UPF 33 S 32.065000 S . pwmt fhi .UPF 34 ATOMIC POSITIONS c r y s t a l 35 Mo 0.000000000 0.000000000 0.000000000 36 S 0.333333333 0.666666666 0.116187965 37 S 0.333333333 0.666666666 0.116187965 38 K POINTS AUTOMATIC 39 24 24 1 0 0 0 The same input parameters were used in the non-self consistent eld calculations needed in order to run the EPW package that solves the Eliashberg equations. The K POINTS setup was naturally changed to CRYSTAL and the exact 576 k-points specied, as is standard in these calculation schemes. The input le for the Quantum ESPRESSO DFPT phonon calculations (ph.in) is show below, that is given as input to the "ph.x" program. We used a 12x12 grid for the phonon 126 APPENDIX B. COMPUTATIONAL PARAMETERS B.1. AB INITIO CALC. Eliashberg function calculations as well as for the input to the Eliashberg solver. For the plot of the phonon dispersion shown in Fig. 5.3 we used a 24x24 grid in order to increase resolution. 1 &inputph 2 p r e f i x = 'MoS2' , 3 f i l d y n = 'MoS2 . dyn ' , 4 amass (1) = 95.940 , 5 amass (2) = 32.065 , 6 outdir = './ ' 7 electron phonon =' interpolated ' , 8 l d i s p = . true . , 9 trans = . true . , 10 f i l d v s c f = ' dvscf ' , 11 nq1=12, 12 nq2=12, 13 nq3=1, 14 tr2 ph = 1.0 d17 15 / The input parameters for the EPW code (executed with "epw.x") is mostly obtained from the inputs to the preceding ab initio calculations as well as the specics of the Wannier construction (discussed next). The remaining parameters were as follows (for the neglected ones the default values were used): 1 f s t h i c k = 0.2 ! eV 2 degaussw = 0.01 ! eV 3 4 degaussq = 0.2 ! meV 5 6 e l i a s h b e r g = . true . 7 nqstep = 500 8 c o n v t h r i a x i s = 1.0 d4 9 10 mp mesh k = . true . 11 nkf1 = 360 12 nkf2 = 360 13 nkf3 = 1 14 15 nqf1 = 180 16 nqf2 = 180 17 nqf3 = 1 The choice to use specically 360x360 k-grid points and 180x180 q-grid points were made based on the convergence tests discussed in the PhD thesis of Dr. Harry Fisher Ref. [154]. The phonon-broadening was also chosen to match the parameter used by Dr. Fisher but the electron smearing ( elec ) was chosen so that the resulting density of states appeared reasonable. The Fermi surface thickness was found to not in uence the calculation results as long as it is chosen bigger than a certain threshold (roughly 4 elec ). 127 B.2. MOS 2 WANNIER APPENDIX B. COMPUTATIONAL PARAMETERS B.2 MoS 2 Wannier The input le for the Wannier construction that supplies the minimal basis is shown below. This le is given as input to the "wannier90.x" program. We used the ab initio results from above to do the Wannier construction. 1 num bands = 23 2 num wann = 3 3 4 begin p r o j e c t i o n s 5 Mo: dxy , dz2 , dx2y2 6 end p r o j e c t i o n s 7 8 num iter = 1000 9 d is n um i te r = 10000 10 dis win min = 0.0 11 dis win max = +6.0 12 d i s f r o z m i n = 0.3 13 dis froz max = 3.75 14 15 bands plot = . true . 16 begin kpoint path 17 G 0.00 0.00 0.00 M 0.00 0.50 0.00 18 M 0.00 0.50 0.00 K 0.33 0.33 0.00 19 K 0.33 0.33 0.00 G 0.00 0.00 0.00 20 end kpoint path 21 bands plot format = gnuplot 22 23 dos = . true . 24 dos kmesh = 360 360 1 25 dos adpt smr = f a l s e 26 dos smr type = ' gauss ' 27 dos smr fixed en width = 0.01 28 29 begin u n i t c e l l c a r t 30 bohr 31 5.900000 0.000000 0.000000 32 2.950000 5.109548 0.000000 33 0.000000 0.000000 25.00000 34 end u n i t c e l l c a r t 35 36 begin atoms frac 37 Mo 0.000000000 0.000000000 0.000000000 38 S 0.333333333 0.666666666 0.116187965 39 S 0.333333333 0.666666666 0.116187965 40 end atoms frac 128 B.3 2d Square Lattice SC-DFT calculations For the evaluation of the bare Coulomb interaction we use = 1:5 A 1 and" 1 = 25 as well asd = 5 A to evaluate" rest (q). The polarization function q (!) is calculated on a 120120 q-grid based on a 120 120k-grid involving 400 frequency points between 0 and 1 eV. The broadening parameter i0 + is set to 10 meV (20 meV and 50 meV result in identical trends, however, with slightly reduced electron-plasmon couplings). The plasmonic Eliashberg function is evaluated on a frequency grid using 800 points between 0 and and 1 eV. All (x) functions are approximated by Gaussian functions with a smearing of 7 meV. The SC-DFT gap equation is solved on an energy grid ranging from the lower to the upper end of the band width. We use 1200 grid points which are distributed logarithmically within a window of0:1 eV around the Fermi level (3=4 of all points). All other points are distributed linearly. B.4 MoS 2 Plasmon calculations For the plasmon calculations from Section 5.4, we apply a centered Monkhorst-Pack 720 720 k-grids and use a broadening of i = 0:0005i. The doping concentration is ad- justed by rigid shifts of the Fermi energy, which change the Fermi functions accordingly. All calculations are carried out forT = 0 K. Furthermore, we restricted the 1 and 2 sum- mations to the partially occupied band only in order to avoid double counting problems within the denition of the total polarization function and not to overload the resulting plots. Appendix C Code Implementation Benchmarks In order to benchmark our single-band SC-DFT implementation, we performed several tests, comparing results from our code to available data and known literature results. 129 C.1. SC-DFT - STEP-LIKE DENSITY OF STATES APPENDIX C. BENCHMARK C.1 SC-DFT - Step-Like Density of States −3.0 −2.5 −2.0 −1.5 −1.0 −0.5 0.0 0.5 1.0 E f [ω E ] 10 15 20 25 30 35 40 T c Eigenvalue method 1000 steps 750 steps 500 steps PRB 88, 014514 (2013) 0.01 0.02 0.03 0.04 0.05 0.06 DOS Figure C.1: SC-DFT Benchmark I. Comparison between ourT c (lines) and data published in Ref. [111] (black dots) for a step-like density of states (gray line) including the eects of electron-phonon coupling and static electron-electron repulsion. We show results obtained by using the iterative method to solve the gap equation with dierent numbers of steps allowed before terminating the loop (dashed lines) and by using the eigenvalue method (solid line). The eigenvalue method appears as the limit of an increasing number of iteration steps. We compare our implementation to results obtained by Akashi and Arita published in Ref. [111]. In analogy to this reference, we use an Einstein phonon with frequency ! E 50 meV with an eective electron-phonon coupling of = 1 and an eective electron- electron repulsion of = 0:5. The underlying density of states is given by a step-like func- tion, as described in Ref. [111] with a band-width of 40 eV. Akashi and Arita determined T c by self-consistently solving the gap equation and setting T c to the highest temperature for which they could nd a non-trivial solution. Such an implementation requires the use of a parameter that limits the number of allowed iterations to obtain a solution to the gap equation before terminating the process. In our implementation we cast the gap equation into an eigenvalue problem and ndT c as the temperature for which the leading eigenvalue is one. In Fig. C.1 we show results we obtained using both methods (showing the dierence due to changes in the parameter) along with the results from Ref. [111]. Besides decreased critical temperatures T c above E f = 0 the agreement between our and the reference data 130 APPENDIX C. BENCHMARK C.2. SC-DFT - MOS 2 MONOLAYER is very good for the case where 500 iterations steps are allowed. Like in Ref. [111] we nd a maximum of T c in the vicinity of E F 0:3! E and similar trends above and below this point. The minor dierences between the two data sets trace back to slightly dierent energy grids, bandwidth cut-os, and dierent tolerance settings. C.2 SC-DFT - MoS 2 Monolayer As a second benchmark, we compare critical temperatures obtained by evaluating the Allen- Dynes equation, using input data from density functional perturbation theory and from our SC-DFT implementation for electron-doped monolayers of MoS 2 . We use the eective frequencies ! log = ! ph and eective electron-phonon couplings ph from Ref. [136] to construct Einstein-phonon models for a range of dierent doping levels. Together with the electronic dispersion of the occupied band of MoS 2 , which we obtain from a corresponding Wannier construction, we can solve the SC-DFT gap equation given in Eq. 3.61 of the main text. To this end, we use similar k meshes and broadenings as in Ref. [136] and x the resulting = 0:15. In Tab. C.1 we list the resulting critical temperatures from the reference (T Allen-Dynes c ) and the solution of the SC-DFT gap equation (T SC-DFT c ). For ph < 2 we nd similar critical temperatures from both approaches. For increased eective electron-phonon couplings, the resulting T c dier, as the Allen-Dynes equation is known to underestimateT c in the strong-coupling regime [171]. Apart from this, small dierences occur due to slightly dierent bandwidths and the applied rigid-shift approximation to describe the doping in the T SC-DFT c data. Considering these circumstances, the agreement between both approaches is good. Table C.1: SC-DFT Benchmark II. Comparison of critical tempera- tures for electron-doped monolayers of MoS 2 , obtained by evaluating the Allen-Dynes equation using density functional perturbation theory input from Ref. [136] and by solving the energy-dependent SC-DFT gap equation. doping a 0:075 0:087 0:100 0:112 0:125 ph 0:821 1:236 1:920 3:089 7:876 ! ph [eV] 0:029 0:025 0:021 0:016 0:009 N 0 [eV 1 ] b 0:440 0:569 0:690 0:857 0:895 T Allen-Dynes c [K] 8:01 15:57 21:12 22:11 15:51 T SC-DFT c [K] 5:40 12:76 22:44 35:07 39:60 a given in additional electrons per unit cell b given in states/spin/eV/unit cell 131 C.3. THREE-BAND MODEL FOR TMDC'S APPENDIX C. BENCHMARK C.3 Three-band model for TMDC's Up to now, there have been basically two theoretical approaches available to study the plasmonic physics in TMDC's. On the one side there are models combining eective kp descriptions of the quadratic electronic bands around the band gap with an evaluation of the dielectric function within the random phase approximation (RPA) [164, 165]. On the other side, there are RPA descriptions based on full density functional theory (DFT) calculations, which include realistic single particle band structures describing the complete Brillouin zone [138,166,192{195]. Here, we add a third approach by utilizing a material-specic low-energy model Hamil- tonian derived from ab initio calculations for the undoped material as the basis for the eval- uation of dynamical response functions in the electron and hole doped situations. Thereby we gain the possibility to accurately calculate the polarization as well as the screening functions for the whole Brillouin zone, which enables us to study plasmons at arbitrary momenta as described in the main text. Figure C.2: (Color online) Plasmon dispersions from our ab initio based model (dots) in comparison with analytic kp-models by Scholz et al. [164] and Kechedzhi et al. [165]. In bothkp-models the authors used a simple 1 q constantly screened Coulomb potential while we take the full background screened Coulomb interaction into account (green dots). When we switch to the same constantly screened Coulomb potential our results match those of the other authors as seen by the red and blue dots respectively. The two kp-models produce dierent plasmon dispersions due to the ways they handle the spin dependencies. In the data by Kechedzhi et al. [165] both spin components and their coupling are included while the data of Scholz et al. [164] include a single spin component only. 132 APPENDIX C. BENCHMARK C.3. THREE-BAND MODEL FOR TMDC'S In Fig. C.2 we compare the resulting plasmon dispersions of our ab initio based model (dots) to the models by Scholz et al. [164] and Kechedzhi et al. [165] (lines) for hole doping with a carrier density of 2 10 13 cm 2 in the case of MoS 2 . In the kp models a simple constantly screened Coulomb interaction of the form U(q)/ 1 q with " (q) = = 5 is used. As long as we use the same constantly screened Coulomb interaction to evaluate the dielectric function we end up with nearly identical plasmon dispersions compared to those derived from thekp-models. However, by including the full ab initio derivedq dependent dielectric function for the background screening in the undoped case, which arises due to the two dimensional geometry and the excluded bands, we nd strongly reduced plasmon energies (green dots). In Fig. C.3 we compare our method to full ab initio results for NbS 2 [138], which is a metal in its ground state. Once again, plasmon dispersions derived form EELS data are shown. Although, the resulting plasmon dispersion are on the eV range (for which our model is actually not set-up), we nd a remarkable agreement with dierences on the order of 100 meV to 200 meV ( 10% to 20%). Hence, we conclude that neglecting the material-specic dielectric function " (q) within the minimal three-band model is a crucial approximation leading to non-realistic plasmonic properties. Figure C.3: (Color online) Plasmon dispersions for undoped NbS 2 from (left) full ab initio calculations from [138] and (right) our ab initio based model. Here we use the macroscopic dielectric function. 133 Bibliography 1 Nobelprize.org. Nobel laureates 1977. https://www.nobelprize.org/nobel_prizes/ physics/laureates/1977/index.html, 1977. Accessed: 2018-06-02. 2 P. W. Anderson. More is dierent. Science, 177(4047):393{396, 1972. 3 Nobelprize.org. Nobel laureates 1933. https://www.nobelprize.org/nobel_prizes/ physics/laureates/1933/, 1933. Accessed: 2018-06-02. 4 E. Pavarini and E. Koch. Emergent Phenomena in Correlated Matter. Forschungszen- trum Julich, 2013. 5 James P. Sethna. Statistical Mechanics: Entropy, Order Parameters and Complexity (Oxford Master Series in Physics). Oxford University Press, 2006. 6 K. S. Novoselov, D. Jiang, F. Schedin, T. J. Booth, V. V. Khotkevich, S. V. Morozov, and A. K. Geim. Two-dimensional atomic crystals. Proceedings of the National Academy of Sciences, 102(30):10451{10453, 2005. 7 Andre K. Geim. Nobel lecture: Random walk to graphene. Rev. Mod. Phys., 83:851{862, Aug 2011. 8 K. S. Novoselov. Nobel lecture: Graphene: Materials in the atland. Rev. Mod. Phys., 83:837{849, Aug 2011. 9 H. Kamerlingh Onnes. The superconductivity of mercury. Comm. Phys. Lab. Univ. Leiden, 119, 1911. 10 W. Meissner and R. Ochsenfeld. Ein neuer eekt bei eintritt der supraleitf ahigkeit. Die Naturwissenschaften, 21(44):787{788, Nov 1933. 11 Michael Tinkham. Introduction to Superconductivity: Second Edition (Dover Books on Physics) (Vol i). Dover Publications, 2004. 12 James S. Schilling. Studies in superconductivity at extreme pressures. Physica C: Superconductivity and its Applications, 460-462:182 { 185, 2007. Proceedings of the 8th 134 BIBLIOGRAPHY BIBLIOGRAPHY International Conference on Materials and Mechanisms of Superconductivity and High Temperature Superconductors. 13 A. P. Drozdov, M. I. Eremets, I. A. Troyan, V. Ksenofontov, and S. I. Shylin. Conven- tional superconductivity at 203 kelvin at high pressures in the sulfur hydride system. Nature, 525(7567):73{76, Aug 2015. 14 Dev Kumar Thapa and Anshu Pandey. Evidence for superconductivity at ambient temperature and pressure in nanostructures, 2018. 15 J. Bardeen, L. N. Cooper, and J. R. Schrieer. Theory of superconductivity. Phys. Rev., 108:1175{1204, Dec 1957. 16 N. F. Berk and J. R. Schrieer. Eect of ferromagnetic spin correlations on supercon- ductivity. Phys. Rev. Lett., 17:433{435, Aug 1966. 17 V. L. Ginzburg. On surface superconductivity. Physics Letters, 13:101{102, November 1964. 18 W. A. Little. Possibility of synthesizing an organic superconductor. Phys. Rev., 134:A1416{A1424, Jun 1964. 19 R. H. Ritchie. Plasma Losses by Fast Electrons in Thin Films. Physical Review, 106(5):874{881, June 1957. 20 Nisha Bhukal, Priya, and R. K. Moudgil. Dispersion of two-dimensional plasmons in GaAs quantum well and Ag monolayer. Physica E: Low-dimensional Systems and Nanostructures, 69:13{18, May 2015. 21 Yu Liu, R. F. Willis, K. V. Emtsev, and Th. Seyller. Plasmon dispersion and damping in electrically isolated two-dimensional charge sheets. Phys. Rev. B, 78(20):201403, November 2008. 22 Long Ju, Baisong Geng, Jason Horng, Caglar Girit, Michael Martin, Zhao Hao, Hans A. Bechtel, Xiaogan Liang, Alex Zettl, Y. Ron Shen, and Feng Wang. Graphene plasmonics for tunable terahertz metamaterials. Nature Nanotechnology, 6(10):630{634, October 2011. 23 Frank H. L. Koppens, Darrick E. Chang, and F. Javier Garc a de Abajo. Graphene Plas- monics: A Platform for Strong Light{Matter Interactions. Nano Letters, 11(8):3370{ 3377, August 2011. 24 B. Wunsch, T. Stauber, F. Sols, and F. Guinea. Dynamical polarization of graphene at nite doping. New Journal of Physics, 8(12):318, 2006. 135 BIBLIOGRAPHY BIBLIOGRAPHY 25 E. H. Hwang and S. Das Sarma. Dielectric function, screening, and plasmons in two- dimensional graphene. Phys. Rev. B, 75(20):205418, May 2007. 26 S. Gangadharaiah, A. M. Farid, and E. G. Mishchenko. Charge Response Function and a Novel Plasmon Mode in Graphene. Physical Review Letters, 100(16):166802, April 2008. 27 A. N. Grigorenko, M. Polini, and K. S. Novoselov. Graphene plasmonics. Nature Photonics, 6(11):749{758, November 2012. 28 Tobias Stauber. Plasmonics in Dirac systems: from graphene to topological insulators. Journal of Physics: Condensed Matter, 26(12):123201, 2014. 29 T. Tudorovskiy and S. A. Mikhailov. Intervalley plasmons in graphene. Physical Review B, 82(7):073411, August 2010. 30 H. Rietschel and L. J. Sham. Role of electron Coulomb interaction in superconductivity. Physical Review B, 28(9):5100{5108, November 1983. 31 A. Bill, H. Morawitz, and V. Z. Kresin. Electronic collective modes and superconduc- tivity in layered conductors. Phys. Rev. B, 68(14):144519, October 2003. 32 Ryosuke Akashi and Ryotaro Arita. Development of Density-Functional Theory for a Plasmon-Assisted Superconducting State: Application to Lithium Under High Pres- sures. Phys. Rev. Lett., 111(5):057006, August 2013. 33 A. Linscheid, A. Sanna, and E. K. U. Gross. Ab initio calculation of a Pb single layer on a Si substrate: two-dimensionality and superconductivity. arXiv:1503.00977 [cond-mat], March 2015. arXiv: 1503.00977. 34 Yasutami Takada. Plasmon Mechanism of Superconductivity in Two- and Three- Dimensional Electron Systems. J. Phys. Soc. Jpn., 45(3):786{794, September 1978. 35 Yasutami Takada. Plasmon Mechanism of Superconductivity in the Multivalley Electron Gas. J. Phys. Soc. Jpn., 61(1):238{253, January 1992. 36 Intel. 14nm technology. https://www.intel.com/content/www/us/en/ silicon-innovations/intel-14nm-technology.html, 2018. Accessed: 2018-06-02. 37 Charles Kittel. Introduction to Solid State Physics. John Wiley & Sons, 1986. 38 Yanqing Wu, Yu ming Lin, Ageeth A. Bol, Keith A. Jenkins, Fengnian Xia, Damon B. Farmer, Yu Zhu, and Phaedon Avouris. High-frequency, scaled graphene transistors on diamond-like carbon. Nature, 472(7341):74{78, Apr 2011. 136 BIBLIOGRAPHY BIBLIOGRAPHY 39 B. Radisavljevic, A. Radenovic, J. Brivio, V. Giacometti, and A. Kis. Single-layer MoS 2 transistors. Nature Nanotechnology, 6:147 EP {, Jan 2011. 40 Qiaoliang Bao and Kian Ping Loh. Graphene Photonics, Plasmonics, and Broadband Optoelectronic Devices. ACS Nano, 6(5):3677{3694, May 2012. 41 F. Javier Garc a de Abajo. Graphene Plasmonics: Challenges and Opportunities. ACS Photonics, 1(3):135{152, March 2014. 42 Helmuth Berger, Jie Shan, Kin Fai Mak, L aszl o Forr o, Liang Zhao, Xiaoxiang Xi, and Zefang Wang. Strongly enhanced charge-density-wave order in monolayer NbSe 2 . Nature Nanotechnology, 10(9):765, September 2015. 43 Donglai Feng, Fangyuan Yang, Liguo Ma, Sang-Wook Cheong, Sejoong Kim, Shiyan Li, Xian Hui Chen, Xiaohai Niu, Xiu Fang Lu, Ya Jun Yan, Yijun Yu, Yong-Heum Cho, Young-Woo Son, and Yuanbo Zhang. Gate-tunable phase transitions in thin akes of 1t-TaS 2 . Nature Nanotechnology, 10(3):270, March 2015. 44 Aaron J. Bradley, Alex Zettl, Alexander Riss, Claudia Ojeda-Aristizabal, Dunghai Lee, Hsin-Zon Tsai, Hyejin Ryu, Mark T. Edmonds, Michael F. Crommie, Miguel M. Ugeda, Seita Onishi, Sung-Kwan Mo, Wei Ruan, Yi Chen, Yi Zhang, Zahid Hussain, and Zhi- Xun Shen. Characterization of collective ground states in single-layer NbSe 2 . Nature Physics, 12(1):92, January 2016. 45 Alberto F. Morpurgo, Davide Costanzo, Helmuth Berger, and Sanghyun Jo. Gate- induced superconductivity in atomically thin MoS 2 crystals. Nature Nanotechnology, 11(4):339, April 2016. 46 Kin Fai Mak, Changgu Lee, James Hone, Jie Shan, and Tony F. Heinz. Atomically Thin MoS 2 : A New Direct-Gap Semiconductor. Phys. Rev. Lett., 105(13):136805, September 2010. 47 Aaron J. Bradley, Diana Y. Qiu, Felipe H. da Jornada, Feng Wang, Michael F. Crommie, Miguel M. Ugeda, Steven G. Louie, Su-Fei Shi, Sung-Kwan Mo, Wei Ruan, Yi Zhang, Zahid Hussain, and Zhi-Xun Shen. Giant bandgap renormalization and excitonic ef- fects in a monolayer transition metal dichalcogenide semiconductor. Nature Materials, 13(12):1091, December 2014. 48 Yuxuan Lin, Xi Ling, Lili Yu, Shengxi Huang, Allen L. Hsu, Yi-Hsien Lee, Jing Kong, Mildred S. Dresselhaus, and Tom as Palacios. Dielectric Screening of Excitons and Trions in Single-Layer MoS2. Nano Lett., 14(10):5569{5576, October 2014. 49 G. Sch onho, M. R osner, R. E. Groenewald, S. Haas, and T. O. Wehling. Interplay of screening and superconductivity in low-dimensional materials. Phys. Rev. B, 94:134504, Oct 2016. 137 BIBLIOGRAPHY BIBLIOGRAPHY 50 Archana Raja, Andrey Chaves, Jaeeun Yu, Ghidewon Arefe, Heather M. Hill, Albert F. Rigosi, Timothy C. Berkelbach, Philipp Nagler, Christian Sch uller, Tobias Korn, Colin Nuckolls, James Hone, Louis E. Brus, Tony F. Heinz, David R. Reichman, and Alexey Chernikov. Coulomb engineering of the bandgap and excitons in two-dimensional ma- terials. Nature Communications, 8:15251, May 2017. 51 A. Steinho, F. Jahnke, G. Sch onho, M. Florian, M. R osner, and T. O. Wehling. Exciton ssion in monolayer transition metal dichalcogenide semiconductors. Nature Communications, 8(1):1166, December 2017. 52 M. A. L. Marques, M. L uders, N. N. Lathiotakis, G. Profeta, A. Floris, L. Fast, A. Con- tinenza, E. K. U. Gross, and S. Massidda. Ab initio. Phys. Rev. B, 72(2):024546, July 2005. 53 E. R. Margine, Feliciano Giustino, and Henry Lambert. Electron-phonon interaction and pairing mechanism in superconducting Ca-intercalated bilayer graphene. Scientic Reports, 6:21414, February 2016. 54 Ryosuke Akashi, Kazuma Nakamura, Ryotaro Arita, and Masatoshi Imada. High- temperature superconductivity in layered nitrides-Li x mncl (m = Ti, Zr, Hf): Insights from density functional theory for superconductors. Phys. Rev. B, 86(5):054513, August 2012. 55 R. Rold an, E. Cappelluti, and F. Guinea. Interactions and superconductivity in heavily doped MoS 2 . Phys. Rev. B, 88(5):054515, August 2013. 56 A. K. Geim and I. V. Grigorieva. Van der Waals heterostructures. Nature, 499(7459):419, July 2013. 57 Giuseppe Grosso and Giuseppe Pastori Parravicini. Solid State Physics. Academic Press, 2000. 58 Qing-Hong Cao. Periodic systems and the bloch theorem. http://www.phy.pku.edu. cn/ ~ qhcao/resources/class/QM/bloch.pdf. Accessed: 2017-02-03. 59 Anthony Paxton. Introduction to the tight binding approximation - implementation by diagonalisation, 3 2009. 60 Nicola Marzari and David Vanderbilt. Maximally localized generalized wannier func- tions for composite energy bands. Phys. Rev. B, 56:12847{12865, Nov 1997. 61 Ivo Souza, Nicola Marzari, and David Vanderbilt. Maximally localized wannier functions for entangled energy bands. Phys. Rev. B, 65:035109, Dec 2001. 138 BIBLIOGRAPHY BIBLIOGRAPHY 62 Arash A. Mosto, Jonathan R. Yates, Young-Su Lee, Ivo Souza, David Vanderbilt, and Nicola Marzari. wannier90: A tool for obtaining maximally-localised Wannier functions. Computer Physics Communications, 178(9):685{699, May 2008. 63 Malte R osner. Electronic Structure of Novel Two-dimensional Materials and Graphene Heterostructures. PhD thesis, Bremen Universitat, 2016. 64 J. C. Slater and G. F. Koster. Simplied lcao method for the periodic potential problem. Phys. Rev., 94:1498{1524, Jun 1954. 65 A P Sutton, M W Finnis, D G Pettifor, and Y Ohta. The tight-binding bond model. Journal of Physics C: Solid State Physics, 21(1):35, 1988. 66 P. Hohenberg and W. Kohn. Inhomogeneous electron gas. Phys. Rev., 136:B864{B871, Nov 1964. 67 W. Kohn and L. J. Sham. Self-consistent equations including exchange and correlation eects. Phys. Rev., 140:A1133{A1138, Nov 1965. 68 Richard M. Martin. Electronic Structure: Basic Theory and Practical Methods (v. 1). Cambridge University Press, 2008. 69 Friedhelm Bechstedt. Many-Body Approach to Electronic Excitations: Concepts and Applications (Springer Series in Solid-State Sciences). Springer, 2014. 70 J. P. Perdew and Alex Zunger. Self-interaction correction to density-functional approx- imations for many-electron systems. Phys. Rev. B, 23:5048{5079, May 1981. 71 A. Sommerfeld and H. Bethe. Elektronentheorie der Metalle (Heidelberger Taschenb ucher) (German Edition). Springer, 1967. 72 Conyers Herring. A new method for calculating wave functions in crystals. Phys. Rev., 57:1169{1177, Jun 1940. 73 E. Wigner and F. Seitz. On the constitution of metallic sodium. Phys. Rev., 43:804{810, May 1933. 74 J. C. Slater. Wave functions in a periodic potential. Phys. Rev., 51:846{851, May 1937. 75 O. Krogh Andersen. Linear methods in band theory. Phys. Rev. B, 12:3060{3083, Oct 1975. 76 E. Wimmer, H. Krakauer, M. Weinert, and A. J. Freeman. Full-potential self-consistent linearized-augmented-plane-wave method for calculating the electronic structure of molecules and surfaces: o 2 molecule. Phys. Rev. B, 24:864{875, Jul 1981. 139 BIBLIOGRAPHY BIBLIOGRAPHY 77 D. R. Hamann, M. Schl uter, and C. Chiang. Norm-conserving pseudopotentials. Phys. Rev. Lett., 43:1494{1497, Nov 1979. 78 Eric L. Shirley, Douglas C. Allan, Richard M. Martin, and J. D. Joannopoulos. Extended norm-conserving pseudopotentials. Phys. Rev. B, 40:3652{3660, Aug 1989. 79 David Vanderbilt. Optimally smooth norm-conserving pseudopotentials. Phys. Rev. B, 32:8412{8415, Dec 1985. 80 Ivo Souza, Tim Wilkens, and Richard M. Martin. Polarization and localization in insulators: Generating function approach. Phys. Rev. B, 62:1666{1683, Jul 2000. 81 Alexander L. Fetter and John Dirk Walecka. Quantum Theory of Many-Particle Systems (Dover Books on Physics). Dover Publications, 2012. 82 R. Heid. Emergent Phenomena in Correlated Matter. Forschungszentrum Julich, 2013. 83 B. M. Santoyo and M Del Castillo-Mussot. Plasmons in three, two and one dimension. Revista Mexicana de Fisica, 39:640, Apr 1993. 84 Enrico Fermi. Nuclear Physics: A Course Given by Enrico Fermi at the University of Chicago. University of Chicago Press, 1974. 85 H. M. Nussenzveig. Causality and Dispersion Relations (Mathematics in Science and Engineering). Academic Press, 1972. 86 M. Born and R. Oppenheimer. Zur quantentheorie der molekeln. Annalen der Physik, 389(20):457{484, 1927. 87 R. P. Feynman. Forces in molecules. Phys. Rev., 56:340{343, Aug 1939. 88 Richard D. Mattuck. A Guide to Feynman Diagrams in the Many-Body Problem: Second Edition (Dover Books on Physics). Dover Publications, 1992. 89 C. Friedrich and A. Schindlmayr. Computing Solids: Models, ab-inito methods and supercomputing, volume 74. J ulich : Forschungszentrum, Zentralbibliothek, 2014. 90 Henrik Bruus and Karsten Flensberg. Many-Body Quantum Theory in Condensed Mat- ter Physics: An Introduction (Oxford Graduate Texts). Oxford University Press, 2004. 91 T. Miyake, C. Martins, R. Sakuma, and F. Aryasetiawan. Eects of momentum- dependent self-energy in the electronic structure of correlated materials. Phys. Rev. B, 87:115110, Mar 2013. 92 Gerald D. Mahan. Many-Particle Physics (Physics of Solids and Liquids). Springer, 2010. 140 BIBLIOGRAPHY BIBLIOGRAPHY 93 D. R. Hartree. The wave mechanics of an atom with a non-coulomb central eld. part iii. term values and intensities in series in optical spectra. Mathematical Proceedings of the Cambridge Philosophical Society, 24(3):426{437, 1928. 94 E. Wigner. On the interaction of electrons in metals. Phys. Rev., 46:1002{1011, Dec 1934. 95 Xinguo Ren, Patrick Rinke, Christian Joas, and Matthias Scheer. Random-phase approximation and its applications in computational chemistry and materials science. Journal of Materials Science, 47(21):7447{7471, June 2012. 96 Wolfgang von der Linden and Peter Horsch. Precise quasiparticle energies and hartree- fock bands of semiconductors and insulators. Phys. Rev. B, 37:8351{8362, May 1988. 97 G. E. Engel and Behnam Farid. Generalized plasmon-pole model and plasmon band structures of crystals. Phys. Rev. B, 47:15931{15934, Jun 1993. 98 Paul Larson, Marc Dvorak, and Zhigang Wu. Role of the plasmon-pole model in the gw approximation. Phys. Rev. B, 88:125205, Sep 2013. 99 Ryosuke Akashi and Ryotaro Arita. Density functional theory for plasmon-assisted superconductivity. Journal of the Physical Society of Japan, 83(6):061016, 2014. 100 Fabio Caruso and Feliciano Giustino. Theory of electron-plasmon coupling in semicon- ductors. Phys. Rev. B, 94:115208, Sep 2016. 101 C. N. Berglund and W. E. Spicer. Photoemission studies of copper and silver: Theory. Phys. Rev., 136:A1030{A1044, Nov 1964. 102 M. L uders, M. A. L. Marques, N. N. Lathiotakis, A. Floris, G. Profeta, L. Fast, A. Con- tinenza, S. Massidda, and E. K. U. Gross. Ab initio theory of superconductivity. i. density functional formalism and approximate functionals. Phys. Rev. B, 72:024545, Jul 2005. 103 Miguel Marques and EKU Gross. Density functional theory for superconductors: ex- change and correlation potentials for inhomogeneous systems. PhD thesis, Bayerischen Julius-Maximilians-Universitat, Wurzburg, 2000. 104 Leon N. Cooper. Bound electron pairs in a degenerate fermi gas. Phys. Rev., 104:1189{ 1190, Nov 1956. 105 N. N. Bogoljubov. On a new method in the theory of superconductivity. Il Nuovo Cimento (1955-1965), 7(6):794{805, Mar 1958. 141 BIBLIOGRAPHY BIBLIOGRAPHY 106 Giovanni AC Ummarino. 13 eliashberg theory. Emergent Phenomena in Correlated Matter: Autumn School Organized by the Forschungszentrum J ulich and the German Research School for Simulation Sciences at Forschungszentrum J ulich 23-27 September 2013; Lecture Notes of the Autumn School Correlated Electrons 2013, 3, 2013. 107 L. N. Oliveira, E. K. U. Gross, and W. Kohn. Density-functional theory for supercon- ductors. Phys. Rev. Lett., 60:2430{2433, Jun 1988. 108 N. David Mermin. Thermal properties of the inhomogeneous electron gas. Phys. Rev., 137:A1441{A1443, Mar 1965. 109 E. K. U. Gross and Stefan Kurth. Density functional theory of the superconducting state. International Journal of Quantum Chemistry, 40(S25):289{297, 1991. 110 Andreas G orling and Mel Levy. Exact kohn-sham scheme based on perturbation theory. Phys. Rev. A, 50:196{204, Jul 1994. 111 Ryosuke Akashi and Ryotaro Arita. Density functional theory for superconductors with particle-hole asymmetric electronic structure. Phys. Rev. B, 88:014514, Jul 2013. 112 W. L. McMillan. Transition temperature of strong-coupled superconductors. Phys. Rev., 167:331{344, Mar 1968. 113 P. Morel and P. W. Anderson. Calculation of the superconducting state parameters with retarded electron-phonon interaction. Phys. Rev., 125:1263{1271, Feb 1962. 114 S.P. Tewari and Paramjeet Kaur Gumber. Eect of plasmons on the critical temper- ature of lanthanum and yttrium based high temperature superconductors. Physica C: Superconductivity, 171(1):147 { 150, 1990. 115 N. N. Lathiotakis, M. A. L. Marques, M. L uders, L. Fast, and E. K. U. Gross. Density functional theory for superconductors. International Journal of Quantum Chemistry, 99(5):790{797, 2004. 116 M Roesner, RE Groenewald, G Schoenho, J Berges, S Haas, and TO Wehling. Plas- monic superconductivity in layered materials. arXiv preprint arXiv:1803.04576, 2018. 117 F. Aryasetiawan, M. Imada, A. Georges, G. Kotliar, S. Biermann, and A. I. Lichtenstein. Frequency-dependent local interactions and low-energy eective models from electronic structure calculations. Phys. Rev. B, 70(19):195104, November 2004. 118 L. V. Keldysh. Coulomb interaction in thin semiconductor and semimetal lms. JETP Letters, 29:716, 1979. 142 BIBLIOGRAPHY BIBLIOGRAPHY 119 M. R osner, E. S a so glu, C. Friedrich, S. Bl ugel, and T. O. Wehling. Wannier function approach to realistic Coulomb interactions in layered materials and heterostructures. Phys. Rev. B, 92(8):085102, August 2015. 120 T. Ritschel, J. Trinckauf, K. Koepernik, B. B uchner, M. v Zimmermann, H. Berger, Y. I. Joe, P. Abbamonte, and J. Geck. Orbital textures and charge density waves in transition metal dichalcogenides. Nature Physics, 11(4):328{331, April 2015. 121 K. Rossnagel. On the origin of charge-density waves in select layered transition-metal dichalcogenides. J. Phys.: Condens. Matter, 23(21):213001, 2011. 122 Erik G. C. P. van Loon, Malte R osner, Gunnar Sch onho, Mikhail I. Katsnelson, and Tim O. Wehling. Competing Coulomb and electron-phonon interactions in NbS 2 . NPJ Quantum Materials, 3(1):32, 2018. 123 Thomas Ayral, Philipp Werner, and Silke Biermann. Spectral Properties of Correlated Materials: Local Vertex and Nonlocal Two-Particle Correlations from Combined gw and Dynamical Mean Field Theory. Phys. Rev. Lett., 109(22):226401, November 2012. 124 Li Huang, Thomas Ayral, Silke Biermann, and Philipp Werner. Extended dynamical mean-eld study of the Hubbard model with long-range interactions. Phys. Rev. B, 90(19):195114, November 2014. 125 A. Einstein. Die plancksche theorie der strahlung und die theorie der spezischen w arme. Annalen der Physik, 327(1):180{190, 1907. 126 R. E. Groenewald, M. R osner, G. Sch onho, S. Haas, and T. O. Wehling. Valley plasmonics in transition metal dichalcogenides. Phys. Rev. B, 93(20):205145, May 2016. 127 Kirsten Andersen, Simone Latini, and Kristian S. Thygesen. Dielectric Genome of van der Waals Heterostructures. Nano Lett., 15(7):4616{4621, July 2015. 128 Yoichiro Nambu. Quasi-Particles and Gauge Invariance in the Theory of Superconduc- tivity. Phys. Rev., 117(3):648{663, February 1960. 129 Fabio Caruso and Feliciano Giustino. Spectral ngerprints of electron-plasmon coupling. Phys. Rev. B, 92(4):045123, July 2015. 130 Philipp Werner and Michele Casula. Dynamical screening in correlated electron systems|from lattice models to realistic materials. J. Phys.: Condens. Matter, 28(38):383001, 2016. 131 Philip B. Allen and Bo zidar Mitrovi c. Theory of Superconducting Tc. In Henry Ehren- reich, Frederick Seitz, and David Turnbull, editors, Solid State Physics, volume 37, pages 1{92. Academic Press, January 1983. 143 BIBLIOGRAPHY BIBLIOGRAPHY 132 Christoph Heil, Samuel Ponc e, Henry Lambert, Martin Schlipf, Elena R. Margine, and Feliciano Giustino. Origin of Superconductivity and Latent Charge Density Wave in NbS 2 . Phys. Rev. Lett., 119(8):087003, August 2017. 133 Tanmoy Das and Kapildeb Dolui. Superconducting dome in MoS 2 and TiSe 2 generated by quasiparticle-phonon coupling. Phys. Rev. B, 91(9):094510, March 2015. 134 Morten N. Gjerding, Mohnish Pandey, and Kristian S. Thygesen. Band structure engi- neered layered metals for low-loss plasmonics. Nature Communications, 8:ncomms15133, April 2017. 135 Yizhi Ge and Amy Y. Liu. Phonon-mediated superconductivity in electron-doped single- layer MoS 2 : A rst-principles prediction. Phys. Rev. B, 87:241408, Jun 2013. 136 M. R osner, S. Haas, and T. O. Wehling. Phase diagram of electron-doped dichalco- genides. Phys. Rev. B, 90(24):245105, December 2014. 137 H. Suderow, V. G. Tissen, J. P. Brison, J. L. Mart nez, and S. Vieira. Pressure induced eects on the fermi surface of superconducting 2H-NbSe 2 . Phys. Rev. Lett., 95:117006, Sep 2005. 138 Kirsten Andersen and Kristian S. Thygesen. Plasmons in metallic monolayer and bilayer transition metal dichalcogenides. Phys. Rev. B, 88(15):155128, October 2013. 139 Pere Mir o, Martha Audired, and Thomas Heine. An atlas of two-dimensional materials. Chem. Soc. Rev., 43(18):6537{6554, August 2014. 140 Franti sek Karlick y, Kasibhatta Kumara Ramanatha Datta, Michal Otyepka, and Radek Zbo ril. Halogenated Graphenes: Rapidly Growing Family of Graphene Derivatives. ACS Nano, 7(8):6434{6464, August 2013. 141 Jasper van Wezel, Roman Schuster, Andreas K onig, Martin Knupfer, Jeroen van den Brink, Helmuth Berger, and Bernd B uchner. Eect of Charge Order on the Plasmon Dis- persion in Transition-Metal Dichalcogenides. Physical Review Letters, 107(17):176404, October 2011. 142 A. K onig, K. Koepernik, R. Schuster, R. Kraus, M. Knupfer, B. B uchner, and H. Berger. Plasmon evolution and charge-density wave suppression in potassium intercalated 2H- TaSe 2 . EPL (Europhysics Letters), 100(2):27002, 2012. 143 A. K onig, R. Schuster, M. Knupfer, B. B uchner, and H. Berger. Doping dependence of the plasmon dispersion in 2H-TaSe 2 . Physical Review B, 87(19):195119, May 2013. 144 Kourosh Kalantar-zadeh, Jian Zhen Ou, Torben Daeneke, Michael S. Strano, Martin Pumera, and Sally L. Gras. Two-Dimensional Transition Metal Dichalcogenides in Biosystems. Advanced Functional Materials, 25(32):5086{5099, August 2015. 144 BIBLIOGRAPHY BIBLIOGRAPHY 145 Kourosh Kalantar-zadeh and Jian Zhen Ou. Biosensors Based on Two-Dimensional MoS2. ACS Sensors, November 2015. 146 J. B. Maurya, Y. K. Prajapati, V. Singh, J. P. Saini, and Rajeev Tripathi. Performance of graphene{MoS2 based surface plasmon resonance sensor using Silicon layer. Optical and Quantum Electronics, 47(11):3599{3611, July 2015. 147 A. Steinho, M. R osner, F. Jahnke, T. O. Wehling, and C. Gies. In uence of Excited Carriers on the Optical and Electronic Properties of MoS2. Nano Letters, 14(7):3743{ 3748, July 2014. 148 R. S. Title and M. W. Shafer. Electron-paramagnetic-resonance studies on arsenic acceptors in natural (2H) and synthetic (3R) MoS 2 crystals. Phys. Rev. B, 8:615{620, Jul 1973. 149 N. Wakabayashi, H. G. Smith, and R. M. Nicklow. Lattice dynamics of hexagonal mos 2 studied by neutron scattering. Phys. Rev. B, 12:659{663, Jul 1975. 150 Davide Costanzo, Sanghyun Jo, Helmuth Berger, and Alberto F. Morpurgo. Gate- induced superconductivity in atomically thin MoS2 crystals. Nature Nanotechnology, 11(4):339{344, January 2016. 151 Wang Qing-Yan, Li Zhi, Zhang Wen-Hao, Zhang Zuo-Cheng, Zhang Jin-Song, Li Wei, Ding Hao, Ou Yun-Bo, Deng Peng, Chang Kai, Wen Jing, Song Can-Li, He Ke, Jia Jin- Feng, Ji Shuai-Hua, Wang Ya-Yu, Wang Li-Li, Chen Xi, Ma Xu-Cun, and Xue Qi-Kun. Interface-induced high-temperature superconductivity in single unit-cell FeSe lms on SrTiO 3 . Chinese Physics Letters, 29(3):037402, 2012. 152 Sinisa Coh, Marvin L Cohen, and Steven G Louie. Large electron{phonon interactions from FeSe phonons in a monolayer. New Journal of Physics, 17(7):073027, 2015. 153 Leiqiang Chu, Hennrik Schmidt, Jiang Pu, Shunfeng Wang, Barbaros Zyilmaz, Taishi Takenobu, and Goki Eda. Charge transport in ion-gated mono-, bi-, and trilayer MoS 2 eld eect transistors. Scientic Reports, 4:7293 EP {, Dec 2014. Article. 154 Harry Fisher. First-principles investigation of electron-phonon interactions in novel superconductors. PhD thesis, Oxford University, UK, 2014. 155 Paolo Giannozzi, Stefano Baroni, Nicola Bonini, Matteo Calandra, Roberto Car, Carlo Cavazzoni, Davide Ceresoli, Guido L Chiarotti, Matteo Cococcioni, Ismaila Dabo, An- drea Dal Corso, Stefano de Gironcoli, Stefano Fabris, Guido Fratesi, Ralph Gebauer, Uwe Gerstmann, Christos Gougoussis, Anton Kokalj, Michele Lazzeri, Layla Martin- Samos, Nicola Marzari, Francesco Mauri, Riccardo Mazzarello, Stefano Paolini, Alfredo Pasquarello, Lorenzo Paulatto, Carlo Sbraccia, Sandro Scandolo, Gabriele Sclauzero, 145 BIBLIOGRAPHY BIBLIOGRAPHY Ari P Seitsonen, Alexander Smogunov, Paolo Umari, and Renata M Wentzcovitch. Quantum espresso: a modular and open-source software project for quantum simula- tions of materials. Journal of Physics: Condensed Matter, 21(39):395502 (19pp), 2009. 156 P Giannozzi, O Andreussi, T Brumme, O Bunau, M Buongiorno Nardelli, M Calandra, R Car, C Cavazzoni, D Ceresoli, M Cococcioni, N Colonna, I Carnimeo, A Dal Corso, S de Gironcoli, P Delugas, R A DiStasio Jr, A Ferretti, A Floris, G Fratesi, G Fugallo, R Gebauer, U Gerstmann, F Giustino, T Gorni, J Jia, M Kawamura, H-Y Ko, A Kokalj, E K u c ukbenli, M Lazzeri, M Marsili, N Marzari, F Mauri, N L Nguyen, H-V Nguyen, A Otero de-la Roza, L Paulatto, S Ponc e, D Rocca, R Sabatini, B Santra, M Schlipf, A P Seitsonen, A Smogunov, I Timrov, T Thonhauser, P Umari, N Vast, X Wu, and S Baroni. Advanced capabilities for materials modelling with quantum espresso. Journal of Physics: Condensed Matter, 29(46):465901, 2017. 157 Youngki Yoon, Kartik Ganapathi, and Sayeef Salahuddin. How good can monolayer MoS 2 transistors be? Nano Letters, 11(9):3768{3773, 2011. PMID: 21790188. 158 Ali Hussain Reshak and S. Auluck. Calculated optical properties of 2H-MoS 2 interca- lated with lithium. Phys. Rev. B, 68:125101, Sep 2003. 159 Christoph Friedrich, Stefan Bl ugel, and Arno Schindlmayr. Ecient implementation of the GW approximation within the all-electron FLAPW method. Physical Review B, 81(12):125102, March 2010. 160 The Juelich FLEUR project, 2014. 161 Z. Y. Zhu, Y. C. Cheng, and U. Schwingenschl ogl. Giant spin-orbit-induced spin split- ting in two-dimensional transition-metal dichalcogenide semiconductors. Physical Re- view B, 84(15):153402, October 2011. 162 G. Kresse and J. Hafner. Ab initio molecular dynamics for liquid metals. Physical Review B, 47(1):558{561, January 1993. 163 G. Kresse and J. Furthm uller. Eciency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set. Computational Materials Science, 6(1):15{50, July 1996. 164 Andreas Scholz, Tobias Stauber, and John Schliemann. Plasmons and screening in a monolayer of MoS 2 . Physical Review B, 88(3):035135, July 2013. 165 K. Kechedzhi and D. S. L. Abergel. Weakly damped acoustic plasmon mode in transition metal dichalcogenides with Zeeman splitting. Physical Review B, 89(23):235420, June 2014. 146 BIBLIOGRAPHY BIBLIOGRAPHY 166 Yufeng Liang and Li Yang. Carrier Plasmon Induced Nonlinear Band Gap Renormal- ization in Two-Dimensional Semiconductors. Physical Review Letters, 114(6):063001, February 2015. 167 Gui-Bin Liu, Wen-Yu Shan, Yugui Yao, Wang Yao, and Di Xiao. Three-band tight- binding model for monolayers of group-VIB transition metal dichalcogenides. Physical Review B, 88(8):085433, August 2013. 168 Friedrich Roth, Andreas K onig, J org Fink, Bernd B uchner, and Martin Knupfer. Elec- tron energy-loss spectroscopy: A versatile tool for the investigations of plasmonic exci- tations. Journal of Electron Spectroscopy and Related Phenomena, 195:85{95, August 2014. 169 R. Rold an, E. Cappelluti, and F. Guinea. Interactions and superconductivity in heavily doped mos 2 . Phys. Rev. B, 88:054515, Aug 2013. 170 G. M. Eliashberg. Interactions between electrons and lattice vibrations in a supercon- ductor. Soviet Phys. JETP, 11, 1960. 171 P. B. Allen and R. C. Dynes. Transition temperature of strong-coupled superconductors reanalyzed. Phys. Rev. B, 12(3):905{922, August 1975. 172 R. F. Frindt. Superconductivity in ultrathin NbSe 2 layers. Phys. Rev. Lett., 28:299{301, Jan 1972. 173 Xiaoxiang Xi, Liang Zhao, Zefang Wang, Helmuth Berger, L aszl o Forr o, Jie Shan, and Kin Fai Mak. Strongly enhanced charge-density-wave order in monolayer NbSe2. Nature Nanotechnology, 10(9):765{769, July 2015. 174 Yijun Yu, Fangyuan Yang, Xiu Fang Lu, Ya Jun Yan, Yong-Heum Cho, Liguo Ma, Xi- aohai Niu, Sejoong Kim, Young-Woo Son, Donglai Feng, Shiyan Li, Sang-Wook Cheong, Xian Hui Chen, and Yuanbo Zhang. Gate-tunable phase transitions in thin akes of 1t-TaS2. Nature Nanotechnology, 10(3):270{276, January 2015. 175 Y. Cao, A. Mishchenko, G. L. Yu, E. Khestanova, A. P. Rooney, E. Prestat, A. V. Kretinin, P. Blake, M. B. Shalom, C. Woods, J. Chapman, G. Balakrishnan, I. V. Grigorieva, K. S. Novoselov, B. A. Piot, M. Potemski, K. Watanabe, T. Taniguchi, S. J. Haigh, A. K. Geim, and R. V. Gorbachev. Quality heterostructures from two- dimensional crystals unstable in air by their assembly in inert atmosphere. Nano Letters, 15(8):4914{4921, July 2015. 176 G abor Zsolt Magda, J anos Pet} o, Gergely Dobrik, Chanyong Hwang, L aszl o P. Bir o, and Levente Tapaszt o. Exfoliation of large-area transition metal chalcogenide single layers. Scientic Reports, 5(1), October 2015. 147 BIBLIOGRAPHY BIBLIOGRAPHY 177 Suljo Linic, Phillip Christopher, and David B. Ingram. Plasmonic-metal nanostructures for ecient conversion of solar to chemical energy. Nature Materials, 10(12):911{921, December 2011. 178 P. Muhlschlegel. Resonant optical antennas. Science, 308(5728):1607{1609, Jun 2005. 179 Ram on Alvarez-Puebla, Luis M. Liz-Marz an, and F. Javier Garc a de Abajo. Light concentration at the nanometer scale. The Journal of Physical Chemistry Letters, 1(16):2428{2434, July 2010. 180 D. Ansell, I. P. Radko, Z. Han, F. J. Rodriguez, S. I. Bozhevolnyi, and A. N. Grig- orenko. Hybrid graphene plasmonic waveguide modulators. Nature Communications, 6(1), November 2015. 181 Jibin Song, Peng Huang, Hongwei Duan, and Xiaoyuan Chen. Plasmonic vesicles of amphiphilic nanocrystals: Optically active multifunctional platform for cancer diagnosis and therapy. Accounts of Chemical Research, 48(9):2506{2515, July 2015. 182 Harry A. Atwater and Albert Polman. Plasmonics for improved photovoltaic devices. Nature Materials, 9(3):205{213, February 2010. 183 Boxiang Song, Yuhan Yao, Roelof E. Groenewald, Yunxiang Wang, He Liu, Yifei Wang, Yuanrui Li, Fanxin Liu, Stephen B. Cronin, Adam M. Schwartzberg, Stefano Cabrini, Stephan Haas, and Wei Wu. Probing Gap Plasmons Down to Subnanometer Scales Using Collapsible Nanongers. ACS nano, 11(6):5836{5843, June 2017. 184 Jon A. Schuller, Edward S. Barnard, Wenshan Cai, Young Chul Jun, Justin S. White, and Mark L. Brongersma. Plasmonics for extreme light concentration and manipulation. Nature Materials, 9(3):193{204, February 2010. 185 Naomi J. Halas, Surbhi Lal, Wei-Shun Chang, Stephan Link, and Peter Nordlander. Plasmons in strongly coupled metallic nanostructures. Chemical Reviews, 111(6):3913{ 3961, June 2011. 186 Serghey V. Vonsovsky and Mikhail I. Katsnelson. Quantum Solid-State Physics (Springer Series in Solid-State Sciences). Springer, 2012. 187 Tom Westerhout, Edo van Veen, Mikhail I. Katsnelson, and Shengjun Yuan. Plasmon connement in fractal quantum systems. Phys. Rev. B, 97:205434, May 2018. 188 Weihua Wang, Thomas Christensen, Antti-Pekka Jauho, Kristian S. Thygesen, Mar- tijn Wubs, and N. Asger Mortensen. Plasmonic eigenmodes in individual and bow-tie graphene nanotriangles. Scientic Reports, 5(1), April 2015. 148 BIBLIOGRAPHY BIBLIOGRAPHY 189 M. R osner, C. Steinke, M. Lorke, C. Gies, F. Jahnke, and T. O. Wehling. Two- dimensional heterojunctions from nonlocal manipulations of the interactions. Nano Letters, 16(4):2322{2327, March 2016. 190 Timothy J. Davis, Daniel E. G omez, and Ann Roberts. Plasmonic circuits for manipu- lating optical information. Nanophotonics, 6(3), January 2016. 191 Enrique Maci a. The role of aperiodic order in science and technology. Reports on Progress in Physics, 69(2):397{441, December 2005. 192 Priya Johari and Vivek B. Shenoy. Tunable Dielectric Properties of Transition Metal Dichalcogenides. ACS Nano, 5(7):5903{5908, July 2011. 193 M. N. Faraggi, A. Arnau, and V. M. Silkin. Role of band structure and local-eld eects in the low-energy collective electronic excitation spectra of 2H-NbSe 2 . Physical Review B, 86(3):035115, July 2012. 194 Pierluigi Cudazzo, Matteo Gatti, and Angel Rubio. Plasmon dispersion in layered transition-metal dichalcogenides. Physical Review B, 86(7):075121, August 2012. 195 Pierluigi Cudazzo, Matteo Gatti, and Angel Rubio. Interplay between structure and electronic properties of layered transition-metal dichalcogenides: Comparing the loss function of 1t and 2h polymorphs. Physical Review B, 90(20):205128, November 2014. 149 Acknowledgments I would like to thank a number of people without who this thesis would not have been possible. Firstly, I am grateful for my supervisor Dr. Stephan Haas, for his guidance throughout the PhD process. I am also very grateful for the opportunity he aorded me to work in Germany every summer from 2014 to 2017. The connections I formed there shaped my PhD in entirety. Therefore, I am extremely grateful for the eort he put in to help me establish various co-operations. On that note, I want to thank all the collaborators I worked with over the years (and will surely work with in the future). At the top of this list is Dr. Tim Wehling and the people from his research group, Dr. Malte R osner, Dr. Gunnar Sch onho and Jan Berges. To Dr. Malte R osner I would like to say a very special, thank you, for all the work we did together. The foundation for this thesis lies in the work he did during his PhD and the body of it is based on the consequent work we did together. It is no exaggeration to say this thesis would not have been if not for the close collaboration I had with him. I'd also like to thank my local collaborators, Dr. Wei Wu and his student Boxiang Song for joined worked we did on nano-particle dimers. Secondly, I would like to thank my committee, Dr. Haas, Dr. Di Felice, Dr. Kresin, Dr. Nakano and Dr. Zanardi, for their time and eort reviewing and refereeing this thesis as well as administering my qualifying exam and thesis defense. Especially, thank you to Dr. Aiichiro Nakano for also guiding me throughout the process of getting my CS masters degree. The skills I learned there undoubtedly aided in this thesis. I'd like to thank the USC Physics and Astronomy department for continuous nancial support throughout my PhD. It goes without saying that I would not have been able to complete this thesis without that crucial component. Special gratitude goes to the ladies in the administrative roles, Betty Byers, Lisa Moeller and Christina Tasulis Williams, for always helping out with 'oce related' problems. Also, thank you to Joe Vandiver for being very helpful in scheduling labs in a way to accommodate my research as much as possible. 150 BIBLIOGRAPHY BIBLIOGRAPHY Then I would also like to thank my LA friends who made the time in between this thesis' work so enjoyable that I could continue the work throughout. Thanks to them, LA has become home (and many of them literally provided me shelter during my, too many, homeless stretches). By name I would like to mention Andrew, Artemis, Emily, Jared, Je, Johann, Jolene, Levon, Lisa, Malte, Michelle, Pancham and Rebecca. A special thanks also to Reni for all her encouragement and support. Last but denitely not least, I thank my parents and sister for their continuous and unconditional love and support in life. 151
Abstract (if available)
Abstract
Screening of the Coulomb interaction is generally suppressed in low dimensional materials and heavily affected by substrates. This in turn affects the interaction between electrons, which plays a major role in many-body effects. In the field of superconductivity, the Coulomb interaction is usually reduced to an approximate static parameter μ (retarded Morel-Anderson Coulomb potential), which describes only the effective Coulomb repulsion, and therefore is responsible for reducing the superconducting transition temperature. ❧ Our calculations of μ* reveal that the Coulomb interactions, renormalized by the reduced layer thickness and the substrate properties, can shift the onset of the electron-phonon driven superconducting phase in monolayer MoS2, but do not significantly affect the critical temperature at optimal doping. ❧ We overcome the inadequate handling of the Coulomb interaction as a static parameter and present an ab initio based dynamical Coulomb description. We apply this method to MoS₂ to account, in a realistic manner, for intrinsic material screening, substrate screening as well as dynamical screening processes. Using this method we can reliably describe the plasmonic excitations and their coupling to the electrons. ❧ Plasmonic excitations, brought on by the dynamical screening of the Coulomb interaction, behave fundamentally different in layered materials in comparison to bulk systems. They form gapless modes, which in turn couple at low energies to the electrons. Thereby they can strongly influence superconducting instabilities. ❧ We show how these excitations can be controlled from the outside via changes in the dielectric environment or in the doping level, which allows for external tuning of the superconducting transition temperature. By solving the gap equation for an effective system, we find, generally, that the plasmonic influence can both strongly enhance or reduce the transition temperature, depending on the details of the plasmon-phonon interplay. We formulate simple experimental guidelines to find plasmon-induced elevated transition temperatures in layered materials. ❧ Lastly, we provide details of how this approach can be used to explore the interplay of plasmons and phonons in real materials (such as MoS₂) and find an optimal ratio between substrate screening and electron doping which maximizes the superconducting critical temperature of the resulting induced phononic and plasmonic superconducting state.
Linked assets
University of Southern California Dissertations and Theses
Conceptually similar
PDF
Plasmonic excitations in quantum materials: topological insulators and metallic monolayers on dielectric substrates
PDF
Superconductivity in low-dimensional topological systems
PDF
Signatures of topology in a quasicrystal: a case study of the non-interacting and superconducting Fibonacci chain
PDF
Light Emission from Carbon Nanotubes and Two-Dimensional Materials
PDF
Optoelectronic, thermoelectric, and photocatalytic properties of low dimensional materials
PDF
Investigation of the superconducting proximity effect (SPE) and magnetc dead layers (MDL) in thin film double layers
PDF
Designing data-effective machine learning pipeline in application to physics and material science
PDF
Neural network for molecular dynamics simulation and design of 2D materials
PDF
Data-driven approaches to studying protein-DNA interactions from a structural point of view
PDF
Healing of defects in random antiferromagnetic spin chains
PDF
Temperature-dependent photoionization and electron pairing in metal nanoclusters
PDF
Analysis and algorithms for distinguishable RNA secondary structures
PDF
Plasmonic excitations in nanostructures
PDF
Hot electron driven photocatalysis and plasmonic sensing using metallic gratings
PDF
Entanglement parity effects in quantum spin chains
PDF
Applications in optical communications: quantum communication systems and optical nonlinear device
PDF
Disordered quantum spin chains with long-range antiferromagnetic interactions
PDF
Photoexcitation and nonradiative relaxation in molecular systems: methodology, optical properties, and dynamics
PDF
Synthesis, characterization, and device application of two-dimensional materials beyond graphene
PDF
Quantum information-theoretic aspects of chaos, localization, and scrambling
Asset Metadata
Creator
Groenewald, Roelof Erasmus
(author)
Core Title
Coulomb interactions and superconductivity in low dimensional materials
School
College of Letters, Arts and Sciences
Degree
Doctor of Philosophy
Degree Program
Physics
Publication Date
10/17/2018
Defense Date
09/12/2018
Publisher
University of Southern California
(original),
University of Southern California. Libraries
(digital)
Tag
ab-initio,OAI-PMH Harvest,plasmons,SCDFT,superconductivity,TMDC's
Format
application/pdf
(imt)
Language
English
Contributor
Electronically uploaded by the author
(provenance)
Advisor
Haas, Stephan (
committee chair
), Di Felice, Rosa (
committee member
), Kresin, Vitaly (
committee member
), Nakano, Aiichiro (
committee member
), Zanardi, Paolo (
committee member
)
Creator Email
groenewa@usc.edu,regroenewald@gmail.com
Permanent Link (DOI)
https://doi.org/10.25549/usctheses-c89-81052
Unique identifier
UC11668911
Identifier
etd-Groenewald-6861.pdf (filename),usctheses-c89-81052 (legacy record id)
Legacy Identifier
etd-Groenewald-6861.pdf
Dmrecord
81052
Document Type
Dissertation
Format
application/pdf (imt)
Rights
Groenewald, Roelof Erasmus
Type
texts
Source
University of Southern California
(contributing entity),
University of Southern California Dissertations and Theses
(collection)
Access Conditions
The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...
Repository Name
University of Southern California Digital Library
Repository Location
USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA
Tags
ab-initio
plasmons
SCDFT
superconductivity
TMDC's