Physical Review B 93: Transport Coefficients of Graphene: Interplay of Impurity Scattering, Coulomb Interaction, and Optical Phonons

Hong-Yi Xie¹,* and Matthew S. Foster¹,²

¹Department of Physics and Astronomy, Rice University, Houston, Texas 77005, USA

²Rice Center for Quantum Materials, Rice University, Houston, Texas 77005, USA

(Received 25 January 2016; revised manuscript received 16 April 2016; published 2 May 2016)

We study the electric and thermal transport of the Dirac carriers in monolayer graphene using the Boltzmann equation approach. Motivated by recent thermopower measurements , we consider the effects of quenched disorder, Coulomb interactions, and electron–optical-phonon scattering. Via an unbiased numerical solution to the Boltzmann equation, we calculate the electrical conductivity, thermopower, and electronic component of the thermal conductivity, and discuss the validity of Mott's formula and of the Wiedemann-Franz law. An analytical solution for the disorder-only case shows that screened Coulomb impurity scattering, although elastic, violates the Wiedemann-Franz law even at low temperature. For the combination of carrier-carrier Coulomb and short-ranged impurity scattering, we observe the crossover from the interaction-limited (hydrodynamic) regime to the disorder-limited (Fermi-liquid) regime. In the former, the thermopower and the thermal conductivity follow the results anticipated by the relativistic hydrodynamic theory. On the other hand, we find that optical phonons become non-negligible at relatively low temperatures and that the induced electron thermopower violates Mott's formula. Combining all of these scattering mechanisms, we obtain the thermopower that quantitatively coincides with the experimental data. DOI: 10.1103/PhysRevB.93.195103

I. Introduction

Electric and thermal transports in monolayer graphene are influenced by various scattering mechanisms, such as quenched impurities, interparticle interactions, and phonons . In weakly disordered graphene, interaction effects can become dominant at nonzero temperature; in the vicinity of zero doping, the Coulomb-interacting massless Dirac carriers form a relativistic electron-hole plasma. In this interaction-limited regime, hydrodynamic theory predicts intriguing non-Fermi-liquid transport properties. First, the electron-hole fluid exhibits a finite and nonvanishing dc electrical conductivity at the Dirac point even in the absence of impurities, entirely due to inelastic electron-hole collisions.

Moreover, Mott's formula and the Wiedemann-Franz law are violated. In a Fermi liquid, these respectively determine the thermoelectric power and the electronic component of the thermal conductivity from the electrical conductivity. Instead, for graphene in the hydrodynamic regime, the thermopower at nonzero doping approaches the thermodynamic entropy per charge, and the thermal conductivity at the Dirac point diverges as the impurity concentration vanishes.

The theory predicts upper bounds for the thermoelectric power and electronic component of the thermal conductivity, limited only by disorder. While violating classical relations between thermoelectric coefficients, the latter are strongly constrained and interrelated by the relativistic hydrodynamics . We emphasize that this violation of the Mott and Wiedemann-Franz relations is different from the usual physics of narrow-gap/gapless semiconductors, for example, the bipolar diffusion process , where separated electron and hole currents are assumed. Strong inelastic electron-hole scattering in ultraclean graphene implies that a composite electron-hole fluid emerges , which cannot be decomposed into valence and conduction band components.

Three very recent experiments have provided substantial evidence for interaction-limited transport in graphene. Measurements of the electronic component of the thermal conductivity near charge neutrality showed large violations of the Wiedemann-Franz law . Nonlocal transport in doped graphene was used to probe the viscosity of the electron fluid . Finally, thermoelectric power measurements showed a substantial deviation from the Mott formula.

In this work, we model the experiment in Ref. using the Boltzmann equation to incorporate carrier-impurity, carrier-carrier, and carrier-optical phonon scattering mechanisms. The thermopower measurements in Ref. were performed on high-mobility graphene encapsulated by hexagonal-boron-nitride. The experiment was done over a large span of dopings, with charge-carrier density n ≡ ρ/(−e) ranging from zero to ± 3.0 × 10¹² cm⁻² (ρ is the charge density and e > 0 is the elementary charge). The measurements were performed at relatively high temperatures (130 K < T < 350 K) in order to fulfill the nondegenerate condition kBT ≫ μ over much of the doping span, while avoiding the electron-hole puddle regime at low temperatures near charge neutrality . Here, μ denotes the chemical potential, determined by the temperature and the fixed charge-carrier density n. In this regime, the measured thermopower is consistently larger than that predicted by Mott's formula , but saturates below the ideal hydrodynamic prediction. This feature suggests that in order to quantitatively characterize the thermoelectric transport in graphene, one should consider additional inelastic scattering mechanisms.

We exclude acoustic phonons, since at low doping the electron–acoustic-phonon scattering is quasielastic and incapable of producing large violations of Mott's formula. As discussed in Ref. , the Bloch-Grüneisen temperature TBG ≡ 2-vakF/kB plays the key role, where va is the acoustic phonon velocity and kF the Fermi wave vector. Assuming the acoustic phonon velocity equal to 2.6 × 10⁶ cm/s, one can estimate the Bloch-Grüneisen temperature as TBG ≈ 70√n K, where the density n is measured in units of 10¹² cm⁻². The experiment in Ref. is performed in the regime T ≫ TBG where the acoustic-phonon scattering is quasielastic. In addition, we disregard the effects of external magnetic fields or spin-flip mechanisms.

In this paper, we consider the inelastic optical-phonon scattering and model graphene by the Hamiltonian:

H = H₀ + Hₒₚₕ + Vᵢₘₚ⁽ˢ⁾ + Vᵢₘₚ⁽ˡ⁾ + Vᵢₙₜ + Vₑ₋ₒₚₕ, (1.1)

where H₀ describes the free Dirac fermions, Hₒₚₕ the optical-phonon bath, Vᵢₘₚ⁽ˢ⁾ and Vᵢₘₚ⁽ˡ⁾ the quenched short-ranged and long-ranged (Coulomb impurity) disorder potentials, respectively, Vᵢₙₜ the interparticle Coulomb interactions, and Vₑ₋ₒₚₕ the electron–optical-phonon coupling. We assume that both the time-reversal symmetry and spin SU(2) rotation symmetry are preserved in the presence of disorder and interactions. We also presume that the particle-hole symmetry as well as the honeycomb lattice space group symmetries (translations, rotations, and reflections) are preserved under disorder average . Concretely, each term in the Hamiltonian (1.1) is constructed as follows.

The short-ranged impurity Hamiltonian Vᵢₘₚ⁽ˢ⁾ takes the form as introduced in Ref. , which incorporates five types of time-reversal-symmetric disorder, all assumed to be zero-mean, short-ranged, and Gaussian-correlated. Five independent parameters {gu, gA, gA3, gm, gv} characterize their statistical fluctuations. In the Boltzmann equation, these parameters appear effectively in certain combinations [*G₀, f, b in Eq. (3.1)]. The term Vᵢₘₚ⁽ˡ⁾ gives the scalar potential due to Coulomb impurities. These are subject to the temperature and density-dependent static screening by the electron-hole plasma .

Coulomb interactions between carriers are encoded in Vᵢₙₜ. Dynamical screening is treated within the random phase approximation . Screening is crucial both in the low-temperature degenerate Fermi-liquid phase, but also in the high-temperature, nondegenerate regime of primary interest here. Different from a single component plasma, graphene ultimately screens better at higher temperatures, due to the proliferation of thermally-excited electron-hole pairs. At intermediate temperatures and finite charge density, the Thomas-Fermi length reaches a maximum. The interaction strength is encoded in the fine structure constant αᵢₙₜ that depends on the dielectric constants of the substrates . In the kinetic theory, dynamical screening suppresses the "collinear" singularity of the Coulomb collision integral, which is due to the linear dispersion of Dirac fermions and the energy conservation (see Appendix A). Note that for simplicity we only consider two-body collision processes that preserve the population of electrons and holes separately. We leave the effects of (three-body or impurity-assisted) electron-hole Auger imbalance relaxation processes to future study.

Three types of in-plane optical-phonon modes allowed by time-reversal symmetry couple to electrons. For simplicity, we consider only the A₁ modes that correspond to the "Kekulé" vibration of the honeycomb lattice and couple the electrons between K and K' valleys . The A₁ phonons have been suggested to be the most relevant optical-phonon branch for influencing electrical transport at relatively low temperature , possessing the lowest excitation energy and the strongest coupling to electrons. We note that in the context of the Boltzmann equation, the collision integral for A₁ phonons can also qualitatively describe the effect of the other optical-phonon branches. Similar to the case of short-ranged impurity scattering, the collision integrals for different optical-phonon branches are distinguished only by the factors (1 ± p̂ · q̂)/2 that enhance the electron forward (+) or backward (−) scattering. Furthermore, we use the single-mode Einstein model (dispersionless) Hₒₚₕ to describe the A₁ phonons; the electron-phonon coupling Vₑ₋ₒₚₕ takes the form introduced in Refs. . Two parameters are present: The A₁-phonon frequency ωA and electron-phonon coupling βA (doping and temperature dependent, see the discussion in Sec. II E). We in addition assume that the phonons are in thermal equilibrium, that is, the phonon kinetics are not involved (no drag effect on electrons) since the optical-phonon dispersion is weak .

All of the scattering mechanisms that we consider are depicted in Fig. 1. In particular, the Coulomb interaction mediates three scattering channels that we label A, B, and C. Channel A describes intraband carrier-carrier scattering. Channels B and C encode interband conduction electron-valence hole ("electron-hole") scattering. These involve different kinematic regions of frequency ω and momentum q transfer across the Coulomb line, as channels A and B have |ω| ≪ vF q ("quasi-static"), while channel C has |ω| ≫ vF q ("optical"). Plasmons appear in channel C.

This paper is organized as follows. In Sec. II, we present the main results of our calculations and interpret the experimental data in Ref. . In Sec. III, we transcribe the Boltzmann equation that is derived via the Schwinger-Keldysh formalism, with the collision integrals for the impurity scattering, Coulomb interaction, and electron–optical-phonon scattering corresponding to the Feynman diagrams depicted in Fig. 1. Then we introduce the orthogonal-polynomial method for solving the linearized Boltzmann equation. Results for impurity-only and interaction-limited transport are discussed in more detail in Sec. III C. The collinear singularity of the Coulomb collision integrals and the RPA dynamical screening are discussed in the Appendix.

II. Main Results

In general one has the linear response relations :

$$\mathbf{J} = \sigma \mathcal{E} + \sigma \alpha_{\infty} (-\nabla_{\mathbf{r}} T),$$ (2.1a)

$${\rm J_{Q} = T \alpha_{\infty} \cal E + (\kappa_{\infty} + T \alpha_{\infty}^{2}) (-\nabla_{\bf r} T)},$$ (2.1b)

where J is the charge current, JQ the heat current, E the electrochemical field, ∇rT the temperature gradient, and σ, κ∞, and α∞ are the electrical conductivity, thermal conductivity, and thermoelectric power, respectively. We use the subscript "∞" to indicate bulk thermoelectric transport coefficients. In a finite (mesoscopic) sample, slow imbalance relaxation (recombination-generation) can give rise to different transport coefficients and/or a spatially inhomogeneous response , but we do not consider this possibility here. The Lorenz ratio is:

$$L \equiv \frac{\kappa_{\infty}}{\sigma T},$$ (2.1c)

for which we discuss the validity of the Wiedemann-Franz law L₀ = π²k²B/(3e²) . Solving the linearized Boltzmann equation (3.6), inserting the distribution function solution into Eq. (3.15), and comparing the result to Eq. (2.1), we obtain the transport coefficients.

A. Quantum Kinetic Equation

We derive the quantum kinetic equation for electron (λ = +1) and hole (λ = −1) distribution functions fλ(p, r, t) via the Schwinger-Keldysh technique , where p is the quasiparticle wave vector, r the position, and t the time. In the presence of an electric driven field, the stationary Boltzmann equation takes the form:

$$\left f_{\lambda}(\mathbf{p}, \mathbf{r}) = \mathfrak{St}_{\lambda},$$ (2.2)

where vF is the Fourier derivative, vF is the Fermi velocity parallel to the wave vector p, e > 0 is the elementary charge, and E is the total electric field. The collision integral Stλ incorporates the three scattering mechanisms in the Hamiltonian (1.1) (see also Fig. 1):

$$\mathfrak{St}_{\lambda} = \mathfrak{St}_{\text{imp},\lambda} + \mathfrak{St}_{\text{int},\lambda} + \mathfrak{St}_{\text{oph},\lambda},$$ (2.3)

where Stimp,λ describes elastic scattering induced by impurities , Stint,λ the inelastic Coulomb scattering between quasiparticles , and Stoph,λ the inelastic scattering of carriers by optical phonons .

Assuming that the distribution functions \( f_{\lambda=\pm1}(p,r) \) are diagonal in valley and spin space , we present the explicit expressions for the collision integrals in Eq. (2.3) in Sec. III A. In the hydrodynamic regime, the response to static fields is dominated by the zero modes of the inelastic carrier-carrier collision integrals, associated with energy and momentum conservation . For this reason, we also neglect the weak off-diagonal components in electron-hole space, which do not directly contribute to the dc response .

B. Benchmark: Impurity-Only Transport

In the presence of only elastic scattering , the linearized Boltzmann equation can be solved exactly (Sec. III C 1). Mott's formula for the thermoelectric power \( \alpha_{\infty} \) manifestly applies , although the integral form must be employed away from Fermi degeneracy.

The short-ranged-impurity-only transport coefficients take simple expressions:

[ \sigma^{(\rm s)}_{\rm imp} = N \frac{e^{2}}{h} \widetilde{g}^{-1}, \quad \alpha^{(\rm s)}_{\infty,{\rm imp}} = 0, \quad \kappa^{(\rm s)}_{\infty,{\rm imp}} = N \frac{\pi^{2}k_{\rm B}^{2}T}{3h} \widetilde{g}^{-1}, \tag{2.4} \]

where \( \widetilde{g} \) is the effective dimensionless short-ranged disorder strength . Note that the Wiedemann-Franz law is manifestly satisfied. The parameter \( N \) is the number of independent two-component Dirac species, equal to four in graphene.

Another analytically solvable limit is the long-ranged-impurity-only case in the absence of screening . Especially, at the charge neutral point \( \mu = 0 \), the Lorenz ratio is independent of temperature:

[ L_{\infty,\mathrm{imp}}^{(1)}(Q_{\mathrm{TF}} \to 0, \mu = 0) = \frac{21}{5} \frac{\pi^{2}k_{\mathrm{B}}^{2}}{3e^{2}}. \tag{2.5} \]

Here \( Q_{\mathrm{TF}} \) denotes the temperature and density-dependent Thomas-Fermi wave vector . Equation (2.5) violates the Wiedemann-Franz law and enhances the Lorenz constant by a factor \( L/L_0 = \frac{21}{5} = 4.2 \).

The transport coefficients due to the combination of short-ranged disorder and screened Coulomb impurities are shown in Figs. 2(i)–2(iii). We compare results obtained by the orthogonal-polynomial algorithm to the exact results evaluated by Eq. (3.44). We used the dimensionless short-ranged impurity strength \( \widetilde{g} \) and Coulomb impurity density \( n_{\mathrm{imp}} \) determined by fitting the low-temperature, density-dependent conductivity data in Ref. . The Thomas-Fermi screening is limited by the fine structure constant:

[ \alpha_{\rm int} = \frac{2e^{2}}{(\kappa_{1} + \kappa_{2})\hbar\nu_{\rm F}}, \tag{2.6} \]

where \( \kappa_{1,2} \) denotes the permittivities of the media above and below the graphene sheet.

Here, we take \( \alpha_{\rm int} = 0.6 \), appropriate for BN encapsulation . We keep the order of the polynomial basis up to \( N = 16 \) in order to recover the analytical result. We observe that in the presence of Coulomb impurities, the Wiedemann-Franz law is in general broken. As shown in Fig. 2(iii), the Lorenz ratio \( L \) is a function of the charge density \( n \) and temperature \( T \) for a fixed \( \alpha_{\rm int} \).

There exist two interesting limits:

• When \( T \to \infty \), the effective long-ranged impurity strength vanishes , so that short-ranged impurity scattering dominates transport and the Wiedemann-Franz law is restored.

• When \( T \to 0 \), the long-ranged impurity becomes dominant. At the charge neutral point \( n = 0 \), the Thomas-Fermi wave vector \( Q_{\mathrm{TF}} \) divided by the temperature becomes a constant . The Lorenz ratio is enhanced relative to the Wiedemann-Franz law by a numerical constant depending on the fine structure constant (\( L/L0 \approx 2.093 \) for \( \alpha{\rm int} = 0.6 \)). The Wiedemann-Franz law is recovered at sufficiently high charge densities \( n \neq 0 \) and/or temperatures.

In Fig. 2(iv), we show the Lorenz ratio \( L \) as a function of the fine structure constant \( \alpha_{\rm int} \) at charge neutrality in the absence of short-ranged impurity scattering. It is clear that for any finite \( \alpha_{\rm int} \), the Wiedemann-Franz law is broken. Especially, for \( \alpha_{\rm int} \to 0 \), we obtain \( L/L_0 = \frac{21}{5} = 4.2 \) , which provides an upper bound for the Lorenz ratio induced solely by impurities.

C. Crossover From Interaction-Limited Regime To Disorder-Limited Regime

We combine Coulomb interactions and short-ranged impurity to verify the predictions of the relativistic hydrodynamic theory . The transport coefficients obtained by the numerical solution of the Boltzmann equation are shown in Fig. 3. Close to charge neutrality \( \mu \sim k_B T \), the conductivity remains finite. This reflects the "minimal" conductivity due to the electron-hole collisions.

In the hydrodynamic (interaction-dominated) regime of primary interest, \( \tau_{\rm in} \ll \tau_{\rm el} \) . Here, \( 1/\tau_{\rm in} \) denotes the inelastic scattering rate due to electron-electron and electron-hole collisions, while \( 1/\tau_{\rm el} \) is the scattering rate due to elastic electron-impurity and (quasi)elastic electron–acoustic-phonon collisions. (In this section, we neglect optical phonons, which are dealt with below.)

Strong inelastic scattering quickly relaxes fluctuations to local equilibrium. Intercarrier scattering is special however, in that it preserves the total energy and momentum of the Dirac fluid . This means that the distribution function for electrons and holes is always close to Fermi-Dirac in some co-moving reference frame, and this translates into strong constraints on kinetic coefficients.

At charge neutrality, the charge flow is decoupled from the momentum flow, and can be relaxed by electron-hole collisions alone. In the interaction-dominated regime, the minimal conductance at the Dirac point is to a first approximation a function only of the dimensionless interaction strength \( \alpha_{\rm int} \) , and is therefore independent of temperature (ignoring logarithmic renormalization effects ).

This is very different from the case of disorder-dominated transport due to Coulomb impurities.

In a disorder-dominated sample, around charge neutrality, scattering off Coulomb impurities leads to a decreasing resistance with temperature, as shown in Fig. 2(i). This can be understood via dimensional analysis, since the resistivity is proportional to the impurity density, and the only other length scale is the thermal de Broglie wavelength: \( \rho(T) \sim n_{\mathrm{imp}} (\hbar v_F / k_B T)^2 \).

In the thermopower experiment , no downturn in resistivity with increasing temperature was observed over the temperature range of interest (130–350 K). Instead, a superlinear rise was seen above 200 K that we attribute to electron–optical-phonon scattering, discussed below. This should be contrasted with earlier high-temperature experiments that observed a decreasing resistance ; the latter can presumably be attributed to disorder-dominated transport .

Away from charge neutrality and at intermediate temperatures, the thermopower shown in Fig. 3(ii) approaches the ideal clean hydrodynamic result:

[ \alpha_{\infty} = \frac{s^}{e n^}, \tag{2.7} \]

which is the thermodynamic entropy per charge; \( s \) denotes the entropy density. At higher densities/lower temperatures, \( \alpha_{\infty} \to 0 \), consistent with the Mott relation . For \( \mu \sim k_B T \), the Lorenz ratio is much larger than that of a Fermi liquid, \( L/L_0 \gg 1 \). Wiedemann-Franz recovers far away from the Dirac point \( \mu \gg k_B T \).

As shown in the inset of Fig. 3(i), the rescaled minimal conductivity \( \alpha^2_{\rm int} \sigma_{\rm min} \) is almost linear in the fine structure constant \( \alpha_{\rm int} \) for \( \alpha_{\rm int} \ll 1 \), which reflects the Coulomb screening effect. In the absence of screening \( \alpha_{\rm int} \gg 1 \), we recover the results of Refs. and . As shown in Fig. 3(iii), the nonmonotonicity of the thermal conductivity \( \kappa_{\infty} \) (or the Lorenz ratio \( L \)) as a function of \( \mu / k_B T \) is simply a consequence of the ideal relativistic thermodynamics . In the interaction-limited regime, the enhancement of the Lorenz ratio diverges as the strength of impurity scattering vanishes. The hydrodynamic enhancement will also dominate over that attributable to Coulomb impurities, Eq. (2.5). For sufficiently weak impurity scattering and in the absence of optical phonons, the hydrodynamic description should generally apply, regardless of the scattering mechanisms that lift the zero modes of the Coulomb collision operator. In reality, both short-ranged and long-ranged impurities are simultaneously present, and the resulting transport coefficients have similar features as shown in Fig. 3.

Our result for the thermopower in the presence of both types of disorder and Coulomb carrier-carrier scattering, but in the absence of optical phonons, is shown in Fig. 5(iii). There it is compared to the experimental results from Ref. . Our numerical results monotonically approach the ideal hydrodynamic limit with increasing temperature, except near charge neutrality where a finite impurity density sends the thermopower to zero as \( n \to 0 \) . The experimental results instead show a saturation of the thermopower midway between the Mott and hydrodynamic bounds. Below we show that the additional inclusion of electron–optical-phonon scattering gives good agreement with the experiment, Fig. 5(i).

D. Optical-Phonon-Limited Transport

The total energy of electrons and holes is no longer conserved in the presence of the optical-phonon bath .

Via Eqs. (3.25), (3.26), (3.39), and (3.32), we calculate the electronic transport coefficients due only to electron–optical-phonon scattering processes. The results are shown in Fig. 4.

We obtain the resistivity as a function of temperature and charge-carrier density that qualitatively coincides with the result in Ref. . The resistivity weakly depends on the charge carrier density, and, moreover, the optical phonons are thermally activated at a temperature about one order of magnitude lower than their frequency. Three temperature regimes can be observed:

1. Collisionless regime (\( T \ll 150 \) K): The resistivity is almost zero since the population of thermally activated phonons is exponentially small when \( T \ll T_A \).

2. Crossover regime (150 K \( \lesssim T \lesssim 400 \) K): The resistance increases superlinearly in temperature.

3. High-temperature regime (\( T \gg 400 \) K): The resistance increases linearly in temperature. For high enough temperatures, the optical phonons play a similar role as impurities, yet the scattering amplitude is enhanced by the Bose-Einstein distribution function \( f_B(T_A / T) \sim T / T_A \). The temperature dependence of the resistivity qualitatively follows the Bose-Einstein distribution function of the optical phonons.

In the crossover regime, the electron–optical-phonon scattering is strongly inelastic. The thermopower due to electron–optical-phonon scattering alone does not follow Mott's formula . Furthermore, the electron-hole imbalance relaxation processes have significant effects at low doping.

E. All Scattering Mechanisms; Comparison To Thermopower Measurements

Finally, we combine all scattering mechanisms to model the data of the experiment in Ref. . In order to interpret the data, we need first to estimate all the effective parameters. Since the graphene sample is encapsulated between two hexagonal-boron-nitride substrates, we estimate the fine structure constant as \( \alpha_{\rm int} = \frac{2e^2}{(\kappa_1 + \kappa_2)\hbar v_F} \approx 0.6 \), where \( \kappa_1 = \kappa_2 \approx 3.8 \) is the dielectric constant of boron nitride .

The dimensionless short-ranged impurity strength \( \widetilde{g} \) and the Coulomb impurity concentration \( n_{\mathrm{imp}} \) are determined by the conductivity data at low temperature and high doping, where inelastic scattering is negligible. According to this analysis, we have \( \widetilde{g} \approx 1.1 \times 10^{-4} \) and \( n_{\mathrm{imp}} \approx 3 \times 10^9 \) cm\(^{-2}\). Finally, the electron–optical-phonon coupling is attained by fitting the electrical conductivity data at high temperatures. Note that to reach a quantitative agreement to the experimental data, we have tuned the optical-phonon frequency to \( T_A = 2200 \) K, which is a little bit higher than the values reported in Refs. and . The reason for this enhancement might be that A phonons are more rigid due to substrate encapsulation or that higher-frequency optical-phonon branches are also involved.

The fitting procedure described above gives an electron–optical-phonon coupling that increases with decreasing temperature, see Ref. for details. This is presumably due to a combination of ultraviolet renormalization and the temperature-dependent Coulomb screening . We leave the theoretical study of the electron–optical-phonon vertex for deeply inelastic energy and momentum transfers to future work.

We have calculated the thermopower for every combination of the scattering sources in Fig. 1 and present the most informative results in Fig. 5. As shown in Fig. 5(i), our theoretical result coincides quantitatively well with the experimental data if we take into account impurities, optical phonons, and Coulomb interactions, yet neglect the electron-hole optical scattering channel C. Figure 5(ii) indicates that the result of the Mott's relation (red dashed line) merely reflects the impurity-only (both short- and long-ranged) thermopower at high doping.

Figure 5(iii) shows the results in the absence of optical phonons, but including short-ranged and Coulomb-impurity scattering, as well as carrier-carrier channels A, B, and C . Although graphene is relatively degenerate for \( n \sim 10^{12} \) cm\(^{-2}\) (\( T_F = 1350 \) K), the Mott relation is not recovered for the measured temperatures. At these high densities, this is due to the pure intraband electron-electron scattering in channel A.

The disorder is so weak in the experiment that we would need very high densities to observe Fermi liquid behavior; in other words, it is possible to be both degenerate and hydrodynamic in a very clean sample. We estimate that at \( T = 170 \) K, the Mott relation would be recovered only at densities above \( n \sim 10^{13} \) cm\(^{-2}\).

Comparing Fig. 5(iii) (results in the absence of optical phonons) to Fig. 5(i), we observe that the optical phonons significantly suppress the thermopower at higher temperatures and drive the system further away from an ideal hydrodynamic fluid. The thermopower \( \alpha \) is well defined and given by Eq. (2.7) in the absence of a mechanism for momentum relaxation.

Away from charge neutrality, however, even within the hydrodynamic regime some such mechanism is necessary to separately define σ and σα in Eq. (2.1a). In general, the ratio α is also sensitive to this mechanism. Here, this role is filled by either disorder or optical phonon scattering. In particular, Coulomb impurities are poorly screened at low temperatures for charge-carrier densities not too large, while optical phonons become important at higher temperatures.

As discussed in Sec. III C 2, the optical-phonon scattering becomes non-negligible when the collision matrix elements satisfy (Mopt)00 (Mimp)00.

For a charge-carrier density, n ∼ 10¹² cm⁻² (TF ≈ TA ≈ 2000 K) and temperature T < 350 K, this leads to T T∗ ∼ TF/ ln (10⁵ α²opt) ∼ 200 K via a simple estimation , based on the parameters in the experiments .

The plasmon pole in the dynamically screened Coulomb interaction can enhance the electron-hole scattering in the Coulomb channel C . This mechanism could strengthen the hydrodynamic response. Comparing Fig. 5(iv) (which includes channel C) to Fig. 5(i) (which neglects it), we conclude that the associated plasmon enhancement is somehow suppressed in the experiments. We propose that this suppression may be due to additional screening by metallic gates that soften the plasmon dispersion, or damping induced by the plasmon–optical-phonon coupling , which is not accounted for in our treatment. Comparing Fig. 5(v) to Fig. 5(i), we observe that these processes also significantly affect the thermopower at lower charge densities.

III. Boltzmann Equation In The Presence Of Impurities, Coulomb Interaction, And Optical Phonons

A. Collision Integrals

The elastic collision integral in Eq. (2.3) gives Fermi's golden rule amplitudes associated to the diagram in Fig. 1(a), and reads

$\mathfrak{S}t_{\text{imp},\lambda}=\mathfrak{S}t_{\text{imp},\lambda}^{(\text{s})}+\mathfrak{S}t_{\text{imp},\lambda}^{(\text{l})},$ (3.1a)

where $St^{(s)}_{\text{imp},\lambda}$ and $St^{(l)}_{\text{imp},\lambda}$ describe the short- and long-ranged impurity scattering, respectively,

$St^{(s)}_{\text{imp},\lambda} = \int_{\bf{q}} \delta({\bf{q}} - {\bf{p}}) \left ,$ (3.1b)

$St^{(l)}_{\text{imp},\lambda} = 2\pi n_{\text{imp}} \int_{\bf{q}} \delta({\bf{q}} - {\bf{p}}) \frac{1 + \hat{\bf{p}} \cdot \hat{\bf{q}}}{2} |U_{\text{eff}}(\omega = 0, |{\bf{p}} - {\bf{q}}|)|^2 .$ (3.1c)

In Eq. (3.1b), the effective short-ranged impurity strengths are $G0 = (2\pi)^2 (2g_A + g{A3})$, $Gf = (2\pi)^2 g_u$, and $G_b = (2\pi)^2 (2g_m + g_v)$ . In Eq. (3.1c), the long-ranged impurity scattering is characterized by the Coulomb impurity number per unit area $n{\text{imp}}$ and the static RPA Coulomb interaction $|U_{\text{eff}}(\omega = 0, {\bf{k}})|^2$ (see Appendix B). The Dirac delta function $\delta({\bf{q}} - {\bf{p}})$ enforces energy conservation. The terms associated to the factors $(1 \pm \hat{\bf{p}} \cdot \hat{\bf{q}})/2$ describe the enhancement of forward (+) and backward (−) scattering. In Eq. (3.1), we have introduced the shorthand notation

$$\int_{\bf{q}}\equiv\int\frac{d^{2}{\bf{q}}}{(2\pi)^{2}}.$$

The Coulomb collision integral is evaluated at the RPA level associated to the three scattering processes depicted in Fig. 1(b),

$St_{\text{int},\lambda} = N \int_{{\bf{p}}2,{\bf{p}}_3,{\bf{p}}_4} \frac{1 + \hat{\bf{p}} \cdot \hat{\bf{p}}_2}{2} \frac{1 + \hat{\bf{p}}_3 \cdot \hat{\bf{p}}_4}{2} (2\pi)^3 \delta^{(3)}({\bf{p}} + {\bf{p}}_4 - {\bf{p}}_2 - {\bf{p}}_3) |U{\text{eff}}({\bf{p}} - {\bf{p}}_2)|^2$

$\times \left\{ f_{\lambda}({\bf{p}}2,{\bf{r}}) [1 - f{\lambda}({\bf{p}}4,{\bf{r}})] f{\lambda}({\bf{p}}3,{\bf{r}}) - [1 - f{\lambda}({\bf{p}}2,{\bf{r}})] f{\lambda}({\bf{p}},{\bf{r}}) f_{\lambda}({\bf{p}}_4,{\bf{r}}) \right\}$ (3.2a)

$+ \delta^{(3)}({\bf{p}} - {\bf{p}}4 - {\bf{p}}_2 + {\bf{p}}_3) |U{\text{eff}}({\bf{p}} - {\bf{p}}2)|^2 \left\{ [1 - f{\lambda}({\bf{p}},{\bf{r}})] f_{\lambda}({\bf{p}}2,{\bf{r}}) f{-\lambda}({\bf{p}}4,{\bf{r}}) [1 - f{-\lambda}({\bf{p}}3,{\bf{r}})] - [1 - f{\lambda}({\bf{p}}2,{\bf{r}})] f{\lambda}({\bf{p}},{\bf{r}}) f_{-\lambda}({\bf{p}}3,{\bf{r}}) [1 - f{-\lambda}({\bf{p}}_4,{\bf{r}})] \right\}$ (3.2b)

$+ \delta^{(3)}({\bf{p}} - {\bf{p}}4 - {\bf{p}}_3 + {\bf{p}}_2) |U{\text{eff}}(-{\bf{p}} - {\bf{p}}2)|^2 \left\{ [1 - f{\lambda}({\bf{p}},{\bf{r}})] f_{\lambda}({\bf{p}}3,{\bf{r}}) f{-\lambda}({\bf{p}}4,{\bf{r}}) - f{\lambda}({\bf{p}},{\bf{r}}) f_{-\lambda}({\bf{p}}2,{\bf{r}}) [1 - f{\lambda}({\bf{p}}3,{\bf{r}})] [1 - f{-\lambda}({\bf{p}}_4,{\bf{r}})] \right\}.$ (3.2c)

Equations (3.2a)–(3.2c) correspond to the Coulomb scattering channels A–C, diagrams (b)i–(b)iii, respectively. The quasiparticle energy and momentum are written in the three-vector form ${\bf{p}} \equiv ({\bf{p}}, p)$, and the three-dimensional Dirac delta functions $\delta^{(3)}(\cdots)$ describe energy and momentum conservation. Channel A, Eq. (3.2a) is electron-electron scattering, while channels B and C, Eqs. (3.2b) and (3.2c) are electron-hole scattering processes. The RPA screened Coulomb interaction takes the form shown in Appendix B. We emphasize that even at charge neutrality, dynamical screening is crucial at finite temperature due to the thermal activation of electron-hole pairs. Interaction-mediated "Auger" imbalance relaxation processes are suppressed because of the linear dispersion of electrons and holes . Due to kinematic constraints, channels A and B act in the "quasistatic" regime $|\omega| \ll v_F q$, while channel C acts in the "optical" regime $|\omega| \gg v_F q$ (Fig. 7). Here, $\omega$ and $q$ are the frequency and momentum transferred across the Coulomb line.

The carrier–optical-phonon scattering is described by the diagrams in Fig. 1(c) and leads to the collision integral

$St_{\text{oph},\lambda} = (2\pi)^2 \beta^2A s_0 \omega_A M \int{\bf{q}} \frac{1 - \hat{\bf{p}} \cdot \hat{\bf{q}}}{2} \left\{ fB(\omega_A) \delta({\bf{p}} - {\bf{q}} - \omega_A) [1 - f{\lambda}({\bf{p}},{\bf{r}})] f_{\lambda}({\bf{q}},{\bf{r}}) - \delta({\bf{p}}

• {\bf{q}} + \omegaA) f{\lambda}({\bf{p}},{\bf{r}}) \right\}$ (3.3a)

$+ \left\{ \delta({\bf{p}} - {\bf{q}} + \omegaA) [1 - f{\lambda}({\bf{p}},{\bf{r}})] f_{\lambda}({\bf{q}},{\bf{r}}) - \delta({\bf{p}}

• {\bf{q}} - \omegaA) f{\lambda}({\bf{p}},{\bf{r}}) \right\}$ (3.3b)

$+ (2\pi)^2 \beta^2A s_0 \omega_A M \int{\bf{q}} \frac{1 - \hat{\bf{p}} \cdot \hat{\bf{q}}}{2} \delta({\bf{p}} + {\bf{q}} - \omega*A) \left\{ f_B(\omega_A) [1

• f*{\lambda}({\bf{p}},{\bf{r}})] - f_{\lambda}({\bf{p}},{\bf{r}}) f_{-\lambda}({\bf{q}},{\bf{r}}) \right\},$ (3.3c)

where $M = 2.0 \times 10^{-23}$ g is the carbon atom mass and $s_0 = 2.62$ Ų the area per carbon atom. Equations (3.3a,3.3b) correspond to the diagrams in Figs. 1(c)i,iii , respectively. We note that the processes (c)ii and (c)iv are absent for acoustic phonon scattering because the acoustic-phonon velocity is much smaller than the Fermi velocity. To compare to the experiment in Ref. , we take the A phonon temperature $T_A \equiv \hbar \omega_A / k_B \approx 2200$ K, larger than in some previous studies .

The coupling strength $\beta_A$ has been suggested to be strongly energy dependent due to renormalization and screening by the Coulomb interactions . We treat $\beta_A$ as a fitting parameter when interpreting the experimental data .

We separate the distribution function $f_{\lambda}({\bf{p}},{\bf{r}})$ into two parts,

$$f_{\lambda}({\bf{p}},{\bf{r}}) \equiv f_{\lambda}^{(0)}({\bf{p}},{\bf{r}}) + \delta f_{\lambda}({\bf{p}},{\bf{r}}),$$ (3.4)

where $f_{\lambda}^{(0)}({\bf{p}},{\bf{r}})$ is the local equilibrium Fermi-Dirac function ($\beta = 1/k_B T$)

$$f_{\lambda}^{(0)}({\bf{p}},{\bf{r}}) = \frac{1}{e^{\beta(\varepsilon_{\bf{p}}

• \mu_{\lambda})} + 1}, \quad \varepsilon_{\bf{p}} = \hbar {\bf{v}}F |{\bf{p}}|, \quad \mu{\lambda} = \lambda \mu,$$ (3.5a)

and $\delta f_{\lambda}({\bf{p}},{\bf{r}})$ is the deviation from the local equilibrium and can be conveniently cast into the form

$$\delta f_{\lambda}({\bf{p}},{\bf{r}}) = \frac{1}{\beta} \left( -\frac{df_{\lambda}^{(0)}}{d\varepsilon_{\bf{p}}} \right) \chi_{\lambda}({\bf{p}},{\bf{r}}).$$ (3.5b)

Via the standard derivation , from Eq. (2.2), we obtain the time-independent linearized Boltzmann's equation for $\chi_{\lambda}$,

$${\bf{v}}F \cdot \left( \lambda e \beta {\cal E} - \frac{p - \lambda \ln z}{T} \nabla{\bf{r}} T \right) = \frac{1}{\hbar \beta} \widetilde{St}_{\lambda},$$ (3.6)

where we have introduced the electrochemical field ${\cal E} \equiv {\bf{E}} + \frac{1}{e} \nabla_{\bf{r}} \mu$, the effective Fermi-Dirac distribution function, and its derivative

$$fk(p,z) \equiv \frac{1}{z^{-\mu} e^{\mu} + 1}, \quad f_k^{\prime}(p,z) \equiv -g{jk} f_k(p,z),$$ (3.7)

which depends on the dimensionless momentum $p = \beta \hbar {\bf{v}}_j |{\bf{p}}|$ and the "fugacity" $z = \exp(\beta \mu)$. On the right-hand side of Eq. (3.6), the linearized collision integral reads

$$\widetilde{St}_{\lambda} = \widetilde{St}_{\text{imp},\lambda} + \widetilde{St}_{\text{int},\lambda} + \widetilde{St}_{\text{oph},\lambda},$$ (3.8)

where the impurity collision integral is

$$\widetilde{St}_{\text{imp},\lambda} = \widetilde{St}_{\text{imp},\lambda}^{(s)} + \widetilde{St}_{\text{imp},\lambda}^{(l)},$$ (3.9a)

with short- and long-ranged components

$$\widetilde{St}_{\text{imp},\lambda}^{(s)} = \int_{\bf{q}} \delta({\bf{p}} - {\bf{q}}) \left \left,$$ (3.9b)

$$\widetilde{St}_{\text{imp},\lambda}^{(l)} = \gamma^2 \int_{\bf{q}} \delta({\bf{p}} - {\bf{q}}) \frac{1 + \hat{\bf{p}} \cdot \hat{\bf{q}}}{2} |U_{\text{eff}}(\omega = 0, |{\bf{p}} - {\bf{q}}|)|^2 \left,$$ (3.9c)

the Coulomb collision integral is

$$\widetilde{St}_{\text{int},\lambda} = \widetilde{St}_{\text{int},\lambda}^{(i)} + \widetilde{St}_{\text{int},\lambda}^{(ii)} + \widetilde{St}_{\text{int},\lambda}^{(iii)},$$ (3.10a)

with components corresponding to channels A–C in Fig. 1(b):

$$\widetilde{St}_{\text{int},\lambda}^{(i)} = 2\pi N \int_{{\bf{k}},{\bf{q}}} \frac{1 + \hat{\bf{p}} \cdot \hat{\bf{p}} - \hat{\bf{q}}}{2} \frac{1 + \hat{\bf{k}} \cdot \hat{\bf{k}} - \hat{\bf{q}}}{2} \delta({\bf{p}} - |{\bf{p}} - {\bf{q}}| - {\bf{k}} + |{\bf{k}} - {\bf{q}}|) |U_{\text{eff}}({\bf{p}} - |{\bf{p}} - {\bf{q}}|, {\bf{q}})|^2 \left,$$ (3.10b)

$$\widetilde{St}_{\text{int},\lambda}^{(ii)} = 2\pi N \int_{{\bf{k}},{\bf{q}}} \frac{1 + \hat{\bf{p}} \cdot \hat{\bf{p}} + \hat{\bf{q}}}{2} \frac{1 + \hat{\bf{k}} \cdot \hat{\bf{k}} - \hat{\bf{q}}}{2} \delta({\bf{p}} - |{\bf{p}} + {\bf{q}}| + {\bf

Collision Matrix

The collision matrix $\hat{M}$ defined by Eq. (3.21) has three parts:

$$\hat{M} = \hat{M}_{\rm imp} + \hat{M}_{\rm int} + \hat{M}_{\rm oph}, \tag{3.27}$$

corresponding to the collision integrals (3.9)–(3.11), respectively.

Impurity Collision Matrix

The impurity collision matrix elements read:

$$\left({\cal M}_{\rm imp}\right)_{mn} = \frac{1}{4\pi} \sum_{\lambda} \int_{0}^{\infty} {\rm d}p\, {\rm f}^{\prime}_{\lambda}(p,z) p^{2} \widetilde{G}(p) \eta_{m}(\lambda,p) \eta_{n}(\lambda,p), \tag{3.28a}$$

with

$$\widetilde{G}(p) = \widetilde{g} + \widetilde{\gamma}^{2} p^{-2} F\left(\frac{q_{\mathrm{TF}}}{2p}\right), \tag{3.28b}$$

where the effective short-ranged disorder strength is $\widetilde{g} \equiv (4g0 + g_f + 3g_b)/(4\pi)$, the effective long-ranged disorder strength is $\widetilde{\gamma}^2 \equiv \gamma^2 \alpha^2{\rm int}/2 = \pi n_{\rm imp} \alpha^2_{\rm int} (\beta - v_F)^2$, and the function

$$F(x) = \int_{0}^{\pi} \mathrm{d}\theta \frac{\sin^2 \theta}{\left(\sin \frac{\theta}{2} + x\right)^2}, \tag{3.29}$$

with $q_{\rm TF}$ the dimensionless Thomas-Fermi wave vector . We note that in the strong-interaction limit $\lim_{\alpha_{\rm int} \to \infty} F(q_{\rm TF}/2p) = 2\pi (p/q_{\rm TF})^2$, so that the Coulomb impurity becomes short-ranged.

Coulomb Interaction Collision Matrix

The Coulomb interaction collision matrix is:

$$\hat{\mathcal{M}}_{\mathrm{int}} = \hat{\mathcal{M}}_{\mathrm{int}}^{\mathrm{A}} + \hat{\mathcal{M}}_{\mathrm{int}}^{\mathrm{B}} + \hat{\mathcal{M}}_{\mathrm{int}}^{\mathrm{C}}, \tag{3.30}$$

where the elements of each component are given by:

$$({\cal M}_{\mathrm{int}}^{A})_{\mathrm{nn}} = \frac{\pi}{2} \sum_{k=\pm1}^{N} \int_{\mathbf{p},\mathbf{k},\mathbf{q}} \frac{1 + \hat{\mathbf{p}} \cdot \widehat{\mathbf{p}-\mathbf{q}}}{2} \frac{1 + \hat{\mathbf{k}} \cdot \widehat{\mathbf{k}-\mathbf{q}}}{2} \delta(p - |\mathbf{p} - \mathbf{q}| - k + |\mathbf{k} - \mathbf{q}|) |\widehat{U}_{\mathrm{eff}}(p - |\mathbf{p} - \mathbf{q}|, q)|^2 \times \widehat{\Xi}_{p, |\mathbf{p} - \mathbf{q}| \lambda, |\mathbf{k} - \mathbf{q}|}^{\lambda, \lambda, \lambda} \cdot [\eta_{n}(\lambda,p) \widehat{\mathbf{p}} - \eta_{n}(\lambda, |\mathbf{p}

• \mathbf{q}|) \widehat{\mathbf{p} - \mathbf{q}} - \eta_{n}(\lambda,k) \widehat{\mathbf{k}} + \eta_{n}(\lambda, |\mathbf{k} - \mathbf{q}|) \widehat{\mathbf{k} - \mathbf{q}}], \tag{3.31a}$$

$$\left({\mathcal{M}}_{\mathrm{int}}^{\mathrm{B}}\right)_{m n} = \frac{\pi}{2} \sum_{\lambda=\pm1} \int_{\mathbf{p},\mathbf{k},\mathbf{q}} \frac{1 + \hat{\mathbf{p}} \cdot \widehat{\mathbf{p}+\mathbf{q}}}{2} \frac{1 + \hat{\mathbf{k}} \cdot \widehat{\mathbf{k}-\mathbf{q}}}{2} \delta(p - |\mathbf{p} + \mathbf{q}| + k - |\mathbf{k} - \mathbf{q}|) |\widehat{U}_{\mathrm{eff}}(p - |\mathbf{p} + \mathbf{q}|, q)|^2 \times \widehat{\Xi}_{p, |\mathbf{p} + \mathbf{q}| \lambda, |\mathbf{k} - \mathbf{q}|}^{\lambda, \lambda, \lambda} \cdot , \tag{3.31b}$$

$$\left({\mathcal{M}}_{\mathrm{int}}^{\mathrm{C}}\right)_{m n} = \frac{\pi}{2} \sum_{\lambda=\pm1} \int_{\mathbf{p},\mathbf{k},\mathbf{q}} \frac{1 - \hat{\mathbf{p}} \cdot \widehat{\mathbf{p} - \mathbf{q}}}{2} \frac{1 - \hat{\mathbf{k}} \cdot \widehat{\mathbf{k} - \mathbf{q}}}{2} \delta(p + |\mathbf{p} - \mathbf{q}| - k - |\mathbf{k} - \mathbf{q}|) |\widehat{U}_{\mathrm{eff}}(p + |\mathbf{p} - \mathbf{q}|, q)|^2 \times \widehat{\Xi}_{p, |\mathbf{p} - \mathbf{q}| \lambda, |\mathbf{k} - \mathbf{q}|}^{\lambda, -\lambda, \lambda, -\lambda} \cdot [\eta_{n}(\lambda,p) \widehat{\mathbf{p}} - \eta_{n}(-\lambda, |\mathbf{p}

• \mathbf{q}|) \widehat{\mathbf{p} - \mathbf{q}} - \eta_{n}(\lambda,k) \widehat{\mathbf{k}} + \eta_{n}(-\lambda, |\mathbf{k} - \mathbf{q}|) \widehat{\mathbf{k} - \mathbf{q}}]. \tag{3.31c}$$

Carrier–Optical-Phonon Collision Matrix

Finally, the carrier–optical-phonon collision matrix reads:

$$(M_{\rm opt})_{mn} = \alpha^2_{\rm opt} \lambda \int_{0}^{\infty} dp f_{\lambda}(p,z) p(p + A) \eta_{n}(\lambda,p) - \alpha^2_{\rm opt} \lambda \int_{A}^{\infty} dp f_{\lambda}(p,z) p(p - A) \eta_{n}(\lambda,p) + \alpha^2_{\rm opt} \lambda \int_{0}^{A} dp f_{\lambda}(p,z) p(A - p) \eta_{n}(\lambda,p) , \tag{3.32}$$

where the effective electron–optical-phonon coupling $\alpha^2_{\rm opt} = \alpha^2_{\rm opt}/(16\pi^2)$.

Orthogonal Polynomials

A computationally efficient method is to choose the basis $\{\eta_n(\lambda,p)\}$ as a set of orthogonal polynomials in two variables $\lambda$ and $p$, taking into account the symmetry properties of the collision matrix and the force vectors. We define the orthonormal condition as:

$$\langle \eta_{m}, \eta_{n} \rangle \equiv \sum_{\lambda=\pm} \int_{0}^{\infty} \mathrm{d}p \, \mathcal{K}_{\lambda}(p,z) \, \eta_{m}(\lambda,p) \, \eta_{n}(\lambda,p) = \delta_{m,n}, \tag{3.33}$$

where the kernel function depends on the fugacity $z$:

$$\mathcal{K}_{\lambda}(p,z) = p \, \mathrm{f}_{\lambda}(p,z), \tag{3.34}$$

and the function $f_{\lambda}(p,z)$ is defined in Eq. (3.7). We note that in general the transport coefficients obtained via Eq. (3.25) are independent of the choice of basis, and moreover, the normalization condition in Eq. (3.33) can be relaxed.

Our objective is to orthonormalize the monomial system $\{p^n, \lambda p^n\}_{n \geq -1}$ by the orthogonal condition (3.33). Note that we have included the negative power $n = -1$ because $\delta f_{\lambda} \sim 1/p$ is the lowest power of $p$ that leads to finite charge and thermal currents in two spatial dimensions according to Eq. (3.15). Via the Gram-Schmidt orthogonalization process, we recursively generate the polynomials in the form:

$$\eta_{-2}(\lambda,p) = \frac{\lambda}{p}, \quad \eta_{-1}(\lambda,p) = \frac{1}{p}, \tag{3.35a}$$

where $\{um\}$ are the monomials $\{p^n, \lambda p^n\}{n \geq 0}$ ordered as follows:

$$\begin{array}{ccccccccccccc} u_0 & u_1 & u_2 & u_3 & u_4 & u_5 & u_6 & \cdots \

\hline p & \lambda & 1 & p^2 & \lambda p & p^3 & \lambda p^2 & \cdots \

\hline \end{array} \tag{3.35c}$$

As discussed in Sec. II, the negative-power basis "$\lambda/p$" is crucial for solving the Boltzmann equation with only disorder or acoustic-phonon scattering processes. The leading positive-power basis "$p$", "$\lambda$", and "$1$" multiplied by the Fermi velocity $v_F$ correspond to the momentum, charge velocity, and energy velocity, respectively, and play the key role in the hydrodynamic description.

In order to calculate the coefficients $\{a_{mn}\}$ in Eq. (3.35b), we introduce the functions:

$$\Omega_{n,\pm}(z) = -\Gamma(n+1) \left, \tag{3.36}$$

where $\Gamma(n)$ is the gamma function and $\mathrm{Li}_n(-z)$ the polylogarithm defined by:

$$\mathrm{Li}n(-z^{\lambda}) \equiv -\frac{1}{\Gamma(n)} \int{0}^{\infty} \mathrm{d}p \, p^{n-1} \, \mathrm{f}_{\lambda}(p,z), \quad n \geq 0. \tag{3.37}$$

The leading coefficients read:

$$a_{00} = \frac{1}{\sqrt{\Omega_{3,+}(z)}}, \quad a_{10} = -a_{11} \frac{\Omega_{2,-}(z)}{\Omega_{3,+}(z)}, \tag{3.38}$$

$$a_{11} = \sqrt{\frac{\Omega_{3,+}(z)}{\Omega_{1,+}(z) \Omega_{3,+}(z) - \Omega_{2,-}^2(z)}},$$

and for $n \geq 0$:

$$\etan(\lambda,p) = \sum{m=0}^{n} a_{nm} u_m, \tag{3.35b}$$

which are important for writing down the force vectors . Higher-order coefficients can be generated numerically.

Substituting Eq. (3.35) into Eq. (3.20), we obtain the force vectors:

$$\begin{array}{l} {\cal F}^{(1)} = \frac{1}{\sqrt{4\pi}} \left^{\sf T}, \

{\cal F}^{(2)} = \frac{1}{\sqrt{4\pi}} \left^{\sf T}, \

\end{array} \tag{3.39}$$

where we have used $\Omega_{0,+}(z) = 1$ and the coefficients $a_{00}, a_{10}, a_{11}$ are given in Eq. (3.38). In Eq. (3.39), only the leading four (three) components of the force vectors ${\cal F}^{(1)}$ are nonzero so that we only need the $4 \times 4$ block of the inverse collision matrix $\hat{M}^{-1}$ to evaluate the transport coefficients in Eq. (3.25).

Thermodynamics

We present some useful thermodynamic relations for the ideal two-component relativistic gas. The charge-carrier density $n$ and the internal energy density $\varepsilon$ are fixed by the Fermi-Dirac function as:

$$n = \sum_{\lambda=\pm1} \lambda \int_{\bf p} f_{\lambda}^{(0)}({\bf p}), \quad \varepsilon = \sum_{\lambda=\pm1} \int_{\bf p} \varepsilon_{\bf p} f_{\lambda}^{(0)}({\bf p}), \tag{3.40}$$

which leads to:

$$n = \frac{N k_{\mathrm{B}}^{2} T^{2}}{4 \pi \hbar^{2} \mathbf{v}_{\mathrm{F}}^{2}} \Omega_{2,-}(z), \tag{3.41a}$$

$$\varepsilon = \frac{N k_{\mathrm{B}}^{3} T^{3}}{6 \pi \hbar^{2} \mathbf{v}_{\mathrm{F}}^{2}} \Omega_{3,+}(z), \tag{3.41b}$$

where $\Omega_{n,\pm}(z)$ are defined in Eq. (3.36). One can use the charge carrier density (3.41a) to determine the fugacity $z = z(n,T)$.

Moreover, the enthalpy density $h$ and entropy density $s$ obey the thermodynamic relations $h = \varepsilon + P$ and $T s = h - n \mu$, where $P$ is the pressure. For the ideal relativistic gas, we exploit scale invariance $h = 3P$ so that:

$$h = \frac{3}{2} \varepsilon, \quad s = \frac{1}{T} \left(\frac{3}{2} \varepsilon - n \mu\right). \tag{3.42}$$

Explicit formulas for all thermodynamic potentials in terms of $n$ and $T$ are useful for analyzing the hydrodynamic description (Sec. II C). Transport coefficients in the interaction-limited regime are expressed in terms of these (irrespective of Fermi degeneracy), see Eq. (3.48).

Transport Coefficients

Impurity-Only Transport

In the presence of only elastic scattering, the linearized Boltzmann equation can be solved by:

$$\chi_{\lambda}({\bf p}) = \Phi_{\lambda}(p) \cos \varphi_{\bf p}, \tag{3.43a}$$

$$\Phi_{\lambda}(p) = \frac{2 \hbar \beta \nu_{\rm F}}{\widetilde{G}(p)} \left, \tag{3.43b}$$

where the dimensionless "scattering rate" $\widetilde{G}(p)$ is defined in Eq. (3.28b). Substituting Eqs. (3.4), (3.5b), (3.43), and (3.43b) into the currents (3.15), we obtain the transport coefficients in the form of Eq. (3.25):

$$\text{L}_{11}^{(\text{imp})} = \sum_{\lambda=\pm1} \int_{0}^{\infty} \text{d}p \, \text{f}_{\lambda}^{\prime}(p,z) \, ^{-1}, \tag{3.44a}$$

$$\text{L}_{12}^{(\text{imp})} = \sum_{\lambda=\pm1} \int_{0}^{\infty} \text{d}p \, \text{f}_{\lambda}^{\prime}(p,z) \, \lambda p ^{-1}, \tag{3.44b}$$

$$\text{L}_{22}^{(\text{imp})} = \sum_{\lambda=\pm1} \int_{0}^{\infty} \text{d}p \, \text{f}_{\lambda}^{\prime}(p,z

where $\delta L^{(h)}_{ij} \sim \alpha^{-2}_{\text{int}} \sim O(1)$ ($i,j \in \{1,2\}$) encode the contributions of polynomial modes $\eta_{n1}(\lambda, p)$ orthogonal to $u_0 = p$; these can be evaluated numerically. Substituting Eq. (3.47) into Eq. (3.25) and exploiting the thermodynamic relations (3.41) and (3.42), we obtain the interaction-limited transport coefficients:

$$ \sigma_{\mathrm{h}} = \frac{(e n \nu_{\mathrm{F}})^2 \tau_{\mathrm{el}}}{\hbar} + \sigma_{\mathrm{min}}, $$

$$ \alpha_{\infty, \mathrm{h}} = \frac{e n \nu_{\mathrm{F}}^2 \tau_{\mathrm{el}}}{T \sigma_{\mathrm{h}}} - \frac{\mu}{e T}, $$

$$ \kappa_{\infty, \mathrm{h}} = \frac{\hbar \nu_{\mathrm{F}}^2 \tau_{\mathrm{el}}}{T \sigma_{\mathrm{h}}} \frac{\sigma_{\mathrm{min}}}{\sigma_{\mathrm{h}}}. $$

Except for a slight discrepancy in the form of the thermal conductivity (discussed below), these match the results of relativistic hydrodynamics .

In Eq. (3.48), $\tau^{-1}_{\text{el}}$ is the elastic scattering rate induced by short-ranged impurities. It is defined by:

$$ \tau_{\text{el}}^{-1} \equiv \frac{k_{\text{B}} T}{2 \hbar} \frac{\Omega_{4,+}(z)}{\Omega_{3,+}(z)} \widetilde{g} \approx \begin{cases} \frac{7 \pi^4}{270 \xi(3)} \frac{k_{\text{B}} T}{\hbar} \widetilde{g}, & |\mu| \ll k_{\text{B}} T, \

\frac{1}{2} \frac{|\mu|}{\hbar} \widetilde{g}, & |\mu| \gg k_{\text{B}} T, \end{cases} $$

where $\zeta(n)$ is the Riemann zeta function. For comparison, we estimate the inelastic-scattering rate due to the Coulomb interactions :

$$ \tau_{\text{ee}}^{-1} \sim \begin{cases} \frac{k_{\text{B}} T}{\hbar} \alpha_{\text{int}}^2, & |\mu| \ll k_{\text{B}} T, \

\frac{(k_{\text{B}} T)^2}{\hbar |\mu|} \alpha_{\text{int}}^2, & |\mu| \gg k_{\text{B}} T. \end{cases} $$

Here, we note that the expression for $|\mu| \gg k_{\text{B}} T$ is the standard Fermi liquid behavior arising from channel A .

The minimal conductivity $\sigma_{\text{min}}$, and the related parameter $\overline{\sigma}_{\text{min}}$ appearing in the thermal conductivity , take the form:

$$ \sigma_{\text{min}} = N \frac{e^2}{h} \delta L_{11}^{(\text{h})}, \quad \overline{\sigma}_{\text{min}} = \sigma_{\text{min}} (1 + \varrho), $$

where $\varrho$ is an enhancement factor:

$$ \varrho = \left( \frac{n}{\beta \hbar} \right)^2 \frac{\delta L_{22}^{(\text{h})}}{\delta L_{11}^{(\text{h})}} - 2 \left( \frac{n}{\beta \hbar} \right) \frac{\delta L_{12}^{(\text{h})}}{\delta L_{11}^{(\text{h})}}. $$

This factor was not taken into account in previous works.

We emphasize four points:

1. As long as the Coulomb interactions dominate the collisions of electrons and/or holes so that $\tau_{\text{ee}}$ is the shortest scattering time, the hydrodynamic description (3.48) applies, where, however, the expression for the scattering rate $\tau^{-1}_{\text{el}}$ and the values of $\delta L^{(\text{h})}_{ij}$ should be determined by the mechanism that lifts the zero modes of the Coulomb collision integrals.

2. As a simplification, if considering the effect of the impurity collision matrix $\hat{M}^{(\text{s})}_{\text{imp}}$ only by its projection onto the zero modes of the Coulomb collision $\hat{M}_{\text{int}}$ , one can show that $\delta L^{(\text{h})}_{12} = \delta L^{(\text{h})}_{22} = 0$, so that $\varrho = 0$ and $\overline{\sigma}_{\text{min}} = \sigma_{\text{min}}$ . In this case, the thermal conductivity (3.48c) precisely recovers the expression in Refs. and .

3. The minimal conductivity $\sigma_{\text{min}}$ dominates the charge conductivity only in the vicinity of the charge neutrality, i.e., for $|\mu| / k_{\text{B}} T \sim \widetilde{g} / \alpha^2_{\text{int}}$ . In this regime, we can show that the enhancement factor $\varrho \sim O(\widetilde{g})$ can be neglected, which is consistent with the conclusion in Refs. and .

4. At high density $|\mu| \gg k_{\text{B}} T$, the impurity scattering starts to dominate when $\tau_{\text{el}} \sim \tau_{\text{ee}}$ , which leads to $|\mu| / k_{\text{B}} T \gg \alpha^2_{\text{int}} / \widetilde{g}$, and the expansion (3.47) is no longer justified.

In the ideal hydrodynamic regime :

$$ \sqrt{\frac{\widetilde{g}}{\alpha_{\text{int}}^2}} \ll \frac{|\mu|}{k_{\text{B}} T} \ll \sqrt{\frac{\alpha_{\text{int}}^2}{\widetilde{g}}}, $$

the thermoelectric power approaches the thermodynamic expression:

$$ \alpha_{\infty, \text{h}} \to \frac{\mathbf{S}}{e n}. $$

The thermal conductivity takes the form:

$$ \frac{\kappa_{\infty, \text{h}}}{T} \to L_{0, \text{h}} \overline{\sigma}_{\text{min}}, $$

where we define the Lorenz ratio of an ideal relativistic gas as:

$$ L_{0, \text{h}} \equiv \left( \frac{\hbar}{e n T} \right)^2. $$

Moreover, the Lorenz ratio $L$ tends to diverge as $|\mu| / k_{\text{B}} T$ approaches the lower bound $\widetilde{g} / \alpha^2_{\text{int}}$:

$$ L \to C \frac{\alpha_{\text{int}}^4}{\widetilde{g}} \frac{h \sigma_{\text{min}}}{N e^2} \left( \frac{k_{\text{B}}}{e} \right)^2, $$

where the constant $C = \frac{189}{1280} \left^2 \approx 90.03$. Manifestly, both Mott's law and the Wiedemann-Franz law are violated in the ideal hydrodynamic regime. By contrast, for $|\mu| / k_{\text{B}} T \gg \alpha^2_{\text{int}} / \widetilde{g}$, we recover the disorder-limited behavior in Eq. (2.4) .

Acknowledgments

We are grateful to Fereshte Ghahari and Philip Kim for sharing their experimental thermopower data before publication and to Kin Chung Fong, Jesse Crossno, and Markus Mueller for stimulating discussions. This research was supported by the Welch Foundation under Grant No. C-1809 and by an Alfred P. Sloan Research Fellowship (No. BR2014-035).

Appendix A: Elliptic Coordinate System for the Coulomb Collisions (3.31)

To evaluate the Coulomb collision matrix , we first perform the integration over the momentum transfer $q$ and for the moment keep the incoming and outgoing momenta $p$ and $k$ constant. It is convenient to solve the energy conservation constraint by parameterizing $q$ in the elliptic (or hyperbolic) coordinate system , where the collinear scattering singularity is shown explicitly. However, as discussed in Appendix B, this singularity is compensated by a line of zeros in the RPA screened Coulomb interaction along the forward scattering direction $\omega = v_{\text{F}} q$.

For channel B , one needs to evaluate an integral in the form:

$$ \mathcal{I}_{\mathbf{B}}(\mathbf{p}, \mathbf{k}) = \int_{\mathbf{q}} (2\pi) \delta(p - |\mathbf{p} + \mathbf{q}| + k - |\mathbf{k} - \mathbf{q}|) G(\mathbf{p}, \mathbf{k}; \mathbf{q}), $$

where $G(\mathbf{p}, \mathbf{k}; \mathbf{q})$ is a general function of $\mathbf{p}$, $\mathbf{k}$, and $\mathbf{q}$. The elliptic coordinates $(\rho, \theta)$ of $q$ are defined by:

$$ \begin{pmatrix} q_x \

q_y \end{pmatrix} = \frac{1}{2} \begin{pmatrix} k_x - p_x \

k_y - p_y \end{pmatrix}

• \frac{1}{2} \begin{pmatrix} k_x + p_x & -(k_y + p_y) \
  
  k_y + p_y & k_x + p_x \end{pmatrix} \begin{pmatrix} \cosh \rho \cos \theta \
  
  \sinh \rho \sin \theta \end{pmatrix}, $$

where $0 \leq \rho < \infty$, $-\pi \leq \theta < \pi$. As shown in Fig. 6(i), the momentum transfer $q$ lies on an ellipse with two foci at the incident momenta $-\mathbf{p}$ and $\mathbf{k}$. Using Eq. (A2), we readily obtain the relations:

$$ |\mathbf{p} + \mathbf{q}| + |\mathbf{k} - \mathbf{q}| = |\mathbf{k} + \mathbf{p}| \cosh \rho, $$

$$ \int_{\mathbf{q}} = \int_{-\pi}^{\pi} \frac{d\theta}{(2\pi)^2} \int_0^{\infty} d\rho \frac{|\mathbf{k} + \mathbf{p}|^2}{4} (\cosh^2 \rho - \cos^2 \theta). $$

Substituting Eq. (A4) into Eq. (A1), we obtain:

$$ \mathcal{I}_{\mathbf{B}}(\mathbf{p}, \mathbf{k}) = \frac{|\mathbf{k} + \mathbf{p}|}{4 \sinh \rho0} \int{-\pi}^{\pi} \frac{d\theta}{2\pi} (\cosh^2 \rho_0 - \cos^2 \theta) G, $$

where $\rho_0$ is fixed by:

$$ \cosh \rho_0 = \frac{k + p}{|\mathbf{k} + \mathbf{p}|}. $$

The phase space of the collinear collision is manifestly divergent since $\rho_0 = 0$ when $k + p = |\mathbf{k} + \mathbf{p}|$.

For channels A and C , we need to evaluate:

$$ \mathcal{I}_{\text{C/A}}(\mathbf{p}, \mathbf{k}) = \int_{\mathbf{q}} (2\pi) \delta(p - k \pm |\mathbf{p} - \mathbf{q}| \mp |\mathbf{k} - \mathbf{q}|) G(\mathbf{p}, \mathbf{k}; \mathbf{q}), $$

where the upper (lower) signs are for channel C (A). The hyperbolic coordinates are defined by:

$$ \begin{pmatrix} q_x \

q_y \end{pmatrix} = \frac{1}{2} \begin{pmatrix} k_x + p_x \

k_y + p_y \end{pmatrix}

• \frac{1}{2} \begin{pmatrix} k_x - p_x & -k_y + p_y \
  
  k_y - p_y & k_x - p_x \end{pmatrix} \begin{pmatrix} \cosh \rho \cos \theta \
  
  \sinh \rho \sin \theta \end{pmatrix}, $$

where $\rho$ and $\theta$ are defined in the intervals in Eq. (A3). As shown in Fig. 6(ii), the momentum transfer $q$ lies on the two branches of a hyperbola (dashed curves) with two foci at the incident momenta $\mathbf{p}$ and $\mathbf{k}$. Via Eq. (A8), we obtain:

$$ |\mathbf{p} - \mathbf{q}| - |\mathbf{k} - \mathbf{q}| = |\mathbf{k} - \mathbf{p}| \cos \theta, $$

$$ \int_{\mathbf{q}} = \int_{-\pi}^{\pi} \frac{d\theta}{(2\pi)^2} \int_0^{\infty} d\rho \frac{|\mathbf{k} - \mathbf{p}|^2}{4} (\cosh^2 \rho - \cos^2 \theta). $$

Substituting Eq. (A9) into Eq. (A7) gives:

$$ \mathcal{I}_{\text{C/A}}(\mathbf{p}, \mathbf{k}) = \frac{|\mathbf{k} - \mathbf{p}|}{8 \pi |\sin \theta_0|} \int_0^{\infty} d\rho (\cosh^2 \rho - \cos^2 \theta_0) \{ G + G \}, $$

where $\theta_0$ is fixed by:

$$ \cos \theta_0 = \pm \frac{k - p}{|\mathbf{k} - \mathbf{p}|}, \quad 0 \leq \theta_0 < \pi, $$

and the sign "+" ("−") is for channel C (A). The phase space of the collinear collision is divergent since $\theta_0 = 0$ when $k - p = \pm |\mathbf{k} - \mathbf{p}|$.

We note that for channel A, Eq. (A11) is consistent with Eq. (A9a), so that the zero momentum transfer condition $q = 0$ resides on the corresponding branch of the hyperbola shown in Fig. 6(ii). For channel C, Eq. (A11) is in general inconsistent with Eq. (A9a), so that $q = 0$ does not reside on this branch. The exception has $p = k$.

Appendix B: RPA Screening of Coulomb Interaction, Cancellation of the Collinear Collision Singularity, and Plasmon Enhancement in Channel C

In this appendix, we employ physical units, but we set $\hbar = 1$ unless noted. At finite temperature, the screened Coulomb interaction takes the form (see Fig. 7):

$$ U_{\text{eff}}(\omega, q) = \frac{Vq}{\epsilon(q, \omega)}, \quad V_q = \frac{2 \pi \alpha{\text{int}} v_{\text{F}}}{q}, $$

where $\alpha_{\text{int}}$ is the fine structure constant and $\epsilon(q, \omega)$ is the dynamical screening function. In the RPA approximation, $\epsilon(q, \omega) = 1 - V_q \chi_0(q, \omega)$, where $\chi_0(q, \omega)$ is the Lindhard function. The real and imaginary parts of $\chi_0(q, \omega)$ take the forms :

$$ \operatorname{Re} \chi0(q, \omega) = -\Xi(T, z) - \frac{q^2}{\sqrt{|\mathbf{v}{\text{F}}^2 q^2 - \omega^2|}} \begin{cases} \Lambda_{-}(q, \omega), & |\omega| \leq \mathbf{v}_{\text{F}} q, \

\Lambda(q, \omega), & |\omega| \geq \mathbf{v}_{\text{F}} q, \end{cases} $$

$$ \operatorname{Im} \chi0(q, \omega) = \frac{q^2}{\sqrt{|\mathbf{v}{\text{F}}^2 q^2 - \omega^2|}} \begin{cases} \Lambda(q, \omega), & |\omega| \leq \mathbf{v}_{\text{F}} q, \

-\Lambda_{+}(q, \omega), & |\omega| \geq \mathbf{v}_{\text{F}} q, \end{cases} $$

where:

$$ \Xi(T, z) = \frac{N k_{\text{B}} T}{2 \pi \nu_{\text{F}}^2} \ln, $$

$$ \Lambda_{\pm}(q, \omega) = \frac{N}{8 \pi} \left, $$

$$ \Lambda(q, \omega) = \frac{N}{8 \pi} \sum_{\lambda, \lambda' = \pm} \lambda G_{\lambda, \lambda'}(q, \omega). $$

Here, $H_{\pm, \lambda}(q, \omega)$ and $G_{\lambda, \lambda'}(q, \omega)$ are defined by:

$$ G_{\lambda, \lambda'}(q, \omega) = \int1^{+\infty} du \frac{\sqrt{u^2 - 1}}{z^{-\lambda'} \exp\left( \frac{|\mathbf{v}{\text{F}} q u + \lambda \omega|}{2 k_{\text{B}} T} \right) + 1}, $$

$$ H_{\lambda, \lambda'}(q, \omega) = \int_{-1}^1 du \frac{\sqrt{1 - u^2}}{z^{-\lambda'} \exp\left( \frac{|\mathbf{v}_{\text{F}} q u + \lambda \omega|}{2 k_{\text{B}} T} \right) + 1}. $$

In numerical calculations, we use the dimensionless form of the Coulomb interaction:

$$ \widetilde{U}_{\text{eff}}(\omega, q) \equiv \frac{1}{\beta (\hbar \nu_{\text{F}})^2} U_{\text{eff}} \left( \frac{1}{\beta \hbar} \omega, \frac{1}{\beta \hbar \nu

$$v_{F}^{2}q^{2} - \omega^{2} v_{F}^{2}|\sin \theta|_{\theta=\theta_{0}} = k p \left[-|k - p| + (k + p) \cosh \rho\right] \sin\left(\frac{\phi_{k} - \phi_{p}}{2}\right) - \text{sgn}(\theta_{0})k p \sinh \rho \sin (\phi_{k} - \phi_{p}), \quad (B9a)$$
$$v_{F}^{2}q^{2} - \omega^{2} v_{F}^{2} \sinh \rho_{\rho=\rho_{0}} = k p \left[|k + p| + (k - p) \cos \theta\right] \sin\left(\frac{\phi_{k} - \phi_{p}}{2}\right) + k p \sin \theta \sin (\phi_{k} - \phi_{p}), \quad (B9b)$$
$$\omega^{2} - v_{F}^{2}q^{2} v_{F}^{2}|\sin \theta|_{\theta=\theta_{0}} = k p \left[|k - p| + (k + p) \cosh \rho\right] \sin\left(\frac{\phi_{k} - \phi_{p}}{2}\right) + \text{sgn}(\theta_{0})k p \sinh \rho \sin (\phi_{k} - \phi_{p}), \quad (B9c)$$

These equations correspond to channels **A**, **B**, and **C**, respectively.

In our calculation, we further simplify the RPA screened interaction (B6). For \(v_{F}q > |\omega|\) (channels **A** and **B**), we take the static limit \(\omega \to 0\) of the denominator and neglect the dielectric enhancement arising from the residue \(-(q,\omega = 0)\):

$$|U_{\alpha}(\omega,q)|^{2} \approx \left(\frac{2\pi\alpha_{m}}{q}\right)^{2} \frac{\nu_{q}^{2}q^{2} - \omega^{2}}{(q + Q\pi)^{2}}, \quad \nu_{q}q > |\omega|.$$

For \(\nu_{q}q < |\omega|\) (channel **C**), up to the leading order in \(\nu_{q}q/|\omega| \ll 1\), we obtain:

$$|U_{d}(\omega,q)|^{2} \approx \frac{(2\pi\omega_{0})^{2}(\omega^{2} - \nu_{0}^{2}q^{2})}{\frac{1}{\tau_{r}^{2}}\left(\sqrt{\omega^{2} - \nu_{0}^{2}q^{2}}(q + \Omega\nu) - \omega \Omega\nu\right)^{2} + \frac{1}{2}\left[\frac{\pi\omega_{0}\nu q^{2} + \omega\omega\omega\omega}{\omega\omega\omega\nu + \omega\omega\omega_{0}}\right]^{2}}, \quad \forall q < |\omega|.$$

As compared in Fig. **8**, Eq. (**810**) provides a good approximation for the RPA interaction (**81**).

The plasmon dispersion [see Fig. 9(ii)] is determined by \(1 - V_{q} \text{Re}\chi_{0} = 0\), which gives for \(v_{F}q/|\omega| \ll 1\):

$$\omega_{\mathrm{p}}(q) = v_{F} \sqrt{\frac{q(q + Q_{\mathrm{TF}})^{2}}{q + 2Q_{\mathrm{TF}}}}. \quad (B11)$$

Channel **C** is strongly enhanced along this plasmon dispersion. However, the presence of metallic gates or plasmon-phonon coupling [58] may significantly broaden the plasmon peak so that the enhancement of channel **C** is suppressed.

---

### References

1. For a review, see, for example, A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. Novoselov, and A. K. Geim, *Rev. Mod. Phys.* **81**, 109 (2009).
2. S. Das Sarma, S. Adam, E. H. Hwang, and E. Rossi, *Rev. Mod. Phys.* **83**, 407 (2011).
3. I. L. Aleiner and K. B. Efetov, *Phys. Rev. Lett.* **97**, 236801 (2006); J. P. Robinson, H. Schomerus, L. Oroszlany, and V. I. Fal'ko, *ibid.* **101**, 196803 (2008).
4. S. Adam, E. H. Hwang, V. M. Galitski, and S. Das Sarma, *Proc. Natl. Acad. Sci. USA* **104**, 18392 (2007).
5. P. M. Ostrovsky, I. V. Gornyi, and A. D. Mirlin, *Phys. Rev. Lett.* **98**, 256801 (2007); J. H. Bardarson, J. Tworzydlo, P. W. Brouwer, and C. W. J. Beenakker, *ibid.* **99**, 106801 (2007); S. Ryu, C. Mudry, H. Obuse, and A. Furusaki, *ibid.* **99**, 116601 (2007); K. Nomura, M. Koshino, and S. Ryu, *ibid.* **99**, 146806 (2007); K. Nomura, S. Ryu, M. Koshino, C. Mudry, and A. Furusaki, *ibid.* **100**, 246806 (2008).
6. S. Adam, E. H. Hwang, and S. Das Sarma, *Physica E* **40**, 1022 (2008).
7. E. H. Hwang and S. Das Sarma, *Phys. Rev. B* **79**, 165404 (2009).
8. C. Jang, S. Adam, J.-H. Chen, E. D. Williams, S. Das Sarma, and M. S. Fuhrer, *Phys. Rev. Lett.* **101**, 146805 (2008).
9. K. Nomura and A. H. MacDonald, *Phys. Rev. Lett.* **96**, 256602 (2006); **98**, 076602 (2007).
10. A. B. Kashuba, *Phys. Rev. B* **78**, 085415 (2008).
11. L. Fritz, J. Schmalian, M. Muller, and S. Sachdev, *Phys. Rev. B* **78**, 085416 (2008).
12. M. Muller and S. Sachdev, *Phys. Rev. B* **78**, 115419 (2008).
13. M. S. Foster and I. L. Aleiner, *Phys. Rev. B* **77**, 195413 (2008).
14. M. Schutt, P. M. Ostrovsky, I. V. Gornyi, and A. D. Mirlin, *Phys. Rev. B* **83**, 155441 (2011).
15. S. A. Hartnoll, P. K. Kovtun, M. Muller, and S. Sachdev, *Phys. Rev. B* **76**, 144502 (2007).
16. M. Muller, L. Fritz, and S. Sachdev, *Phys. Rev. B* **78**, 115406 (2008).
17. M. S. Foster and I. L. Aleiner, *Phys. Rev. B* **79**, 085415 (2009).
18. B. N. Narozhny, I. V. Gornyi, M. Titov, M. Schutt, and A. D. Mirlin, *Phys. Rev. B* **91**, 035414 (2015).
19. E. H. Hwang and S. Das Sarma, *Phys. Rev. B* **77**, 115449 (2008).
20. D. K. Efetov and P. Kim, *Phys. Rev. Lett.* **105**, 256805 (2010).
21. T. Sohier, M. Calandra, C.-H. Park, N. Bonini, N. Marzari, and F. Mauri, *Phys. Rev. B* **90**, 125414 (2014).
22. S. S. Kubakaddi, *Phys. Rev. B* **79**, 075417 (2009).
23. W. S. Bao, S. Y. Liu, and X. L. Lei, *J. Phys. Condens. Matter* **22**, 315502 (2010).
24. Y. M. Zuev, W. Chang, and P. Kim, *Phys. Rev. Lett.* **102**, 096807 (2009).
25. P. Wei, W. Bao, Y. Pu, C. N. Lau, and J. Shi, *Phys. Rev. Lett.* **102**, 166808 (2009).
26. N. F. Mott and E. A. Davis, *Electronic Processes in Noncrystalline Materials* (Clarendon, Oxford, 1971), p. 47.
27. N. W. Ashcroft and N. D. Mermin, *Solid State Physics* (Saunders College, Philadelphia, 1976), p. 255.
28. H. Julian Goldsmid, *Introduction to Thermoelectricity* (Springer-Verlag, Berlin, 2010), p. 34.
29. J. Crossno, J. K. Shi, K. Wang, X. Liu, A. Harzheim, A. Lucas, S. Sachdev, P. Kim, T. Taniguchi, K. Watanabe, T. A. Ohki, and K. C. Fong, *Science* **351**, 1058 (2016).
30. D. A. Bandurin, I. Torre, R. Krishna Kumar, M. Ben Shalom, A. Tomadin, A. Principi, G. H. Auton, E. Khestanova, K. S. Novoselov, I. V. Grigorieva, L. A. Ponomarenko, A. K. Geim, and M. Polini, *Science* **351**, 1055 (2016).
31. F. Ghahari, H.-Y. Xie, T. Taniguchi, K. Watanabe, M. S. Foster, and P. Kim, *Phys. Rev. Lett.* **116**, 136802 (2016).
32. A. Lucas, J. Crossno, K. C. Fong, P. Kim, and S. Sachdev, *Phys. Rev. B* **93**, 075426 (2016).
33. A. Principi and G. Vignale, *Phys. Rev. Lett.* **115**, 056603 (2015).
34. M. Muller, J. Schmalian, and L. Fritz, *Phys. Rev. Lett.* **103**, 025301 (2009).
35. L. Levitov and G. Falkovich, *arXiv:1508.00836v1* [cond-mat.mes-hall].
36. I. Torre, A. Tomadin, A. K. Geim, and M. Polini, *Phys. Rev. B* **92**, 165433 (2015).
37. U. Briskot, M. Schutt, I. V. Gornyi, M. Titov, B. N. Narozhny, and A. D. Mirlin, *Phys. Rev. B* **92**, 115426 (2015).
38. V. V. Cheianov, V. I. Fal'ko, B. L. Altshuler, and I. L. Aleiner, *Phys. Rev. Lett.* **99**, 176801 (2007).
39. J. Martin, N. Akerman, G. Ulbricht, T. Lohmann, J. H. Smet, K. von Klitzing, and A. Yacoby, *Nat. Phys.* **4**, 144 (2008).
40. B. Wunsch, T. Stauber, F. Sols, and F. Guinea, *New J. Phys.* **8**, 318 (2006).
41. E. H. Hwang and S. Das Sarma, *Phys. Rev. B* **75**, 205418 (2007).
42. M. R. Ramezanali, M. M. Vazifeh, R. Asgari, M. Polini, and A. H. MacDonald, *J. Phys. A: Math. Theor.* **42**, 214015 (2009).
43. We assume that the out-of-plane phonons are suppressed since the graphene sample in Ref. [31] is encapsulated between two boron-nitride substrates.
44. J. L. Manes, *Phys. Rev. B* **76**, 045430 (2007).
45. D. M. Basko, *Phys. Rev. B* **78**, 125418 (2008).
46. The acoustic-phonon drag effect should be also negligible compared to the direct electronic component due to phonon-phonon interactions (see Refs. [22] and [23]).
47. G. D. Mahan, *Many-Particle Physics*, 3rd ed. (Kluwer Academic/Plenum, New York, 2000), p. 177.
48. J. Rammer, *Quantum Field Theory of Non-Equilibrium States* (Cambridge University Press, New York, 2007); J. Rammer and H. Smith, *Rev. Mod. Phys.* **58**, 323 (2007).
49. We argue that no valley or spin imbalance develops, so that the distribution functions remain diagonal in valley and spin spaces. Both elastic impurity scattering and inelastic Coulomb scattering preserve the valley and spin symmetry. Although optical phonons involve intervalley scattering, a non-valley-diagonal component of the distribution function is not expected.
50. Q. Shao, G. Liu, D. Teweldebrhan, and A. A. Balandin, *Appl. Phys. Lett.* **92**, 202108 (2008).
51. F. T. Vasko and V. Ryzhii, *Phys. Rev. B* **76**, 233404 (2007).
52. The nonmonotonicity in doping at low temperature could be an artifact of the semiclassical Boltzmann equation and might be smeared out by multiple scattering processes that are not captured by the Boltzmann equation approach.
53. D. M. Basko and I. L. Aleiner, *Phys. Rev. B* **77**, 041409 (2008).
54. We estimate \((M_{opt})_{00} \sim \alpha_{opt}^{2}(T_{F}/T)^{4}e^{-T_{F}/T}\) and \((M_{imp})_{00} \sim (g + \pi\alpha_{int}^{-2}v_{F}^{2}n_{imp}/T_{F}^{2})(T_{F}/T)^{4}\), where \(g\), \(\alpha_{int}\), and \(n_{imp}\) take the values extracted from experiment [31].
55. G. Catelani and I. L. Aleiner, *Zh. Eksp. Teor. Fiz.* **127**, 372 (2005) [JETP **100**, 331 (2005)]; G. Catelani, *Phys. Rev. B* **75**, 024208 (2007).
56. K. Flensberg and Ben Yu-Kuang Hu, *Phys. Rev. B* **52**, 14796 (1995).
57. L. S. Levitov, A. V. Shtyk, and M. V. Feigelman, *Phys. Rev. B* **88**, 235403 (2013).
58. F. J. G. de Abajo, *ACS Photonics* **1**, 135 (2014).
59. E. M. Lifshitz and L. P. Pitaevskii, *Physical Kinetics*, Course of Theoretical Physics Vol. 10 (Elsevier, Singapore, 2008).
60. P. B. Allen, *Phys. Rev. B* **17**, 3725 (1978).
61. H. Fehske, R. Schneider, and A. Weisse (Eds.), *Computational Many-Particle Physics*, Lecture Notes in Physics Vol. 739 (Springer, Berlin, Heidelberg 2008), Chap. 8.
62. G. E. Uhlenbeck, G. W. Ford, and E. W. Montroll, *Lectures in Statistical Mechanics* (American Mathematical Society, Providence, 1963).
63. L. D. Landau and E. M. Lifshitz, *Statistical Physics Part 1*, Course of Theoretical Physics Vol. 5, 3rd ed. (Elsevier, Singapore, 2007), p. 178.
64. See the diagram (b)i in Fig. 1. If we artificially suppress channel A, one should obtain \(\tau_{ee}^{-1} \sim e^{-\mu/k_{B}T}k_{B}T^{-1}\alpha_{int}^{2}\) for \(\mu \ll k_{B}T\). Here we assume finite screening.
65. S. Sachdev, *Phys. Rev. B* **57**, 7157 (1998).
66. M. Polini and G. Vignale, *arXiv:1404.5728* [cond-mat.mes-hall].