Theory and computation of hot carriers generated by surface plasmon polaritons in noble metals

Hot carriers (HC) generated by surface plasmon polaritons (SPPs) in noble metals are promising for application in optoelectronics, plasmonics and renewable energy. However, existing models fail to explain key quantitative details of SPP-to-HC conversion experiments. Here we develop a quantum mechanical framework and apply first-principles calculations to study the energy distribution and scattering processes of HCs generated by SPPs in Au and Ag. We find that the relative positions of the s and d bands of noble metals regulate the energy distribution and mean free path of the HCs, and that the electron–phonon interaction controls HC energy loss and transport. Our results prescribe optimal conditions for HC generation and extraction, and invalidate previously employed free-electron-like models. Our work combines density functional theory, GW and electron–phonon calculations to provide microscopic insight into HC generation and ultrafast dynamics in noble metals.

S urface plasmon polaritons (SPPs) are electron collective excitations generated by light at the interface between a metal and a dielectric 1 . Modelling of SPPs is dominated by approaches that employ classical electromagnetism and account empirically for the properties of the materials supporting the SPP. Yet, the inherently quantum mechanical nature of SPPs 2,3 becomes manifest in their scattering and decay processes in bulk and nanoscale materials. A key example is the decay of SPPs to electron-hole pairs, a process whose crucial importance is twofold: first, it is a main energy loss mechanism of SPPs in metals, currently limiting the applicability of plasmonics 4,5 ; second, it leads to generation of hot electrons and holes with application in several branches of applied physics 6 . Recent experiments reported the extraction of such SPP-induced hot carriers (HCs) in metals before they thermalize-typically over nanometer lengths and sub-picosecond times-by injecting the HCs over a Schottky barrier into a semiconductor or oxide [7][8][9] , or by transferring them into surface adsorbates to perform chemical reactions with large activation barriers 6,10,11 . The vast literature on HCs in noble metals generated by intense light pulses [12][13][14] (instead of SPPs), together with recent calculations on SPP-to-HC conversion 15,16 , help interpret these recent experiments on SPPinduced HCs.
Lack of predictive and accurate quantum mechanical approaches to study HCs generated by SPPs or photons has led to ambiguity in the microscopic interpretation of experiments involving HCs in Au and Ag, two mainstay materials in plasmonics 4 . In particular, recent studies indicate that SPP decay excites electronic transitions from occupied states close to the Fermi energy (E F ) 10 , implying that most of the SPP energy goes into hot electrons rather than hot holes, and that nanocrystals with o10 nm diameter are necessary to obtain HCs with significant energy 17 (for example, 1-2 eV away from E F ). The mean free paths (MFPs) of HCs in noble metals appear to be short 10 , though a quantitative estimate of the MFPs and their dependence on crystal direction is not available. In addition, it is often assumed that HC scattering and thermalization are dominated by Auger processes 10,14 , since thermalization induced by phonons occurs on a slower time scale 14 . As shown in this work, all these conclusions must be revised.
Here, we develop a quantum mechanical framework to study the energy distribution of HCs generated by SPPs and photons in Au and Ag, and employ ab initio calculations of the electron-phonon (e-ph) and electron-electron (e-e) interactions to study the MFPs and relaxation times of HCs within 5 eV of the Fermi energy. Our approach is free of experimental input, and combines density functional theory (DFT) 18 , the GW (where G is the Green function, W is the screened Coulomb potential, and GW is the diagram employed for the electron self energy) method 19 , and ab initio e-ph calculations 20 . Our ability to use extremely fine grids for Brillouin zone (BZ) integrations allows us to resolve HC scattering with unprecedented accuracy 21 . We find that the interband transition threshold (between d and s states) defines two regimes for HC generation and transport. The decay of SPPs with energy lower than the interband threshold leads to generation of long-lived HCs with long isotropic MFPs of up to 40 nm and energy within 1-2 eV of E F . On the other hand, decay of SPPs with energy higher than the interband threshold leads to generation of short-lived hot holes in d states with anisotropic and short (sub 5 nm) MFPs, and hot electrons with only o1 eV energy above E F . The regime characterized by SPP energy below the interband threshold is better suited to employ HCs in applications requiring long MFPs, and allows one to optimize HC generation by tuning the SPP energy. These results represent an important first step to understand SPP decay and energy loss, and to control HC generation and transport in bulk and nanostructured noble metals.

Results
Theory of SPP-electron coupling. Coupling of SPPs or photons to electrons in materials can be described in the framework of many-body perturbation theory 22 . The lowest order Feynman diagram for this coupling process, shown in Table 1, describes a boson (here, SPP or photon) coupling to the electron gas through the electronic polarizability w and a coupling matrix element g. This diagram is analogous to the lowest order coupling of phonons to electrons 23,24 . The decay rate G(q p ,E p ) for a boson of momentum :q p and energy E P to electron-hole pairs is proportional to the imaginary part of the self-energy 22 : where À E n,k and E n 0 ;k þ q p are, respectively, the quasiparticle energies of the hot hole in a state with band n and crystal momentum -k, and hot electron in a state with band n 0 and crystal momentum k þ q p , as generated in the boson decay process. In addition, f n,k and f n 0 ;k þ q p are Fermi occupations, Z is a small broadening, and the sum is contributed by electronic transitions from occupied to empty states that differ by :q p in crystal momentum and E p in energy. If the matrix elements have a weak band and k dependence, the decay rate of a SPP or photon is proportional to the imaginary part of the electronic polarizability, Imw(q p ,E P ).
Here, we employ a quantity proportional to Imw(q p ,E P ), the finite momentum joint density of states (FM-JDOS) J(q p ,E P ), defined as the number of states per unit energy separated by momentum :q p and energy E p in the quasiparticle bandstructure: Photon-electron g PH n;n 0 ;k SPP-electron g SPP n;n 0 ;k 2 ¼ e" h=m ð Þ 2 1=ð2E eff E p AL z Þ ½ e TM q p Á p n;n 0 ;k q p À Á 2 SPP, surface plasmon polariton. Square modulus of the coupling matrix element, g j j 2 , and Feynman diagram for the self-energy of the boson in the lowest-order interaction for the three bosons considered here-plasmons, photons and SPPs-with electrons in a material. These interactions share a common diagram (right column) consisting of a polarizability bubble w multiplied by two coupling vertices g. The wiggly lines on the two sides of the self-energy loop g 2 w indicate the incoming and outgoing bosons.
The FM-JDOS is a measure of the phase space available for the decay of a SPP or photon with momentum :q p and energy E p . Greater FM-JDOS values generally would correspond to greater SPP decay and HC generation rates.
The coupling matrix elements in equation (1) can be derived with a range of theoretical approaches. We first discuss the coupling of electrons to plasmons-that is, collective electronic excitations in the absence of light-since theories developed for plasmons have been occasionally applied to SPPs 25 . The so-called Landau damping 26 treatment of the electron gas derives the decay rate G p of a bulk plasmon in the weak coupling limit via the zeros of the complex macroscopic dielectric function e M of the metal close to the real energy axis, that is, E M ðq p ; E p þ i' G p ðq p ; E p ÞÞ ¼ 0. Due to the bosonic nature of plasmons, the decay rate G(q p ,E p ) has the form given in equation (1), with coupling matrix elements 26 g PL j j 2 ¼ 2pe 2 E p =q 2 p . In spite of its common use, Landau damping is rarely discussed in terms of the diagram in Table 1 and the coupling g PL , though this viewpoint elucidates its physical origin.
The coupling matrix elements for SPPs and photons can be derived by second quantization of the vector potential 27 (see Methods). For photons, the coupling is given by the well-known dipole matrix element jg PH n;n 0 ;k j 2 ¼ ðe' =mÞ 2 1=ð2EE p VÞ Â Ã je q p Á p n;n 0 ;k ðq p ¼ 0Þj 2 , where V is the volume, e q p is the polarization vector, and p n;n 0 ;k ðq p Þ the transition dipole matrix element. For SPPs at a metal-dielectric interface, we apply a vector potential second quantization procedure 28,29 (see Methods) to obtain the SPP-electron coupling jg SPP n;n 0 ;k j 2 given in Table 1. The SPP-electron coupling has a form analogous to that for photons, with two important differences. First, the volume V is replaced by the SPP-field volume V SPP ¼ A Á L z , where A is the metal-dielectric interface area and L z the decay length of the SPP in the metal. Since L z is of order 10 nm at visible wavelength for noble metals 1 and thus much shorter than the light absorption depth, the electric field associated with the SPP is concentrated to a small volume V SPP at the metal-dielectric interface 1 , resulting in enhanced local generation of HCs from SPPs compared with light. Second, the SPP wavevector q p is in general comparable to the size of the BZ, while for photons q p E0. The finite SPP wavevector can significantly alter the phase space available for SPP decay, as dictated by the q p dependence of the FM-JDOS. We remark that additional momentum can be transferred by phonons and defects assisting the SPP scattering process.
For both SPP and photon decay, the FM-JDOS regulates the phase space for HC generation. Since the experimental conditions and the material-specific SPP dispersion relation dictate the SPP energy E p and wavevector q p 1,30,31 , we study the FM-JDOS as a function of q p and E p and arbitrarily treat them for the purpose of this work as independent variables; HC generation from light is obtained as the specific limit q p ¼ 0. We note that the SPP decay rate G p (q p ,E p ) could also be computed by evaluating the coupling matrix elements explicitly for fixed experimental conditions and dispersion relation 1 . This approach would enable first-principles calculations of plasmonic losses and the MFPs of SPPs in materials 4,5 .
HC generation. To study HC generation from SPPs at a metal-dielectric interface (see Fig. 1a), we compute the FM-JDOS for Au and Ag, as shown in Fig. 1b, for SPP energies up to 5 eV and increasing SPP wavevectors along the G-L direction in the BZ (that is, for SPPs propagating along the [111] crystallographic direction). Our calculations are carried out using GW bandstructures (see Methods).
Detailed features of the bandstructures of Au and Ag, shown in Fig. 2, regulate HC generation and transport in these materials. The main features of the Au and Ag bandstructures are the presence of relatively dispersionless, occupied d bands with large associated electronic DOS, straddled by a free-electron-like s band. The top of the d bands is located at energy E int below the Fermi energy, where E int is the threshold for interband (d to s) electronic transitions 32,33 . From our GW calculations, we find E int values of 2.3 eV in Au and 3.7 eV in Ag (see Fig. 2), in excellent agreement with experiment 5,32,33 . The free-electron-like s band spans an energy window between E int below and over Energy (eV) 4 eV above E F , and hence dominates the electronic properties of Au and Ag near E F . The FM-JDOS in Fig. 1b follows similar trends in Au and Ag, with deviations due to the different E int values in the two materials. We first discuss the FM-JDOS with q p ¼ 0, which quantifies the phase space for HC generation for both SPPs with very small wavevectors and photons. The FM-JDOS curve with q p ¼ 0 (G) has a zero value for SPP energies up to E int , indicating the absence of possible SPP decay into HCs in this energy range. For SPP energies greater than E int , the FM-JDOS increases linearly with energy, indicating that a large phase space opens up for HC generation. In this regime, SPP damping is strong, and HC generation is intense and dominated by formation of hot d holes collecting most of the energy from the SPP, as discussed below.
For small SPP wavevectors along the G-L direction up to B0.2 L, we find a peak in the FM-JDOS corresponding to SPP decay through intraband transitions within the s band. The energy position of the peak increases linearly with SPP wavevector, and for large enough wavevectors merges with the FM-JDOS above E int . In Au, the intraband peak is below E int for q p values up to B0.2 L, and in Ag for q p values up to B0.3 L due to the higher E int value in Ag. The presence of an intraband peak associated with s-s transitions can be understood by examining the energy difference between states with momentum k and k þ q p in the free-electron-like s band: For small wavevector q p , the term quadratic in q p can be neglected, while the term linear in q p increases by B' 2 p 2 ffiffi ffi 3 p = 10ma 2 ð Þ % 1 eV (a is the lattice constant) for an increase of q p by 1/10 the G-L distance. The computed FM-JDOS confirms these trends, as seen by the increase in the energy of the intraband peak by B1 eV in going from 0.1 L to 0.2 L and from 0.2 L to 0.3 L. In this small wavevector regime, generation of HCs is appreciable only if the SPP energy and wavevector are matched to the intraband FM-JDOS peak, namely if E p (q p )E(' 2 k Á q p )/m for a set of k points in the BZ. When this condition is met, intense HC generation from intraband transitions can occur, resulting in an almost equal apportioning of the SPP energy between hot electrons and holes, as discussed below.
For SPP wavevectors larger than B0.3 L, the FM-JDOS resembles that for q p ¼ 0 except from a flat tail that extends to low SPP energy, as shown by the dashed lines in Fig. 1. This tail is due to generation of HCs from intraband transitions with large transferred momentum from the decay of SPPs with wavevector comparable to the size of the BZ. The FM-JDOS for q p 40.5 L is unchanged and nearly identical to the 0.5 L case.
To summarize, we find two HC generation regimes. For SPPs with energy E p oE int , optimal generation can be achieved using SPPs with relatively small wavevector (less than B1/5 the BZ size) that are matched with the FM-JDOS peak, while HC generation is relatively weak for SPPs with larger momentum. For SPPs with energy E p 4E int , HC generation is intense regardless of the SPP momentum. We find identical trends for SPPs with wavevector in the G-X high-symmetry direction, that is, propagating along the [100] direction. Removal of the electron gas approximation used in simplified treatments of SPPs in noble metals 17 is necessary to capture the delicate interplay found here between the energetics of the s and d states and intraband (s-s) and interband (d-s) transitions.
Energy distribution of HCs. Figure 2 shows the energy distribution of HCs generated by specifically assumed SPPs with energy increasing in small steps from B1 to 5 eV and a range of wavevectors in the G-L direction, along with the GW bandstructures. From now on throughout the manuscript, all HC energies are referenced to the Fermi energy (E F ), and the energy of the hot hole is À E, that is, the direction of increasing energy in the bandstructure is downward for holes and upward for electrons. The energy distribution histograms identify the HCs generated in the two regimes discussed above. For SPPs with energy lower than E int , hot holes and electrons are generated with a roughly uniform energy distribution as a result of the s-s intraband transitions. The HC energy distributions shown in Fig. 2 in this energy range pertain to SPPs with energy and wavevectors matched to the FM-JDOS intraband peaks (see Fig. 1). While the details of the energy distribution in this regime are sensitive to the SPP energy E p , the overall trends indicate that HCs are generated with an average energy of BE p /2.
For SPPs with energy greater than E int , we find a change to a different HC generation regime. At the onset of the interband transitions for SPP energy of 2.3 eV in Au and 3.7 eV in Ag, on SPP decay the hot holes absorb most of the SPP energy, creating a spike in the hot hole population at the top of the d band. On the other hand, hot electrons have only modest energies in this regime. At higher SPP energies of B5 eV, the HC population become approximately uniformly distributed over a wide energy range, although a peak in the hot d hole population is still present in Ag due to its higher E int . These trends are common to both Au and Ag, and we have verified that Cu, a material less commonly used in plasmonics 5 and not discussed here, follows identical HC generation and energy distribution trends as Au and Ag.
HC scattering processes. While ample experimental 13,32 and theoretical 14 data exists on HCs in noble metals, characterization of the MFPs and relaxation times of HCs using ab initio theory has been limited by the absence of accurate e-ph calculations including all phonon modes over the entire BZ. First-principles calculations combining e-ph and e-e interactions are necessary to characterize the anisotropic and energydependent MFPs of HCs in materials, as shown in our recent work on HCs in semiconductors 21 . We compute the scattering rate (and its inverse, the relaxation time) for the e-ph and e-e interactions from the imaginary part of the respective selfenergies, ImAE e À ph and ImAE e À e , and resolve these quantities for each electronic state with band n and k-point on very fine BZ grids (see Methods and ref. 21). The total relaxation times t n,k , combining the e-ph and e-e interactions, are obtained as t n;k ¼ 2=' ImAE e À e n;k þ ImAE e À ph n;k h i À 1 . Figure 3 shows the total relaxation times t n,k for HCs in Au and Ag with energy up to 5 eV. The total relaxation times decay rapidly away from E F following a volcano-shaped trend, with the peak of the volcano centred close to E F . In Ag, we find an additional peak in the relaxation times at the top of the d bands, in agreement with recent photoemission experiments 34 . Analysis of the e-ph and e-e scattering rates (see the bottom panels in Fig. 3) highlights the differences between the energy windows spanned by the s and d states, and elucidates the origin of the relaxation time trends. In the energy range spanned by the s states, the e-e scattering rates (that is, ImAE e À e ) form a parabola with minimum value of zero at E F as predicted by Fermi liquid theory 14 , while the e-ph scattering rates are relatively constant and exhibit a minimum at the onset of the d states and a maximum 1-2 eV above E F . The e-e scattering rates become greater than the e-ph B2 eV away from E F , thus indicating that Auger and impact ionization processes included in the e-e term dominate HC scattering only 2 eV or more away from E F , while e-ph dominates HC scattering within 2 eV of E F . Combining the ARTICLE two scattering mechanisms leads to total relaxation times with a broad maximum centred at E F and a rapid decay 1-2 eV away from E F as the parabolic e-e scattering rate becomes dominant (see Fig. 3). In the energy window spanned by the d bands, the large increase in the electronic DOS causes very strong e-ph scattering, and the localized nature of the d states leads to deviations from the free-electron parabolic trends for e-e scattering. The total relaxation times show that hot holes arising from d states lose energy on a sub-5-fs time scale, thus making such hot d holes challenging to extract before thermalization. Finally, our calculations yield comparable time scale and scattering rate for e-ph and e-e scattering, in contrast with the common heuristic assumption that scattering by phonons occurs on a much slower time scale 14 . We thus conclude that previous models including only e-e scattering 14 and missing the important e-ph component are incomplete for understanding HC relaxation times and photoemission linewidths in noble metals.
To validate our results, we compare our relaxation times at the Fermi energy in Au and Ag with those inferred from carrier transport measurements. Kopitzki 35 combined the Drude model with room temperature resistivity measurements, and obtained semiempirical relaxation times of 29 fs for Au and 41 fs for Ag, in close agreement with those obtained by Ashcroft and Mermin 36 with the same approach. Averaging our relaxation times over a small energy window around E F yields relaxation times of 24 fs for Au and 27 fs for Ag. While the agreement is excellent for Au and good for Ag, our data shows a large spread in the relaxation times at E F , a feature not captured by the Drude model and stemming from the Fermi surface anisotropy. It is therefore puzzling that the Drude relaxation times 33,35 are widely used in plasmonics 4,5 , given that their physical meaning is questionable when applied to non-equilibrium situations involving SPPs. We conclude that further investigation on this point is necessary, and remark that ours is the first truly ab initio determination of the relaxation times on the Fermi surface of noble metals, to be discussed in detail elsewhere.

Discussion
The data obtained so far allow us to formulate optimal conditions to generate HCs and extract or utilize them before thermalization. Figure 5 summarizes HC generation and transport in Au and Ag by defining two regimes for HC generation from intraband (E p oE int ) and interband (E p 4E int ) transitions induced by SPP decay. The energy distribution, relaxation time and MFP of the generated HCs are dramatically different in these two regimes. For E p oE int , the HCs consist of s holes and s electrons with average energies up to E int /2 (that is, 1-2 eV) generated by weakly damped SPPs due to the small FM-JDOS values. HC generation can be optimized by matching the energy and momentum of the SPP to the FM-JDOS intraband peak described above and shown in Fig. 1. The HCs in this energy range possess relatively long MFPs of B10-40 nm, and can thus be extracted before thermalization over B50 nm lengths in ideally pure samples. This regime is thus optimal for HC generation and extraction. For hot electrons, we predict that the (110) surface would be best suited to extract HCs due to the longer MFPs in the [110] direction. For hot holes with B1-2 eV energy, the (100) and (110) surfaces are predicted to enable optimal extraction. Our results further suggest that Ag may be better suited than Au for HC generation due to the wider energy window for intraband transitions, thus motivating studies of HC generation from SPPs in Ag. On the other hand, the strong SPP damping regime with E p 4E int is non-optimal for extracting HCs, since in spite of strong HC generation the SPP energy is mostly transferred to short-lived hot d holes with B1 nm MFP, while hot electrons possess only low energy of o1 eV. These conditions prevent HC extraction in this regime, unless a nanometer-thick thin film or nanostructured metal is employed.
Our data further show that nanostructuring the metal to sub-10-nm size is not required to obtain energetic HCs, and not needed to extract them unless HCs with 42 eV energy are desired. In addition, for SPP with energy lower than E int , the SPP energy is equally distributed between hot electrons and holes, and energy loss is dominated by e-ph rather than Auger scattering. These findings address the misconceptions emerged in recent work as outlined above, and show that optimal HC generation is possible by carefully tuning the SPP energy and wavevector at noble metal-dielectric interfaces. Finally, we note that our approach is for now limited to SPP at metal-dielectric interfaces, while localized surface plasmons in nanostructures need to be treated differently due to their localized nature 39 . We believe that the method presented in this work can be extended to study HC generation from localized surface plasmon modes, which will be the subject of future investigation.
In summary, we establish a theoretical framework to study SPP damping and HC generation and transport on the same footing using many-body perturbation theory. Our first-principles calculations can accurately describe SPP-induced generation and ultrafast scattering of HCs in noble metals of use in plasmonics, photocatalysis, photovoltaics and optoelectronics. Our work highlights the interplay of the s and d bands in noble metals, and prescribe optimal experimental conditions for generation and extraction of HCs. Our approach paves the way to first-principles calculations of SPP losses in materials 4 .

Methods
First-principles calculations. We carry out ab initio calculations on Au, Ag and Cu in the face-centred cubic structure with lattice parameters of 7.72, 7.71 and 6.82 bohr, respectively. The ground-state electronic structure is computed within the local density approximation (LDA) 40 of DFT using the QUANTUM ESPRESSO code 41 . We use norm-conserving pseudopotentials (which include semicore s and p states) to describe the core-valence interaction 42 , together with a plane-wave basis set with kinetic energy cutoff of 60 Ry for Ag and Au, and 240 Ry for Cu. The e-e (that is, GW) and e-ph self-energies are computed separately, and then combined together.
The real and imaginary parts of the GW self-energy 19 are computed using the Berkeley-GW code 43 . The real part ReðAE GW n;k Þ is obtained with a generalized plasmon-pole calculation 19 on a 12 Â 12 Â 12 k-point grid and then interpolated using Wannier functions (see below). Kinetic energy cutoffs of 50 Ry and 120 Ry are used, respectively, for the screened and bare Coulomb interactions, together with B1,000 empty bands and a 8 Â 8 Â 8 q-point grid for the dielectric screening. The DFT eigenvalues E n;k are corrected with GW self-energies to obtain the quasiparticle bandstructures, E n;k ¼ E n;k þ ReðAE GW n;k Þ, used in all calculations in this work. The imaginary part of the GW self-energy, ImðAE GW n;k Þ, is computed using fullfrequency GW calculations 43 . Here ImðAE GW n;k Þ are evaluated on-shell at the LDA eigenvalues E n;k , and then plotted versus the corresponding GW eigenvalues E n,k . We use 20 Â 20 Â 20 k-point grids to converge ImðAE GW n;k Þ within B10 meV. In addition, we use kinetic energy cutoffs of 20 Ry and 120 Ry, respectively, for the screened and bare Coulomb interactions, together with B100 empty bands and a 20 Â 20 Â 20 q-point grid for the dielectric screening.
The imaginary part of the e-ph self-energy, ImðAE e À ph n;k Þ, is computed using the EPW code 20 . We compute the self-consistent potential and Kohn-Sham states on a 12 Â 12 Â 12 k-point grid using DFT, and lattice-dynamical properties with DFPT 44 on a 4 Â 4 Â 4 q-point grid. The e-ph matrix elements are first computed on these coarse grids, and then obtained on significantly finer grids using an interpolation procedure based on Wannier functions as implemented in the EPW code 20,24 . Maximally localized Wannier functions 45 are obtained starting from five d orbitals on the metal atoms and two s orbitals, each in a tetrahedral interstitial point of the face-centred cubic lattice, for a total of seven wannierized bands (we skip the semicore states). The Wannier interpolated and DFT bandstructures agree within 5 meV in an energy window up to 5 eV above the Fermi energy. The Blochto-Wannier rotation matrices are then used to interpolate the GW bandstructures, which are used to compute ImðAE e À ph n;k Þ. The fine grids for the calculation of ImðAE e À ph n;k Þ consist in a 40 Â 40 Â 40 k-point grid, and up to 1 million random phonon q-points in the BZ. Such fine grids allow us to converge ImðAE e À ph n;k Þ within 1 meV. Further details of our approach to compute the e-ph self-energy are discussed in ref. 21.
The FM-JDOS calculations are carried out using Wannier interpolated GW bandstructures, a 100 Â 100 Â 100 k-point grid and a small Lorentzian broadening ZE25 meV. The energy distribution of the generated HCs in Fig. 2 are obtained by counting the momentum-and energy-conserving transitions in equation (2) as a function of occupied band (for hot holes) and empty band (for hot electrons) quasiparticle energies.
The MFPs are obtained using velocities computed from the slope of the GW bandstructure along the given high-symmetry direction, together with total relaxation times t n,k combining the e-ph and e-e relaxation times with Matthiessen's rule. The total relaxation times are obtained from ImðAE e À ph n;k Þ computed on grids of over 100 k points along each high-symmetry direction, and ImðAE GW n;k Þ computed on a 20 Â 20 Â 20 k-grid and then averaged over k points to generate an energy-dependent ImðAE GW n;k Þ evaluated at each energy E n,k for which a calculation of ImðAE e À ph n;k Þ is performed. ARTICLE Derivation of the coupling matrix elements. The derivation of the plasmonelectron matrix element in Table 1 within the Landau damping approximation is given in ref. 26 and will not be discussed. Here we derive the coupling matrix elements for the photon-electron and SPP-electron interactions given in Table 1.

Energy
We use SI units throughout the derivations. The photon case is well known 27 and discussed here as a starting point for the SPP case. The second quantized vector potentialÂ r; t ð Þ for a photon with momentum :q p , energy E q p ¼ ' o p and polarization unit vector e q p can be obtained by associating the amplitude with a destruction operatorâ q p and its complex conjugate with a creation operatorâ y q p . We use the vector potential with unit amplitude since we are interested in the coupling matrix element, and thus the transition rate per unit of incident power. It reads 27 where f q p is a normalization constant. Using the Coulomb gauge rA q p ¼ 0, the second quantized electric field can be obtained fromÊ ¼ À @Â=@t.
The normalization constant f q p is determined by equating the energy of the classical and quantized fields. The classical (cycle averaged) energy hU q p i for a photon with momentum :q p in a homogeneous and isotropic material with where V is the volume of the material. This energy needs to be equal to its quantum counterpart, ' o q p ðn q p þ 1=2Þ. In addition, hU q p i can be obtained by averaging the quantized electric fieldÊ q p ¼ À @Â q p =@t over the quantum state jn q p i with n q p photons. We get: where in the last equation we used the time derivative of equation (4) the electric field and the commutation rules for boson creation and annihilation operators, (5) and (6) yields a normalization constant , and a vector potential: The photon-electron interaction Hamiltonian is H int ¼ e mÂ q p Áp, as obtained by expanding the kinetic energy jp þ eÂ q p j 2 = 2m ð Þto first order in the vector potential. The second quantized fields associated with the electronic Bloch states are: whereĉ n;k andĉ y n 0 ;k 0 are, respectively, fermion annihilation and creation operators satisfying the anticommutation ruleĉ n;kĉ y n 0 ;k 0 þĉ y n 0 ;k 0ĉn;k ¼ d n;n 0 d k;k 0 . The second quantized interaction Hamiltonian can be written as: Crystal momentum conservation requires k 0 ¼ ± q p , where the plus and minus sign apply for the first and second integral above, respectively. The interacting Hamiltonian becomes: e q p Á p n;n 0 ;k q p e À ioq p tĉy n 0 ;k þ q pĉ n;kâq p þ h:c: where h.c. is the hermitian conjugate, and we defined the momentum matrix element as the integral p n;n 0 ;k ðq p Þ ¼ 1 O R dru Ã n 0 ;k þ q p r ð Þðe q p ÁpÞu n;k r ð Þ over the unit cell volume O, as commonly used in the study of optical absorption 46 . We obtain the photon decay rate G q p by applying Fermi's golden rule. Photon absorption occurs with a transition rate G abs q p between an initial state j0; 1 q p 4, consisting of the ground electronic state at room temperature and one photon with momentum q p , and a final state |S;0i with an electron in the Bloch state |n 0 ,k þ q p i, a hole in the Bloch state |n,ki and no photons. The absorption rate is given by: Only the first term in equation (8) gives a non-zero contribution, and we get: X n;n 0 ;k ' e 2 2Vm 2 eoq p ! eq p Á p n;n 0 ;k q p 2 fn;k 1 À fn0;k þ q p d ' oq p þ En;k À En0;k þ q p The inverse process of photon emission from an initial state |S;0i to a final state j0; 1 q p i has a rate G em q P given by: The net rate G q p ¼ G abs q p À G em q p is the decay rate of a photon with momentum :q p and energy E p ¼ ' o q p into a hot electron-hole pair. Writing the delta function as a Lorentzian with broadening Z-0, namely d x ð Þ ¼ 1 p Im 1 x À iZ , we obtain the decay rate as: This expression is identical to equation (1) with the squared coupling matrix elements jg PH n;n 0 ;k j 2 in Table 1. Since for photons q p E0 compared with the size of the BZ, the squared coupling matrix elements are: 2EE p V e q p Á p n;n 0 ;k q p ¼ 0 2 as shown in Table 1. We note that for slowly varying matrix elements in the BZ, the photon decay rate is proportional to the zero-momentum JDOS commonly used in the theory of the optical properties of solids 46 , that is, The derivation of the coupling matrix elements jg SPP n;n 0 ;k j 2 for SPPs follows similar steps to the photon case discussed above. Following Nkoma et al. 28 , we study the case of a SPP generated at the interface between a metal with frequency-dependent dielectric function E M o ð Þ and a dielectric material with frequency-independent dielectric constant E D . A similar second quantization procedure was also developed by Elson et al. 29 . It is well known that the intensity of the SPP wave decays exponentially in the direction perpendicular to the dielectric-metal interface 1,30 . To study HC generation, we focus here on the SPP field in the metal, and its associated vector potential. Since SPP possess transverse magnetic polarization, their wavevector can be taken as real in the xy plane-here the plane of the metaldielectric interface, with area A-and purely imaginary in the z direction normal to the interface. With the metal in the z40 space and the x-axis parallel to the propagation direction of the SPP wave, the SPP wavevector is q p ¼ (q p,x ,0,iq p,z ). While in the main text, we define q p ¼ (q p,x ,0,0) as the propagation wavevector in the xy plane, for the purpose of this derivation q p is defined differently. Here q p,x is a real positive number since the electron-induced dissipation effects we aim to find would be encoded in its imaginary part, and q p,z is related to the decay length L z of the SPP field intensity by 30 q p,z ¼ 1/(2L z ). The relative magnitude of the components of q p are determined by the dielectric properties of the metal and dielectric 1,30 , which impose q p;x =q p;z ¼ ffiffiffiffiffiffiffiffiffiffiffiffi ffi E D =E M p (here, the real parts of E D and E M should be employed for consistency). Similarly, the components of the unit polarization vector e TM q p are related by 30 The vector potential of the SPP wave in the metal, with unit amplitude as in the photon case, can be written in second quantized form as: A q p r; t ð Þ ¼ f q p e TM q pâ q p e i q p Ár À opt ð Þ þ e TM q p Ãâ y where, as in the photon case, the normalization constant f q p can be determined by energy considerations. We write the classical electric field associated with the SPP in the metal as E ¼ E 0 e i q p Ár À oq p t À Á . Using a result in ref. 28 , the cycle-averaged energy in the SPP field-which accounts for contributions from both the field in the dielectric and the metal-can be written as hU q p i ¼ AE eff E 0 j j 2 =ð4q p;z Þ, where we defined an effective dielectric constant: We use the correspondence principle for a quantum state jn q p i with n q p SPP quanta. Similar to the photon case, the total energy in the SPP mode ' o q p ðn q p þ 1=2Þ can be equated to its classical counterpart, while the square modulus of the electric field jE q p j 2 can be replaced with its quantum mechanical average hn q p jÊ 2 q p jn q p i, withÊ q p ¼ À @Â q p =@t. We obtain the two equations: By equating equations (12) and (13), and using q p,z ¼ 1/(2L z ), we obtain the ARTICLE NATURE COMMUNICATIONS | DOI: 10.1038/ncomms8044 normalization constant f q p ¼ ½' =ð2AL z E eff o q p Þ 1=2 , and thus the SPP vector potential: e TM q pâ q p e i q p Ár À opt ð Þ þ e TM q p Ãâ y q p e À i q p Ár À opt ð Þ h i This vector potential is similar to the photon case, equation (7). The key differences are the decay in the z direction for the SPP field, which introduces an effective SPP volume AL z in the vector potential, the appearance of an effective dielectric constant for the medium and the fact that the polarization is transverse magnetic with a direction imposed by the properties of the metal and dielectric material.
Since the form of the SPP vector potential is completely analogous to the one for photons, the decay rate of SPPs to hot electron-hole pairs can be carried out following the same steps as in the photon case discussed above, with the substitutions V-AL z , e q p ! e TM q p , and allowing for a finite momentum 1 q p with magnitude q . This leads to the SPP-electron coupling matrix element in Table 1