Examining electron-boson coupling using time-resolved spectroscopy

Nonequilibrium pump-probe time domain spectroscopies can become an important tool to disentangle degrees of freedom whose coupling leads to broad structures in the frequency domain. Here, using the time-resolved solution of a model photoexcited electron-phonon system we show that the relaxational dynamics are directly governed by the equilibrium self-energy so that the phonon frequency sets a window for"slow"versus"fast"recovery. The overall temporal structure of this relaxation spectroscopy allows for a reliable and quantitative extraction of the electron-phonon coupling strength without requiring an effective temperature model or making strong assumptions about the underlying bare electronic band dispersion.


I. INTRODUCTION
Pump-probe spectroscopies offer the exciting opportunity to perturb and measure electrons and collective modes in solids on their intrinsic ultrafast time scales. Moreover, they provide a tool to gain insights into the behavior of matter pushed out of its thermodynamic equilibrium state, thereby probing microscopic details of complex many-body systems beyond effective thermodynamic variables. The opportunity, however, comes with a challenge: When the equilibrium concepts of entropy and temperature are overcome, we are also in need to develop a language for true nonequilibrium behavior. Moreover, the value of ultrafast spectroscopies must be enhanced by the ability to relate measurements out of equilibrium to microscopically derivable quantities in an effort to understand dynamics and, ultimately, nonequilibrium itself. In this paper, we address these issues by showing a full nonequilibrium solution of a generic electron-phonon coupled model system, and by relating the temporal spectroscopic response to microscopic quantum-many body details of the system. We utilize this concept as a basis to directly access the intrinsic decay rates of quasiparticles in the time domain without assumptions about bare band quantities that are needed in the frequency domain.
Angle-resolved photoemission (ARPES) has been widely used to determine the low-energy electronic structure of a variety of materials. 1,2 Many-body interaction effects are manifest via the self energy Σ giving rise to renormalizations of the electron mass and an overall broadening of the spectral function. As an example, it was thought that ARPES data for the cuprates could thereby be studied to extract the boson coupling that is responsible for high temperature superconductivity. [3][4][5][6][7][8][9][10][11][12][13][14][15][16][17][18][19][20] These investigations have provided important insights into the underlying bosons and interactions responsible for kinks in the electronic dispersion.
However, extracting the self-energy and a meaningful electron-boson coupling constant from frequency-dependent ARPES data can be challenging. In fact there is currently no reliable way to extract a frequency-dependent electron-phonon coupling. In addition, the interpretation of ARPES data in terms of a self-energy naturally rests upon the assumption of an underlying bare electronic band, which becomes renormalized by the coupling to bosonic modes. Similarly, the widths of momentum distribution curves (MDCs) need to be supplemented with a bare Fermi velocity in order to determine scattering rates. Estimates of the dimensionless electron-boson coupling strength λ therefore vary strongly in the literature, from very small effective scattering rates 21 to extremely strong coupling 18 .
As an alternative, a time domain approach using the relaxation of excited electrons back to equilibrium as a direct technique to measure electron-phonon coupling was proposed by Allen. 22 The basic quantities for this method are the deviations from the initial equilibrium values of the electron and phonon occupations in momentum space. A semi-classical Boltzmann approximation that averages over momenta and energy to extract a single energy relaxation time scale was used to extract estimates for electron-phonon couplings from pump-probe photoemission data in simple elemental metals, 23,24 but it was also applied to analyze pump-probe data for high-temperature superconductors. [25][26][27][28][29] However, in correlated materials, such as transition metal oxides like the cuprates, or in multiband systems with multiple disconnected Fermi surface pieces, such as the iron pnictides, the momentum and energy dependences must be resolved, since they are crucial ingredients for these materials' interesting properties. 30,31 In addition, a full quantum mechanical formulation is required to treat arbitrarily strong fields and strong energy and momentum dependent scattering, particularly in materials where the balance between various competing phases may be tipped with the application of strong electric fields. Examples include possible light-induced superconductivity in a stripe-ordered cuprate, 32 a dramatically enhanced conductivity in an oxide heterostructure, 33 coherent phonon oscillations serving as a marker for a photoinduced phase transition in VO 2 , 34 , the theoretical proposal of dynamical band flipping, 35 and pump-probe photoemission providing a novel classification scheme for charge-density wave insulators in the time domain. 36 In this paper, we combine the strength of the time domain approach with the power of ARPES self-energy microscopics. We present a microscopic quantum nonequilibrium approach with full momentum and energy resolution that can be applied when effective temperature phenomenology and the semiclassical Boltzmann equation fail, and that can be used in conjunction with standard many-body and material specific computational approaches. Using our model system, we show how the time (t), momentum (k), and energy (ω) resolution provided by state-of-the-art tr-ARPES [25][26][27][28][37][38][39][40][41][42][43][44][45][46][47] experiments in conjunction with the dependence of extracted decay times on the initial equilibrium sample temperature provides a direct method to measure the dominant phonon energy and electron-phonon coupling. Hence, we show that relaxational spectroscopy in the time domain allows the possibility to extract detailed information about the system in equilibrium and lends itself to a simple interpretation that does not involve any assumptions about the bare electronic band dispersion and goes beyond the traditional two-temperature model phenomenology.
The paper is organized as follows. In Section I we outline the electron-phonon coupled model system and the calculation of double-time Green functions on the Keldysh contour and the time-resolved photoemission response. In Section III we present our results for the photoexcited electronic structure (III A) and the extraction of the quasiparticle self-energy and the implied electron-phonon coupling from the temporal relaxation of the electrons (III B). Section IV contains a discussion of the results and their implications. Additional considerations regarding the tradeoff between time and energy resolution, the breakdown of an effective temperature description, the connection between relaxation rates and the self-energy, and the additional effect of including electron-electron scattering are provided in a Supplementary Material.

A. Introductory remarks
For the formulation of the nonequilibrium problem, we employ the well-established Keldysh Green function technique. The single-particle properties of the electrons are described by the double time Green function G(t, t ) on the complex Keldysh contour. These encode the electrons' initial equilibrium state, their real time evolution during and after a perturbation via a pump laser pulse, and their memory of the initial state. The doubletime Green functions allow us to extract the relevant information (occupation of states, tr-ARPES response) about the excited electron system in a microscopic way.
The pump laser field of arbitrary temporal shape is included as a space-independent and time-dependent vector potential A(t) in a nonperturbative way via Peierls substitution. 48,49 This approach has been formulated earlier by Jauho and Wilkins 49 which we extend numerically to longer time scales and arbitrary field strengths. The coupling of the electrons to a nondispersive optical phonon mode is modeled by the well-established and controlled diagrammatic Migdal approximation. 50 The resulting integro-differential equations of mo-tion for the Green functions 51 are solved using massively parallel numerical integration. The code is parallelized easily since the equations for different momentum points decouple in the discretized Brillouin zone. An additional Fourier transform of the lesser component of the resulting interacting Green function is performed to extract the tr-ARPES spectra. 52 In this step, a shape function is introduced to take the probe pulse width into account.
Additionally, the momenta need to be time-shifted depending on the pump field to ensure gauge invariance. 53

B. Model
We study the Holstein Hamiltonian 54 where the individual terms are the electronic kinetic energy with dispersion (k), the phonon energy for dispersionless Einstein phonons of frequency Ω, and the electron-phonon interaction, respectively. Here c i is the bosonic phonon annihilation (creation) operator on lattice site i, and g is the momentum-independent electron-phonon coupling strength.
The starting bare nonequilibrium Green function 53 where n F (ω) = 1/(1 + exp(ω/T )) is the Fermi function at temperature T , t and t are times on the Keldysh contour C, and θ c (t, t ) is the Heaviside function on the contour, which equals 1 if t is later on the contour than t and 0 otherwise. (k) is the single-particle energy dispersion for a square lattice with nearest (V nn ) and next nearest neighbor (V nnn ) hopping, (k) = −2V nn (cos k x + cos k y ) + 4V nnn cos k x cos k y − µ, measured with respect to the chemical potential µ. Here, we use the convention thath = c = e = 1, and we work in the Hamiltonian gauge, i.e. the scalar potential is set to zero.
We incorporate the electron-phonon interactions in the Migdal limit, which is appropriate for weak coupling. In this case the electronic self-energy due to electron-phonon coupling is The bare phonon Green's function D 0 (t, t ) is 55 where n B (Ω) = 1/(exp(Ω/T ) − 1) is the Bose function.

C. Method
With the self-energy above, we solve the Dyson equation for the interacting Green func- This can be done by casting the Dyson equation as a matrix equation. As an alternative, the Dyson equation can be rewritten by expanding the integral with Langreth rules, which gives equations of motion for the Matsubara, retarded, real-imaginary, and lesser Green's functions. 51,56 This leads to a set of Volterra integro-differential equations that can be solved via standard numerical integration. For the results presented here, we have used the Volterra scheme. Its main advantages are that the equations are manifestly causal, and the required memory is reduced by the use of symmetries.
The pulse that is of direct interest to pump-probe experiments is, by nature, a propagating light pulse; this implies an oscillating field without a zero-frequency component. We model the pump via an oscillating vector potential along the (11) or Brillouin zone diagonal direction with a Gaussian profile, The tr-ARPES intensity for a probe pulse of width σ pr is computed from 52 where the integral is along real times, and each component k i of the momentum argument of where A i is the corresponding i-th component of the vector potential) to restore gauge invariance. 53 In the Supplementary Material, we illustrate the importance of the temporal probe pulse width σ pr for the effective energy resolution of spectral features.

D. Setup and model parameters
We consider 2D tight-binding electrons in a partially filled metallic band with parameters For the electric pump and probe fields, we choose parameters that allow us to compute the full time evolution on the Keldysh contour by the solution of the integro-differential equations with a reasonable computational effort. Specifically we choose a Gaussian field envelope of width σ = 4 femtoseconds (fs) and a center frequency of 0.5 eV, which is smaller than UV pump pulses around 1.5 eV, which are currently typically used in tr-ARPES experiments.
The smaller frequency affects the transient behavior quantitatively. However, it allows us to use a coarser time step discretization and therefore to observe the relaxation of the excited electrons for a longer time, while it does not affect our conclusions about relaxation effects.
The delay time is defined with respect to the center of the Gaussian pulse envelope. The electric field has a peak strength E max = 0.4 V/a 0 . While this field strength is at the high end of what can be currently achieved experimentally, this choice is motivated by the fact that it leads to relatively strong, clearly visible spectral redistributions whose temporal decay can be followed easily. Again, the details of the pump do affect the results, but not the conclusions. This is explicitly demonstrated in the Supplement by analyzing the relaxation times for different maximal pump field strengths.
The subsequent tr-ARPES probe pulse has a Gaussian envelope width of σ pr = 16.5 fs. This choice is motivated by the maximum contour length that can be reached with reasonable computational effort, since the earliest and latest tr-ARPES observation delay times for which one can extract spectra have to fulfill a 3 σ pr minimum distance to the actual earliest and latest times on the contour. In the Supplement, we also show the equilibrium tr-ARPES spectra for smaller σ pr , which enables a better time resolution but necessarily leads to worse energy resolution.
In In equilibrium at zero temperature, kinematics prevents electron relaxation unless the quasiparticles have sufficient energy to emit phonons, which is why this window is expected.
In Figure 1f, we show that raising the initial equilibrium temperature leads to smaller spectral changes compared to Figure 1d at the same time delay. Correspondingly, the k-integrated spectral changes shown in Figure 1g are significantly smaller than the ones in 1e (shown on the same vertical scale), although the phonon window effect is still visible.
We note here that this result is quite reminiscent of electronic relaxation measured in The breakdown of the "hot electron" picture becomes evident when we integrate the tr-ARPES spectra in the momentum window around the Fermi momentum between k 1 and k 7 (see Figure 1). The energy distribution curves in Figure 2c are not simply featureless decaying exponentials, as they would be if the "hot electron" model was valid, but show distinct higher energy peaks related to the pump characteristics. Although it is possible to extract an effective temperature from the linear slope of the logarithmic intensities in a very small energy window (< 0.1 eV) near the Fermi level (see Figure S3 in Supplement), this "temperature" does not provide a meaningful description of the photoexcited system.
In particular, it fails to account for the actual excess energy in the electronic subsystem due to the spectral features at higher energies. In fact, once this high spectral weight has dissipated, a fit to an effective Fermi distribution shows that the increase in the effective temperature is very small, much smaller than that obtained from forcing a fit a shorter time delays (see Supplement).
These nonthermal features are even more pronounced for the strong pump. This observation clearly shows that photoexcitation is not necessarily the same as heating, and that a full nonequilibrium modeling beyond effective temperature phenomenology is often essential for correctly describing pump-probe spectroscopies. Since we include only electron-phonon scattering in our simulation and neglect electron-electron scattering, the discrepancy between quasithermal behavior and the observed transient spectra is not too surprising. The inclusion of electron-electron scattering leads to a closer resemblance of spectra with two-temperature phenomenology (see additional calculations including electron-electron scattering in the Supplement), but this is not guaranteed to be the case. Moreover, even if it is the case, the single relaxation rate from two-temperature fits is not simply connected to energy dependent population decay rates, as we show in the Supplement. Therefore, since tr-ARPES experiments probe the full energy-and momentum-resolved population dynamics, the data analysis should go beyond effective temperature decay rates.

B. Extracting equilibrium system properties by pump-probe photoemission
We now focus on a quantitative analysis of energy-dependent relaxation rates for the weaker pump and the higher initial temperature T = 290 K. Figure 3a shows the pump- not depend too strongly on temperature, then it is weakly dependent on the distribution function, and as a result, we will see the nonequilibrium system relax back to equilibrium at the longer time scales via the equilibrium relaxation rates. These conditions are usually fulfilled for sufficiently high temperature. Care has to be taken when the temperature becomes so low compared to the lowest relevant bosonic mode energy that relaxation via the equilibrium self-energy is completely blocked at low binding energies. In this case, relaxation pathways of higher order in the electron-boson coupling become relevant for the long-time relaxation, corresponding to effects of higher order in the deviation from the equilibrium distribution.
We now show how more detailed knowledge about the electron-phonon scattering processes allows us to quantitatively extract the electron-phonon mode coupling strength from a comparison of the tr-ARPES relaxation rates Γ tr (ω) and the self-energy Σ(ω). To this end, we treat our tr-ARPES relaxation rates as "experimental data" and assume that we know the line shape of Im Σ(ω). Within the Migdal approximation to the Holstein model, the relaxation rate is proportional to g 2 and of the form Γ(ω) = g 2 φ(ω) (see Supplement).
The function φ(ω) reflects electron-phonon scattering processes and takes into account the available scattering phase space related to the electronic band structure via the local electronic density of states N (ω). This is the reason why an asymmetry in the band structure results in an electron-hole asymmetry in the relaxation rates.
The measured tr-ARPES relaxation rates Γ tr are naturally convolved with the probe pulse shape and therefore include an additional energy averaging, since there is a tradeoff between time resolution and energy resolution. Therefore we apply an energy averaging procedure to φ(ω) (see Supplement), which yields φ aver . The electron-phonon coupling strength is then extracted from g 2 = Γ tr (ω)/φ aver (ω). Figure 4b shows this ratio as a function of energy, which serves as an a posteriori estimate of g 2 given the numerical simulation data. The extracted mean value including an error estimate from the statistical standard deviation is g 2 = 0.021 ± 0.004 eV 2 , to be compared with the exact input value g 2 = 0.02 eV 2 . This agreement highlights the fact that a nonequilibrium measurement can serve as an accurate tool for extracting equilibrium system properties, a feat which may be difficult to accomplish by pure equilibrium measurements.

IV. DISCUSSION
The central result of this paper is that equilibrium ideas can be used to understand the nonequilibrium dynamics of a pumped system. This is facilitated by the connection between relaxation rates for the decay of photoexcited electrons and the equilibrium electron-phonon scattering rates. Importantly, both the dependence of the rates on the binding energy and on the initial equilibrium temperature affect the result, but not the conclusion. We have shown  (e) Intensity difference ∆I ω between photoexcited (89 fs) and equilibrium tr-ARPES data. (f),(g) EDCs and spectral changes (same as in (d),(e)) for a higher initial equilibrium temperature T = 290 K, showing that the time for recovery shortens for increased temperature. for a weak pump field (same as in Figure 1) with maximum field strength E max = 0.4 V/a 0 . b, Ten times stronger pump field (E max = 4 V/a 0 ). For a lattice constant a 0 = 3.4Å, these peak field strengths are ≈ 12 MV/cm ("weak pump") and 118 MV/cm ("strong pump"), respectively. c, Energy distribution curves integrated over momenta along the Brillouin zone diagonal (same momentum window k 1 to k 7 close to k F as indicated in Figure 1) for various pump-probe delay times.  For the interpretation of time-resolved photoemission experiments, it is of central importance to include the influence of σ pr , the temporal width of the probe pulse, 1 which is illustrated in the pump-probe setup in Figure S1. A narrower probe pulse implies a better time resolution but a worse energy resolution, and vice versa. This is illustrated in Figure   S2, which shows the probe pulse width dependence of tr-ARPES intensity maps and corresponding EDCs. Apparently, the necessary energy resolution to resolve the phonon kink and the sharpening of spectral changes to be inside the phonon range has to be achieved by trading off intrinsic time resolution. Figure S3 shows fits of Fermi functions with effective temperatures and chemical potentials to the low-energy structure of the tr-ARPES energy distribution curves. This figure serves to show that one can make such a fit in a small energy window and extract decay constants for the time-dependent temperatures and chemical potentials. However, these do not provide an adequate description of the excited system, which is not at thermal equilibrium. Notably, the time-dependent spectra deviate considerably from the fits already at slightly higher binding energies than the fit window. The figure therefore demonstrates the expected breakdown of an effective temperature description and calls for a more advanced data analysis. Additionally, we note that even if a Fermi function fit seems to work (see, e.g., Ref.

S2. BREAKDOWN OF EFFECTIVE TEMPERATURE DESCRIPTION
2), it is easily shown that the extracted decay time for the effective temperature is not related in a simple fashion to energy-dependent decay times for the spectral intensities.
This is due to the fact that the intensity at a given binding energy depends nonlinearly on the effective temperature in this case, namely through a Fermi function. Therefore it is doubtful which information can be gained from the temperature decay time.

S3. RELAXATION RATES
In Figure 1 of the main text, we showed how the relaxation process inside the phonon window is strongly temperature dependent. For low temperatures, the time scales become considerably longer and thus lead to a slower relaxation. Figure S4 shows the correspond-ing equilibrium relaxation times. Obviously, the times due to electron-phonon scattering are basically independent of temperature outside the phonon window. Inside the phonon window, however, they are strongly temperature dependent. The time scales for our model parameters are well below a picosecond for the higher equilibrium temperature T =290 K, but become as large as tens of picoseconds for T = 116 K.
The relaxation times shown in Figure S4 are obtained from the imaginary part of the equilibrium retarded self-energy according to where Im Σ R (ω) = −πg 2 N (ω + Ω) (n B (Ω) + n F (Ω + ω)) and N (ω) = L −1 k δ(ω − (k)) (the sum is over the two-dimensional Brillouin zone, with L the number of momentum points summed over) is the bare electronic density of states. 3 The relaxation rate Γ(ω) = 1/τ (ω) is thus of the form g 2 φ(ω), which is used for the extraction of g 2 from the tr-ARPES data in the main text. Here, the equilibrium self-energy is taken at the initial temperature of the pump-probe simulation.
For the relaxation rates shown in Figure 4 of the main text, we apply an additional averaging procedure to account for the energy resolution due to the finite width of the probe pulse. The averaged self-energy is obtained by convoluting its inverse with a resolution function: For the probe pulse width σ pr = 16.5 fs, the corresponding energy resolution is given by σ res = √ 2h/σ pr ≈ 57 meV. Note that this is in contrast to Matthiessen's rule for scattering times, which states that inverses of times should be averaged instead of the times themselves. The difference is that Matthiessen's rule applies for averages of scattering processes for various sources of scattering, all for the same particle. Here, however, we are interested in averages of time scales for different particles. The reason why one should average the relaxation times here, instead of the rates, is the following: In the long-time limit, the slowest decaying components determine the effective relaxation rate. To give those components with the longest time scales the largest weight, one has to average times, not rates. The same averaging procedure is applied to yield φ aver , which is used to extract g 2 in Figure 4b of the main text.
Finally, we demonstrate that the relaxation times for particular energies do not significantly depend on the pump field strength as long as one does not measure times too close to the application of the pump. To this end, we show in Figure S5 a comparison of the time scales for exponential fits to the time traces of the energy-dependent tr-ARPES intensities for the weak and strong pumps (see Figure 2 of main text). The times extracted for weak and strong pump pulses agree reasonably well. This overall agreement of the time scales demonstrates that the interpretation of the nonequilibrium dynamics in terms of equilibrium electron-phonon scattering rates is valid even beyond the weak pump regime.

S4. INCLUSION OF ELECTRON-ELECTRON SCATTERING
Here we present the effects of an additional electron-electron scattering channel by adding a second-order perturbative local Coulomb repulsion of strength U in the local approximation to the electron-phonon self-energy: Σ U (t, t ) = U 2 G 0 (t, t )G 0 (t, t )G 0 (t , t), where G 0 = k G 0 k is the bare local electronic Green function. Note that this representation is equivalent to summing the full momentum-dependent Σ U k (t, t ) over momenta. As Coulomb repulsion values representative of weak to moderate coupling we choose U = 0.25 eV and 0.5 eV for our calculations, which correspond to W/8 and W/4, where W is the bare bandwidth of the electrons.
The analysis of these results shows that a treatment of electron-electron scattering on the same level as we treat electron-phonon scattering supports our main results, namely the phonon window effect for slow versus fast relaxation ( Figure S6