Predicting charge transport in the presence of polarons: The beyond-quasiparticle regime in SrTiO$_{\mathrm{3}}$

TL;DR

This paper introduces a first-principles method to accurately predict charge transport in SrTiO$_{3}$ with strong electron-phonon interactions, revealing a transition from quasiparticle to incoherent transport regimes with temperature.

Contribution

It develops a cumulant diagram-resummation technique to model strong electron-phonon interactions and predicts charge transport beyond the quasiparticle approximation in polar materials.

Findings

Accurately predicts electron mobility in SrTiO$_{3}$ between 150-300 K.

Identifies a transition from band-like to incoherent charge transport with increasing temperature.

Provides a broadly applicable method for materials with strong electron-phonon interactions.

Abstract

In materials with strong electron-phonon ($e$-ph) interactions, the electrons carry a phonon cloud during their motion, forming quasiparticles known as polarons. Predicting charge transport and its temperature dependence in the polaron regime remains an open challenge. Here, we present first-principles calculations of charge transport in a prototypical material with large polarons, SrTiO$_{3}$. Using a cumulant diagram-resummation technique that can capture the strong $e$-ph interactions, our calculations can accurately predict the experimental electron mobility in SrTiO$_{3}$ between 150$-$300 K. They further reveal that for increasing temperature the charge transport mechanism transitions from band-like conduction, in which the scattering of renormalized quasiparticles is dominant, to a beyond-quasiparticle transport regime governed by incoherent contributions due to the interactions between the electrons and their phonon cloud. Our work reveals long-sought microscopic details of charge transport in SrTiO$_{3}$, and provides a broadly applicable method for predicting charge transport in materials with strong $e$-ph interactions and polarons.

Equations

$$G_{n\mathbf{k}}^R(t)=-i\theta (t)e^{-i{\epsilon }_{n\mathbf{k}}t}{e}^{C_{n\mathbf{k}}\left(t\right)},$$

$$C_{n\mathbf{k}}\left(t\right)={\int }_{-\infty }^{\infty }d\omega \frac{{\beta }_{n\mathbf{k}}\left(\omega \right)}{{\omega }^2}\left(e^{-i\omega t}+i\omega t-1\right),$$

$$\begin{array}{rl}{\sigma }_{\alpha \beta }(\omega )=& \frac{\pi \hslash e^2}{V_{\text{uc}}}\int d{\omega }^{\prime }\frac{f({\omega }^{\prime })-f({\omega }^{\prime }+\omega )}{\omega }\\ & \times \sum _{n\mathbf{k}}{v}_{n\mathbf{k}}^{\alpha }{v}_{n\mathbf{k}}^{\beta }{A}_{n\mathbf{k}}({\omega }^{\prime })A_{n\mathbf{k}}({\omega }^{\prime }+\omega ),\end{array}$$

$$n_c=\sum _{n\mathbf{k}}{\int }_{-\infty }^{\infty }d\omega f(\omega )A_{n\mathbf{k}}(\omega ).$$

$$\Phi (\omega )=\frac{\pi \hslash e^2}{V_{\text{uc}}}\sum _{n\mathbf{k}}{v}_{n\mathbf{k}}^{\alpha }{v}_{n\mathbf{k}}^{\beta }{\left|A(n\mathbf{k},\omega )\right|}^2\left(-\frac{\partial f(\omega )}{\partial \omega }\right),$$

$${\tau }_{\text{tr}}^{\ast }=-\frac{2}{\pi {\sigma }^{dc}}{\int }_0^{\infty }\frac{1}{{\omega }^{\prime }}\frac{\partial \sigma ({\omega }^{\prime })}{\partial {\omega }^{\prime }}d{\omega }^{\prime }.$$

$$C_{n\mathbf{k}}\left(t\right)=P{\int }_{-\infty }^{\infty }d\omega \frac{{\tilde{\beta }}_{n\mathbf{k}}\left(\omega \right)}{{\omega }^2}\left(e^{-i\omega t}-1\right)-it{\Sigma }_{n\mathbf{k}}({\epsilon }_{n\mathbf{k}}),$$

$$\begin{array}{r}{\int }_{-\infty }^{\infty }\frac{\left(e^{-i\omega t}+i\omega t-1\right)}{{\omega }^2}d\omega =-\pi t\\ P{\int }_{-\infty }^{\infty }\frac{{\tilde{\beta }}_{n\mathbf{k}}(\omega )}{\omega }d\omega =-\mathrm{Re}{\Sigma }_{n\mathbf{k}}\left({\epsilon }_{n\mathbf{k}}\right).\end{array}$$

$$C_{n\mathbf{k}}^s\left(t\right)=\sum _{{\omega }_l\ne 0}\frac{{\tilde{\beta }}_{n\mathbf{k}}({\omega }_l)}{{\omega }_l^2}\left(e^{-i{\omega }_{l}t}-1\right)\Delta \omega ,$$

$$G_{n\mathbf{k}}^R(t)=-i\theta (t)e^{-it({\epsilon }_{n\mathbf{k}}+{\Sigma }_{n\mathbf{k}}({\epsilon }_{n\mathbf{k}}))}{e}^{C_{n\mathbf{k}}^s\left(t\right)}.$$

$$A_{n\mathbf{k}}(\omega )=A_{n\mathbf{k}}^0(\omega )\ast A_{n\mathbf{k}}^s(\omega )/2\pi ,$$

$$A_{n\mathbf{k}}^0(\omega )=\frac{-\mathrm{Im}\Sigma ({\epsilon }_{n\mathbf{k}})/\pi }{(\omega -{\epsilon }_{n\mathbf{k}}-\mathrm{Re}\Sigma ({\epsilon }_{n\mathbf{k}}){)}^2+{\mathrm{Im}\Sigma ({\epsilon }_{n\mathbf{k}})}^2}$$

$${\Sigma }_{n\mathbf{k}}(\omega ,T)=\sum _{m,\nu \mathbf{q}}{\left|g_{mn,\nu }(\mathbf{k},\mathbf{q})\right|}^2\times$$

$$\left[\frac{N_{\nu \mathbf{q}}+f_{m\mathbf{k}+\mathbf{q}}}{\omega -{\epsilon }_{m\mathbf{k}+\mathbf{q}}+{\omega }_{\nu \mathbf{q}}+i\eta }+\frac{N_{\nu \mathbf{q}}+1-f_{m\mathbf{k}+\mathbf{q}}}{\omega -{\epsilon }_{m\mathbf{k}+\mathbf{q}}-{\omega }_{\nu \mathbf{q}}+i\eta }\right],$$

$$g_{mn,\nu }\left(\mathbf{k},\mathbf{q}\right)=\sqrt{\frac{\hslash }{2{\omega }_{\nu \mathbf{q}}}}\sum _{\kappa \alpha }\frac{{\mathbf{e}}_{\nu \mathbf{q}}^{\kappa \alpha }}{\sqrt{M_{\kappa }}}⟨m\mathbf{k}+\mathbf{q}\left|{\partial }_{\mathbf{q}\kappa \alpha }V\right|n\mathbf{k}⟩,$$

Full text

I introduction

Understanding charge transport in complex materials is a grand challenge of fundamental and technological relevance. The interactions between electrons and phonons (the quanta of lattice vibrations) set an intrinsic limit for the conductivity and typically control charge transport near room temperature. When electron-phonon ($e$-ph) interactions are weak, charge transport is well described by the scattering of quasiparticles (QPs) Ziman (2007), leading to the well-known band-like conduction regime. As the $e$-ph interactions become stronger, the electrons are dressed by a cloud of phonons, forming composite charge carriers known as polarons Fröhlich (1954); Holstein (1959); Devreese and Alexandrov (2009). In the limit of strong $e$-ph coupling, the electrons are self-trapped by the lattice distortions, and the conduction mechanism becomes the thermally activated hopping of localized polarons Emin (2012).

Many oxides and organic crystals exhibit $e$-ph coupling strengths intermediate between the band-like and polaron hopping limits Iadonisi et al. (2006); Fratini et al. (2016). In this so-called “large polaron” regime, the charge transport mechanisms and their temperature dependence are not well understood. The transition from band-like to hopping conduction for increasing $e$-ph coupling strength is also unclear, and recent work on the Holstein model uncovered an incoherent transport regime at intermediate coupling Mishchenko et al. (2015). Yet, predictive calculations and microscopic understanding of charge transport in the intermediate $e$-ph coupling regime remain an open challenge.

Strontium titanate (SrTiO3), which is stable in the cubic phase above 105 K, is a prototypical material with intermediate $e$-ph coupling in which large polaron effects are clearly seen in experiments van Mechelen et al. (2008); Devreese et al. (2010); Chen et al. (2015); Wang et al. (2016). Charge transport in SrTiO3 is a decades-old problem Frederikse and Hosler (1967); Wemple et al. (1969); Verma et al. (2014); Zhou et al. (2016); Cain et al. (2013); Lin et al. (2017), yet its underlying microscopic mechanisms are still debated Collignon et al. (2019). The electron mobility in cubic SrTiO3 exhibits a roughly $T^{-3}$ temperature dependence above 150 K Cain et al. (2013); Lin et al. (2017), which is commonly attributed to the scattering of electron QPs with phonons Frederikse and Hosler (1967); Wemple et al. (1969); Verma et al. (2014); Zhou et al. (2016). Different phenomenological models based on QP scattering have been proposed that can fit the experimental transport data Verma et al. (2014); Mikheev et al. (2015). However, the carrier mean free paths extracted from experiment in SrTiO3 fall below the interatomic distance Lin et al. (2017), violating the Mott-loffe-Regel (MIR) criterion for the applicability of the QP scattering picture Hussey et al. (2004). While there is consensus that large polarons are present in SrTiO3 van Mechelen et al. (2008); Devreese et al. (2010); Chen et al. (2015); Wang et al. (2016); Devreese (2016), charge transport in this regime cannot yet be predicted, and detailed microscopic understanding has remained elusive Collignon et al. (2019).

First-principles calculations based on lowest-order $e$-ph scattering plus the Boltzmann transport equation (BTE) Bernardi (2016) can accurately predict the conductivity in simple metals and semiconductors Mustafa et al. (2016); Zhou and Bernardi (2016); Liu et al. (2017); Ma et al. (2018); Lee et al. (2018). We have recently shown Zhou et al. (2018a) that when this approach is applied to SrTiO3, one can obtain an accurate temperature dependence for the electron mobility if all the phonons (including the soft modes) are taken into account. However, the absolute value of the computed electron mobility is an order of magnitude greater than experiment Zhou et al. (2018a). It is clear that QP scattering alone cannot explain charge transport in SrTiO3, consistent with the MIR limit violation; a fully quantum mechanical framework is needed to predict charge transport in this beyond-QP regime.

Here we show first-principles calculations of charge transport in cubic SrTiO3 using a finite-temperature retarded cumulant approach that includes higher-order $e$-ph interactions and goes beyond the QP scattering picture. Our calculations can accurately predict the experimental electron mobility in SrTiO3 between 150$-$300 K, and further shed light on its microscopic origin. We show that the weight of the QP peak in the electron spectral function is strongly renormalized, with significant weight transfer to the incoherent phonon satellites. While the renormalized QPs control transport at low temperature, the incoherent contributions from the phonon satellites and large momentum states are significant at room temperature, indicating a transport regime beyond the QP scattering paradigm. Consistent with these trends, our analysis shows that near room temperature the scattering rate (extracted from the optical conductivity) breaks the Planckian limit of $k_{\rm B}T$ for semiclassical transport Hartnoll (2015), clearly indicating a beyond-QP charge transport regime. Our work opens new avenues for computing charge transport in complex materials with large polarons and beyond the QP scattering regime.

II Electron spectral function

Central to our approach for computing charge transport is the electron spectral function, $A(\omega)$, which can be seen as the density of states of a single electron. Due to the interactions with the phonons, the spectral function consists of a QP peak representing a single-electron-like excitation and an incoherent part including both phonon satellite peaks and a background contribution Damascelli et al. (2003); Martin et al. (2016). While the spectral weights of the QP peak and incoherent part may vary, they always add up to one due to the sum rule $\int\!d\omega A(\omega)=1$, which amounts to a conservation of the electron. To investigate the dynamics of the electrons and their interactions with the phonons, we compute the spectral function with a finite-temperature retarded cumulant approach that can account for higher-order $e$-ph interactions (see Appendix A), and use it directly to predict transport, without relying on QP scattering approaches.

We focus on the electron spectral function in cubic SrTiO3 above 110 K. Leveraging our recently developed approach Zhou et al. (2018a), which combines density functional theory (DFT), its linear response extension Baroni et al. (2001) and the temperature dependent effective potential (TDEP) method Hellman et al. (2011), we compute the band structure, lattice vibrations and $e$-ph interactions for all phonon modes, including the soft modes due to the lattice anharmonicity (see Appendix B). Using these quantities, we calculate the spectral function using a retarded cumulant formalism Kas et al. (2014); Kas and Rehr (2017). The retarded cumulant approach is based on an exponential ansatz for the retarded Green’s function of an electronic state $\left|n\bm{k}\right\rangle$ (with band index $n$ and momentum $\bm{k}$) in the time domain:

[TABLE]

where $t$ is time, $\theta(t)$ is the Heaviside step function, and $\varepsilon_{n\bm{k}}$ is the non-interacting electron energy from DFT. The $e$-ph interactions are included in the cumulant $C_{n\bm{k}}\left(t\right)$; an approximate expression derived in Ref. Kas et al. (2014) computes $C_{n\bm{k}}\left(t\right)$ using the off-shell lowest-order $e$-ph self-energy, $\Sigma_{n\bm{k}}(\omega)$ Bernardi (2016) as input:

[TABLE]

where $\beta_{n\bm{k}}(\omega)\equiv|\mathrm{Im}\Sigma_{n\bm{k}}\left(\omega+\varepsilon_{n\bm{k}}\right)|/\pi$. The spectral function can be obtained for each state from the retarded Green’s function in the frequency domain, using $A_{n\bm{k}}(\omega)=-\mathrm{Im}G_{n\bm{k}}^{R}(\omega)/\pi$. Our scheme to compute the retarded Green’s function at finite-temperature (see Appendix A) further allows us to compute the spectral function as a function of temperature. The retarded cumulant approach includes higher-order $e$-ph Feynman diagrams beyond the Migdal approximation Hedin (1999), and it produces accurate spectral functions (see below) that can capture the strong $e$-ph interactions. On the other hand, we have verified that the lowest-order Dyson-Migdal approximation generates spectral functions with large errors in the QP spectral weight and satellites energies, consistent with recent results at zero temperature Nery et al. (2018).

The computed electron spectral functions for the three lowest conduction bands in cubic SrTiO3 at 110 K are combined in a color map and shown in Fig. 1(a). Each state exhibits a rather sharp QP peak at low energy and broader phonon satellite peaks at higher energies ($\sim$60 meV or more above the QP peak). By tracking the low-energy QP peaks, we map the energy-momentum dispersion of the QPs. Figure 1(b) shows that the interacting QPs exhibit a heavier effective mass than in the DFT band structure calculations, in which the $e$-ph interactions are not included. The mass enhancement is a factor of 1.8$-$2.6 for different bands and directions, and it increases only slightly with temperature. Taking the lowest bands along $\Gamma-$M as an example, the DFT effective mass is roughly 0.75$m_{e}$ ($m_{e}$ is the electron mass), compared to a QP effective mass of 1.4$m_{e}$ at 110 K and a slightly heavier mass of 1.6$m_{e}$ at 300 K. The mass enhancement is thus roughly a factor of 2, in excellent agreement with experimental results at light doping van Mechelen et al. (2008); Wang et al. (2016).

The interactions with the surrounding phonons not only make the electron QPs heavier, they also significantly reduce the QP spectral weight to a value well below one, transferring weight to the higher-energy incoherent phonon satellites. The QP peak even disappears at large electron momenta, leaving a spectral function made up entirely by the incoherent background. These trends are analyzed in detail below and in Fig. 2, focusing on how the spectral function of the lowest conduction band changes as a function of temperature and momentum.

The spectral function at the conduction band minimum at $\Gamma$ exhibits a main QP peak, two main satellites (also known as phonon sidebands or replicas), and weaker additional satellites at higher energy [see Fig. 2(a)]. The two main satellite peaks are at an energy $\omega_{\text{LO-1}}=98$ meV and $\omega_{\text{LO-2}}=57$ meV above the main QP peak; these values correspond, respectively, to the energies of the two longitudinal optical (LO) modes with long-range $e$-ph interactions that exhibit the strongest coupling with electrons Zhou et al. (2018a). Note also the presence in Fig. 2(a) of weak phonon sideband peaks at energies of $2\omega_{\text{LO-1}}$ and $\omega_{\text{LO-1}}+\omega_{\text{LO-2}}$. These higher-order replicas, which are known to occur in the strong coupling limit of the Holstein model Berciu (2006), are akin to the higher harmonics observed in phonon Floquet states Hübener et al. (2018); they are a signature of strong coupling with the LO modes.

The main phonon sidebands are associated with the polaron plus one-phonon continuum, and are a hallmark of the large polaron regime Mahan (2000). Note that our calculations are performed on lightly $n$-doped SrTiO3, with the chemical potential lying below the lowest QP peak; therefore, as expected, the phonon satellites appear at energy higher than the QP peak as they correspond to the excitation of an “electron-like” QP plus one LO phonon Damascelli et al. (2003); Nery et al. (2018). On the other hand, recent angle-resolved photoemission measurements on heavily $n$-doped samples, in which the chemical potential is above the QP peak, revealed a phonon sideband roughly $\sim$100 meV below the QP peak, which corresponds to the excitation of a “hole-like” QP plus one LO phonon Chen et al. (2015); Wang et al. (2016); Damascelli et al. (2003). The energy difference between the QP peak and phonon sideband observed in experiments is in very good agreement with the $\omega_{\text{LO-1}}=98$ meV energy difference between the QP peak and the most intense satellite peak found in our computed spectral function. We have verified that at high doping our approach also gives satellites with energies lower than the QP peak, consistent with experiment. The fact that the satellite position can be higher or lower than the QP peak depending on doping is well known Damascelli et al. (2003); Zhou et al. (2018b), but since inverse photoemission measurements are challenging, the low doping regime we compute here is rarely probed experimentally.

To compute the spectral weights of the QP peak and incoherent part, we integrate the spectral function up to an energy $\omega$, obtaining the spectral weight $N(\omega)=\int_{-\infty}^{\omega}A(\omega^{\prime})d\omega^{\prime}$ given in Fig. 2(b). We find that the spectral weight of the QP peak is $\sim$0.4 at 110 K, and thus much less than the unit value of the weak $e$-ph interaction limit. As the temperature increases from 110 K to 300 K, both the QP peak and the phonon satellites are broadened and smeared out [see Fig. 2(a)], but the QP spectral weight changes only slightly, primarily due to an overlap between the QP peak and the phonon satellites near 300 K. At all temperatures between 110$-$300 K, the QP weight is strongly renormalized to a value of $\sim$0.4, implying that (pictorially) only half of the electron resides in the QP state, while the other half contributes to the incoherent dynamical interactions with the phonons.

Figure  2(c) reveals the disappearance of the QP peak at large enough momentum $\bm{k}$ by showing how the spectral function changes as we increase $\bm{k}$ along the \Gamma$$-$$M Brillouin zone line. We find that the QP spectral weight decreases with increasing momentum [see Fig. 2(d)] and that the QP peak ultimately disappears at \bm{k}=0.12$$M, leading to a fully incoherent spectral function at larger momenta. These so-called end points of the QP peak, which have been predicted in both the Fröhlich Mishchenko et al. (2000) and Holstein models Goodvin and Berciu (2010), are yet another signature of the strong $e$-ph interactions. The decrease of the QP spectral weight and the disappearance of the QP peak at large momentum have a significant impact on transport, as we discuss below.

III Electron mobility

Large polaron transport is commonly believed to be the band-like conduction of QPs with enhanced effective mass. However, this simplified picture neglects the fact that the QP weight can drop to values much smaller than one (here, to roughly 0.4, as discussed above), and that the contribution to transport from the incoherent part of the spectral function can be significant. Temperature also plays a primary role. At low temperature, the electrons occupy the low-energy QP states, and there are only few LO phonons due to their relatively high energy. As the temperature increases, thermal fluctuations push the electrons to higher energies, exciting the electrons outside of the QP peak into the incoherent regime. In addition, the number of LO phonons grows rapidly with temperature, leading to strong dynamical interactions between the electron and its phonon cloud. The incoherent contributions are thus expected to significantly influence transport at higher temperatures.

To investigate these points quantitatively in SrTiO3, we compute the conductivity directly from the spectral functions $-$ therefore including both the QP and incoherent contributions $-$ using the Kubo formula Mahan (2000). In the absence of current-vertex corrections, the conductivity can be expressed as Economou (2006); Basov et al. (2011)

[TABLE]

where $\bm{v}_{n\bm{k}}$ is the band velocity of the electronic state $\left|n\bm{k}\right\rangle$, $f(\omega)$ is the Fermi-Dirac distribution, $V_{\text{uc}}$ is the unit cell volume, and $\alpha$ and $\beta$ are Cartesian directions. We also compute the DC conductivity, using $\sigma^{dc}=\sigma(\omega\rightarrow 0)$, and the electron mobility as $\mu=\sigma^{dc}/n_{c}e$, where the carrier concentration $n_{c}$ is computed as

[TABLE]

We perform mobility calculations in the lightly doped regime, with electron concentrations ranging between 10^{17}$$-$$10^{18} cm*-3*, and find that the computed mobility is nearly independent of the chosen concentration above 150 K, consistent with experimental data Cain et al. (2013); Lin et al. (2017).

Our computed electron mobility in SrTiO3 as a function of temperature is shown in Fig. 3(a) and compared with experimental data taken from Refs. Cain et al. (2013); Lin et al. (2017). Both the absolute value and the temperature dependence of the computed mobility are in excellent agreement with experiments. Our computed mobility at 300 K is about 8 cm2/V s, versus an experimental value of $\sim$5 cm2/V s Cain et al. (2013); Lin et al. (2017). For comparison, the mobility obtained using lowest-order $e$-ph scattering plus the BTE is an order of magnitude higher than experiments Zhou et al. (2018a), as also shown in Fig. 3(a).

The order-of-magnitude mobility drop from the BTE to the cumulant approach in SrTiO3 is due to both the QP renormalization and the incoherent contributions. The cumulant approach clearly shows that the QP peak is strongly renormalized with significant weight transfer to higher energy incoherent satellites, and the QP peak even disappears as momentum $\bm{k}$ goes beyond the end points [see Fig. 2 and Fig. 3(b,c)]. However, the BTE inherently assumes that the QP states are well defined for all bands and momentum values, so it fails to capture the contributions from the incoherent part of the spectral function, thus placing all the weight in the QP peak and significantly overestimating the mobility. Our results show unambiguously that the established lowest-order $e$-ph plus BTE approach Bernardi (2016); Mustafa et al. (2016); Zhou and Bernardi (2016); Liu et al. (2017); Ma et al. (2018); Lee et al. (2018) is accurate only in the case of weak $e$-ph interactions. Including higher-order $e$-ph processes and incoherent polaron effects via the cumulant approach greatly improves the computed mobility in materials with strong $e$-ph coupling.

Correctly taking into account the contributions from all the phonon modes, including the soft modes Zhou et al. (2018a), is also essential for predicting the electron mobility and its temperature dependence in SrTiO3. Neglecting the soft mode in the calculations leads to significant errors in the computed mobility and its temperature dependence (see Fig. 5 in Appendix C), both in the BTE and in the cumulant approach. These results show that the soft modes, and not just the LO phonons as is widely believed, also contribute to charge transport in the large polaron regime.

Lastly, note also that in the weak $e$-ph coupling limit, in which the spectral function consists only of a sharp QP peak and no incoherent contributions, the expression we use for $\sigma^{dc}$ reduces to the conductivity obtained from the relaxation time approximation (RTA) of the BTE Economou (2006); Allen (2015). For a material with weak $e$-ph interactions and negligible polaron effects, one thus expects that the cumulant and BTE approaches predict the same transport results. We perform this sanity check for GaAs, showing that its mobility curves computed with the cumulant and BTE-RTA are in close agreement (see Fig. 6 in Appendix D). In the case of SrTiO3, one the other hand, it is apparent that the BTE cannot predict the mobility correctly.

IV Coherent and incoherent contributions

Our approach for computing the conductivity in the presence of polarons further allows us to resolve the coherent and incoherent contributions to transport and to uncover different transport regimes as a function of temperature. We write the DC conductivity as $\sigma^{dc}=\int\Phi(\omega)d\omega$, where the integrand

[TABLE]

is referred to as the transport distribution function (TDF), and is employed here to quantify the contributions to the DC conductivity as a function of electron energy $\omega$. We analyze the TDF in the 150$-$300 K range in Fig. 3(b,c) , and find that both the coherent QP peak and the incoherent part contribute to transport. At temperatures below $\sim$200 K, the TDF spans primarily the QP peaks [see Fig. 3(b)], implying that transport is dominated by the scattering of QPs with a strongly renormalized spectral weight. As the temperature increases, the high-energy tail of the TDF extends over the incoherent contributions [see Fig. 3(c)], which become important and contribute by as much as 40% to the conductivity at 300 K, as discussed below. Therefore, transport near room temperature is governed not only by the weight-renormalized QPs, but also by the incoherent phonon satellites above the QP peaks and by the polaron states at large momenta beyond the end points, where the QP peaks disappear. The picture of QP scattering is inadequate to describe transport at room temperature in SrTiO3, and a more complex picture emerges in which transport is an interplay between the QP renormalization and the contributions from the incoherent phonon sidebands and from the polaron states beyond the end points, all of which are consequences of the dynamical interactions between the electrons and the phonon cloud.

These conclusions are supported by experiments on the optical conductivity. Recent experiments in SrTiO3 show that the Drude peak at low frequency in the optical conductivity, which is associated with the coherent band-like transport of QPs, loses weight for increasing temperatures up to 300 K van Mechelen et al. (2008); Devreese et al. (2010). Figure 4(a) compares the low-energy optical conductivities at 150 K and 300 K, both computed with Eq. (3) and normalized to possess the same integral, consistent with the optical sum rule. The optical conductivities exhibit a Drude-like peak centered at zero frequency, and an incoherent shoulder structure consisting of phonon sidebands plus a broad background. We find a significant weight transfer from the Drude peak to the incoherent shoulder as the temperature increases from 150 K to 300 K, in agreement with experiments van Mechelen et al. (2008). The Drude peak is sharp at 150 K, but it broadens rapidly as the temperature increases [see Fig. 4(b)]. These trends confirm the transition seen in our transport results from a renormalized QP regime at low temperature to an incoherent, beyond-QP regime near room temperature.

We also extract an effective transport relaxation time, $\tau_{\text{tr}}^{*}$, from the optical conductivity through the extended Drude analysis of Ref. Deng et al. (2014):

[TABLE]

Figure 4(c) shows the inverse of $\tau_{\text{tr}}^{*}$, namely the effective scattering rate characterizing the width of the Drude peak. We find that this effective scattering rate increases rapidly with temperature, reaching values much greater than the QP scattering rate extracted from the QP peak of the spectral function in Fig. 2(a). Due to the uncertainty principle, in a semiclassical transport regime the scattering rate cannot exceed the so-called Planckian limit of $k_{B}T$ Hartnoll (2015). We find that the effective scattering rate in SrTiO3 exceeds the Planckian limit $k_{B}T$ above 250 K, highlighting the beyond-quasiparticle nature of charge transport in SrTiO3 near room temperature. The breaking of the Planckian limit in SrTiO3 at room temperature is consistent with very recent results obtained by Mishchenko et al. using a model Hamiltonian Mishchenko et al. (2019).

Finally, one expects the incoherent contributions to transport to be significant at temperatures where the scattering rate exceeds the Planckian limit, namely above $\sim$250 K in our calculations. We quantity the incoherent contribution to transport by defining the conductivity ratio $\sigma^{dc}_{\rm inc}/\sigma^{dc}$, where the incoherent contribution to the conductivity $\sigma^{dc}_{\rm inc}$ is obtained by integrating the TDF at energies greater than all QP peaks ($\omega>32$ meV in our calculations). Figure 2(d) shows that above 250 K the beyond-QP incoherent contributions amount to up to $\sim$40% of the total DC conductivity, confirming that the Planckian limit breaking is associated with a fully quantum mechanical transport regime beyond the QP scattering paradigm. While the Planckian limit breaking has been typically associated with strange metals and other strongly correlated phases of matter Emery and Kivelson (1995); Gunnarsson et al. (2003); Bruin et al. (2013), our results highlight that strong $e$-ph interactions can also lead to this quantum mechanical transport regime.

V Conclusion

In summary, we developed a broadly applicable approach for computing charge transport in the large polaron regime in materials with intermediate $e$-ph coupling strength. Our calculations on SrTiO3 unveil a transition from band-like transport of strongly weight-renormalized QPs at low temperature to an incoherent transport regime beyond the QP picture near room temperature. Our approach can shed new light on broad classes of materials with polaron effects, ranging from perovskites Miyata et al. (2017) and transition metal oxides Moser et al. (2013); Cancellieri et al. (2016) to high-$T_{c}$ superconductors Lanzara et al. (2001); Legros et al. (2019).

Acknowledgements.

J.-J.Z. has benefited from discussion with N.-E. Lee. This work was supported by the Joint Center for Artificial Photosynthesis, a DOE Energy Innovation Hub, supported through the Office of Science of the U.S. Department of Energy under Award No. DE-SC0004993. M.B. acknowledges support by the National Science Foundation under Grant No. ACI-1642443, which provided for code development, and Grant No. CAREER-1750613, which provided for theory and method development. This work was partially supported by the Air Force Office of Scientific Research through the Young Investigator Program, Grant FA9550-18-1-0280. This research used resources of the National Energy Research Scientific Computing Center, a DOE Office of Science User Facility supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231.

Appendix A Finite temperature retarded cumulant approach

We implement and employ a retarded cumulant approach Kas et al. (2014); Kas and Rehr (2017) to compute the electron spectral function at finite temperature. The retarded cumulant approach [see Eq. (1) and (2)] was recently employed to study the electron spectral function near the band edge in insulators at zero temperature Nery et al. (2018). However, the separation scheme used in Ref. Nery et al. (2018) to compute the cumulant and the spectral function is not well behaved at finite temperature Kas et al. (2014), where $\beta_{n\bm{k}}(0)>0$. Here, we develop a new scheme to compute the retarded Green’s function at finite temperature. We rewrite Eq. (2) as

[TABLE]

where $\tilde{\beta}_{n\bm{k}}\left(\omega\right)\equiv\beta_{n\bm{k}}(\omega)-\beta_{n\bm{k}}(0)$ and $P$ denotes the Cauchy principal value of the integral. The following two relations are used to derive Eq. (7) :

[TABLE]

We evaluate $\tilde{\beta}_{n\bm{k}}\left(\omega\right)$ on a discrete frequency grid, $\omega_{l}=l\Delta\omega$ with $l$ an integer. We label the first term in Eq. (7) as $C_{n\bm{k}}^{s}\left(t\right)$ and compute it from $\tilde{\beta}_{n\bm{k}}\left(\omega_{l}\right)$,

[TABLE]

Substituting Eq. (7) and (8) into Eq. (1), the retarded Green’s function becomes

[TABLE]

We then evaluate the spectral function by Fourier transforming to the frequency domain, and obtain

[TABLE]

where $*$ denotes a convolution operation; the two quantities entering this formula are

[TABLE]

and $A^{s}_{n\bm{k}}(\omega)$, which is the Fourier transform of $e^{C^{s}_{n\bm{k}}\left(t\right)}$.

Note that $C^{s}_{n\bm{k}}\left(t\right)$ defined in Eq. (8) is periodic with a period of $\mathcal{T}=2\pi/\Delta\omega$; we compute $e^{C^{s}_{n\bm{k}}\left(t\right)}$ on a discrete time grid $t_{j}=j\Delta t$ with $j\in[-N,N]$ and $\Delta t=\mathcal{T}/(2N+1)$. It is seen from Eq. (8) that $C_{n\bm{k}}^{s}\left(-t_{j}\right)=C_{n\bm{k}}^{s}\left(t_{j}\right)^{*}$ and $C_{n\bm{k}}^{s}\left(0\right)=0$, which allows us to compute $C_{n\bm{k}}^{s}\left(t_{j}\right)$ only for $j\in[1,N]$ and obtain $e^{C^{s}_{n\bm{k}}\left(t_{j}\right)}$ on the full time grid. We then obtain $A^{s}_{n\bm{k}}(\omega)$ on a frequency grid with the same step $\Delta\omega$ as $\tilde{\beta}_{n\bm{k}}\left(\omega\right)$ via a discrete Fourier transform of $e^{C^{s}_{n\bm{k}}\left(t_{j}\right)}$. The value of $A^{s}_{n\bm{k}}(\omega)$ is guaranteed to be real. We have carefully checked the convergence of $A^{s}_{n\bm{k}}(\omega)$ with respect to the range and step size of both the time and frequency grids. Specifically, we employed $\omega$ values ranging from $-0.8$ to $0.8$ eV with a small step $\Delta\omega$ of 0.05 meV and $N=32000$ for the time grid in our calculations.

Appendix B Computational details

We carry out ab initio calculations on cubic SrTiO3 with a lattice parameter of 3.9 Å. The ground-state electronic structure is computed within the Perdew-Burke-Ernzerhof generalized gradient approximation Perdew et al. (2008) of density functional theory (DFT) using the Quantum ESPRESSO code Giannozzi et al. (2009). We employ fully-relativistic norm-conserving pseudopotentials, which include spin-orbit coupling (SOC) Hamann (2013); van Setten et al. (2018), together with a plane-wave basis set with a kinetic energy cutoff of 85 Ry. The anharmonic lattice-dynamical properties of cubic SrTiO3 are computed using the temperature-dependent effective potential (TDEP) method Hellman et al. (2011) with 4 $\times$ 4 $\times$ 4 supercells, together with a scheme we recently developed Zhou et al. (2018a) to accurately include the long-range dipole-dipole contributions in the interatomic force constants.

We use our in-house developed perturbo code 111The code employed in this work will be released in the future at http://perturbo.caltech.edu to compute the lowest-order $e$-ph self-energy at temperature $T$ for the Bloch state $\left|n\bm{k}\right\rangle$ with band $n$ and crystal momentum $\bm{k}$,

[TABLE]

where $\varepsilon_{n\bm{k}}$ is the DFT band electron energy, $\omega_{\nu\bm{q}}$ is the energy of a phonon with branch index $\nu$ and wave vector $\bm{q}$, $f_{n\bm{k}}$ and $N_{\nu\bm{q}}$ are, respectively, the Fermi and Bose occupation numbers evaluated at temperature $T$, and $\eta$ is a small broadening. The key quantities are the $e$-ph coupling matrix elements, $g_{mn,\nu}(\bm{k},\bm{q})$, defined as

[TABLE]

where $\partial_{\bm{q}\kappa\alpha}V\equiv\sum_{p}e^{i\bm{q}\bm{R}_{p}}\partial_{p\kappa\alpha}V$, with $\partial_{p\kappa\alpha}V$ the variation of the Kohn-Sham potential for a unit displacement of the atom $\kappa$ (with mass $M_{\kappa}$ and in the unit cell at $\bm{R}_{p}$) in the $\alpha$-direction, and $\bm{e}_{\nu\bm{q}}$ is the phonon displacement eigenvector. We obtain the $e$-ph matrix elements with the method described in our recent work Zhou et al. (2018a). Briefly, we first compute the electronic wave functions $\left|n\bm{k}\right\rangle$ on an 8 $\times$ 8 $\times$ 8 $\bm{k}$-point grid using DFT, construct Wannier functions for the Ti-$t_{2g}$ orbitals from the three lowest conduction bands using the Wannier90 code Mostofi et al. (2014), compute $\partial_{\bm{q}\kappa\alpha}V$ on an 8 $\times$ 8 $\times$ 8 $\bm{q}$-point grid with density functional perturbation theory Baroni et al. (2001), and compute anharmonic phonon energies $\omega_{\nu\bm{q}}(T)$ and eigenvectors $\bm{e}_{\nu\bm{q}}(T)$ with TDEP. Using these ingredients, we first evaluate the $e$-ph matrix elements in Eq. (13) on coarse grids. Wannier interpolation together with a long-range $e$-ph correction for the polar modes Giustino et al. (2007); Sjakste et al. (2015); Verdi and Giustino (2015) is then employed to interpolate $g_{mn,\nu}(\bm{k}\bm{q})$ for $\bm{k}$- and $\bm{q}$-points on very fine Brillouin zone grids, which are needed to converge the $e$-ph self-energy in Eq. (12). An approach we recently developed Zhou and Bernardi (2016) is employed to converge the imaginary part of the self-energy, $\mathrm{Im}\Sigma_{n\bm{k}}(\omega)$, which is computed off-shell on a fine energy $\omega$ grid, and the real part of the self-energy, $\mathrm{Re}\Sigma_{n\bm{k}}$, which is evaluated on-shell at the electron energy $\varepsilon_{n\bm{k}}$. The spectral functions $A_{n\bm{k}}(\omega)$ are then obtained using the finite temperature retarded cumulant approach described in Appendix A, using $\mathrm{Im}\Sigma_{n\bm{k}}(\omega)$ and $\mathrm{Re}\Sigma_{n\bm{k}}(\varepsilon_{n\bm{k}})$ as input.

We compute spectral functions $A_{n\bm{k}}(\omega)$ on fine $\bm{k}$-point grids with up to 1003 points in the Brillouin zone, which are then employed in the conductivity and mobility calculations using Eq. (3) and (4). Note that the spectral functions depend on both the temperature $T$ and chemical potential $\mu$. For each temperature and carrier concentration $n_{c}$ considered in our calculations, we determine self-consistently the chemical potential $\mu(T,n_{c})$ using Eq. (4); this means that the spectral functions computed using the chemical potential $\mu(T,n_{c})$ give a consistent carrier concentration $n_{c}$.

Appendix C Soft mode contribution to transport

To highlight the contribution from the soft modes to charge transport in SrTiO3, we compare the electron mobility computed with and without the soft mode contribution, as shown in Fig. 5. Our recent calculations using the BTE show that the ferroelectric soft mode plays an important role in transport, especially at low temperatures below 250 K, and that including the soft mode contribution is critical to obtaining an accurate temperature dependence of the electron mobility Zhou et al. (2018a). The results from the cumulant approach exhibit the same trend $-$ neglecting the soft mode contribution leads to a significant overestimate of the mobility below 250 K and to an inaccurate temperature dependence. Although the phonon satellites in the spectral function is primarily due to the strong coupling between the electron and LO phonons, the ferroelectric soft mode can impact the width of the QP peak, thus contributing to charge transport.

Appendix D Computed electron mobility in GaAs: BTE vs cumulant approach

The method we developed in this work, which uses the retarded cumulant approach plus the Kubo formula to compute charge transport, is general and can be applied broadly in materials with $e$-ph interactions ranging from weak to strong. To illustrate this point, we apply our approach to GaAs, a high mobility semiconductor in which large-polaron effects are very weak. We have recently computed the temperature dependent electron mobility in GaAs using the BTE approach within the RTA Zhou and Bernardi (2016). With the same computational settings for the DFT and DFPT calculations as in Ref. Zhou and Bernardi (2016), we compute the electron mobility at temperatures between 250$-$400 K using the cumulant approach. Figure 6(a) compares the electron mobility from the cumulant and BTE-RTA calculations with experimental data. Due to the weak $e$-ph coupling, the spectral function in GaAs consists of a sharp QP peak with near-unit spectral weight [see Fig. 6(b)], plus a weak phonon sideband at energy $\omega_{\rm LO}=36$ meV above the QP peak, where $\omega_{\rm LO}$ is the LO phonon energy in GaAs. The cumulant approach and the BTE-RTA give mobility results in close agreement with each other, as is expected in materials with weak $e$-ph coupling. We have thus shown that our implementation of the cumulant plus Kubo approach without current-vertex corrections reduces to the RTA solution of the BTE in the weak $e$-ph coupling limit Economou (2006); Allen (2015).

For linear-response theory to be consistent with the BTE, the Kubo formula with current-vertex corrections should reduce in the weak $e$-ph coupling limit to the full solution of the BTE, e.g. the solution obtained in practice with an iterative method Ma et al. (2018). Our calculation on GaAs in the absence of current-vertex corrections shows good agreement with the BTE-RTA, and similarly we expect that including the current-vertex corrections would give results consistent with the iterative BTE solution for GaAs. In the case of SrTiO3, the BTE-RTA and the iterative solution of BTE produce very similar results as we have shown previously Zhou et al. (2018a). Therefore, the current-vertex corrections are not essential for SrTiO3 and are neglected in our calculations.