s41467-021-25864-1.pdf

Type: Document | Status: ready
← All sources

ARTICLE Carrier transport theory for twisted bilayer graphene in the metallic regime Gargee Sharma 1,2,10, Indra Yudhistira 1,3,10, Nilotpal Chakraborty 1,4,5, Derek Y. H. Ho 1,4, M. M. Al Ezzi 1,3, Michael S. Fuhrer 6,7, Giovanni Vignale 1,4,8 & Shaf fique Adam 1,3,4,9 ✉ Understanding the normal-metal state transport in twisted bilayer graphene near magic angle is of fundamental importance as it provides insights into the mechanisms responsible for the observed strongly correlated insulating and superconducting phases. Here we provide a rigorous theory for phonon-dominated transport in twisted bilayer graphene describing its unusual signatures in the resistivity (including the variation with electron density, tempera- ture, and twist angle) showing good quantitative agreement with recent experiments. We contrast this with the alternative Planckian dissipation mechanism that we show is incom- patible with available experimental data. An accurate treatment of the electron-phonon scattering requires us to go well beyond the usual treatment, including both intraband and interband processes, considering the finite-temperature dynamical screening of the electron- phonon matrix element, and going beyond the linear Dirac dispersion. In addition to explaining the observations in currently available experimental data, we make concrete predictions that can be tested in ongoing experiments. https://doi.org/10.1038/s41467-021-25864-1 OPEN 1 Centre for Advanced 2D Materials, National University of Singapore, 6 Science Drive 2, 117546 Singapore, Singapore. 2 School of Basic Sciences, Indian Institute of Technology Mandi, Mandi 175005, India. 3 Department of Physics, National University of Singapore, 2 Science Drive 3, 117551 Singapore, Singapore. 4 Yale-NUS College, 16 College Avenue West, 138527 Singapore, Singapore. 5 Max-Planck-Institut für Physik komplexer Systeme, Nöthnitzer Straße 38, Dresden 01187, Germany. 6 ARC Centre of Excellence in Future Low Energy Electronic Technologies, Monash University, Monash, VIC 3800, Australia. 7 School of Physics and Astronomy, Monash University, Monash, VIC 3800, Australia. 8 Department of Physics and Astronomy, University of Missouri, Columbia, MO 65211, USA. 9 Department of Materials Science and Engineering, National University of Singapore, 9 Engineering Drive 1, 117575 Singapore, Singapore. 10These authors contributed equally: Gargee Sharma, Indra Yudhistira. ✉email: shaffi[email protected] NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications 1 1234567890():,; T he seminal observation of superconductivity and correlated insulating states in twisted bilayer graphene (tBG) 1–3 has generated tremendous excitement in the physics community4. At present, there is no consensus on the mechanism responsible for these observations. It was anticipated almost 15 years ago5 that when two sheets of graphene are stacked on top of each other with a slight relative rotation, a large wavelength moiré superlattice potential emerges. By reducing the twist angle in these moiré systems, the Bloch period can be increased by two orders of magnitude thereby bridging the lengthscales between naturally occurring lattices in materials and optical traps of cold atoms6. The addition of a moiré potential signi ficantly modi fies the underlying electronic structure, including both a reduction in the Fermi velocity at low energy and a reduction of the bandwidth of the lowest energy band. Both these effects enhance the impor- tance of electron –electron interactions 7. These properties can be understood as follows: In the absence of any coupling between the layers, the original Dirac-like bands are just folded onto the smaller moiré Brillouin zone as determined by sym- metry, but not modi fied; it is the interlayer coupling that causes level repulsion between the folded moiré bands. The moiré band closest to charge neutrality remains Dirac-like at low energy (the sublattice symmetry protecting the Dirac cones is not broken by the moiré potential) but with a reduced Fermi velocity, and the first moiré band as a whole gets squeezed by the level repulsion. This reduced bandwidth is quanti fied by the separation in the energy of the two van Hove singularities (VHS) that are found at the midpoint between the two original (relatively rotated in the Brillouin zone) Dirac cones, but pushed closer in energy by the level repulsion. Numerical ab initio studies soon con firmed the predictions of this long-wavelength continuum picture 8; however, the experimental situation remained controversial for a while (see e.g., ref. 9). Since then, the continuum model has been largely con firmed experimentally (see e.g., refs. 10,11). Taking the continuum model to its logical conclusion, Bis- tritzer and MacDonald predicted 12 that the Fermi velocity would vanish at a family of so-called “magic angles”. Their original work assumed that the lattices remained rigid. More recent work including lattice relaxation effects 13,14 suggests that only the first and largest magic angle ( θM~1.06∘) is stable, and that the rigid lattice continuum approximation breaks down for smaller angles. It should be emphasized that within the continuum model, strictly speaking, the bandwidth or 2 εVHS remains finite at the magic angle. However, experimentally, at least in local spectro- scopy measurements (e.g., refs. 15–18), an alternate de finition of magic angle is possible, i.e., when εVHS = 0. These would occur at angles below the original magic angle and in the regime where lattice relaxation effects are dominant (and it is not clear, in this case, what the electronic structure would look like). Given the observation of strongly correlated physics in other twisted 2D materials19–22, it seems that the vanishing bandwidth is more germane than the vanishing Fermi velocity, although, at present, the relation between the two has not been established. In this work, we establish yet another special angle, θcr ~1.15∘, the angle at which the Fermi velocity equals the phonon velocity. We show that at this angle, the phonon contribution to the resistivity strictly vanishes, and the experimentally measured resistivity would increase by several orders of magnitude for small deviations in angle on either side of θcr. It has become normative in this quickly evolving field to attribute factor of ~5 changes in the resisitivity 23 as evidence for superconductivity, and our work suggests more caution. By construction, each magic angle must be accompanied by two critical angles θcr (above and below θM), and therefore its effects should be robust to lattice relaxation effects (as we demonstrate explicitly below). We demonstrate that the Fermi velocity of the linear bands and the VHS at the edges of the moiré Brillouin zone have distinct effects on the resistivity, and these could therefore be used in transport experiments to dis- entangle the importance of each in the correlated regime. The present work is not about the observed superconductivity or correlated insulators. As we explain here, there is a geometric enhancement of the electron –phonon coupling in such moiré systems24,25 that would favor a phonon mechanism for super- conductivity; however, in a separate paper 26 we show that plas- mons are also strongly enhanced and that superconductivity can arise from a purely electronic mechanism. Similarly, at present, it is unclear if the correlated insulator is a Mott insulator (see e.g., refs. 27–35) or a Wigner crystal 36,37. Furthermore, it is possible for long-range interactions to signi ficantly distort the noninteracting bands away from charge neutrality due to the formation of inhomogeneous electrostatic potentials (see e.g., refs. 38,39) although the experiments seem to suggest otherwise. Rather, this work is about the carrier transport theory in the metallic regime (including at the van Hove singularity at higher carrier densities). We find that the role of phonons in tBG is perhaps as interesting as that of electrons: the same moiré potential that gives rise to the flat electronic bands, also results in enhancement of the electron–phonon coupling. Soon after the first experiments1–3, we predicted 40 that charged impurities would always dominate the resistivity at the lowest carrier densities and temperatures, but that gauge phonons would dominate for most of the experimental window. This crossover is also present in monolayer graphene, but occurs at a temperature of ~ 500 K, while for tBG the crossover happens at ~5 K. As we show in Methods, available experimental data largely confirm our earlier predictions. By now there have been two experimental transport studies focusing on the metallic regime. The first is from the MIT group41 and the second is a UCSB-Columbia collaboration 42. While the two experiments are largely consistent with each other, they arrive at very different conclusions on the dominant scat- tering mechanisms at play. Ref. 41 argues for a Planckian mechanism to explain their data, which implies a scattering rate ℏτ−1 = CkBT, where C ≲ 143. Here C = 1 is the Planckian bound set by holography and believed to be relevant for strange metals 44. They argue that the linear-in-temperature behavior persisting well below the Bloch –Gruneisen temperature and the saturation of resistivity at higher temperature are both inconsistent with the 05 0 1 0 0 0 0.05 0.1 2468 0.07 0.075 0.08 0.085 Fig. 1 Illustration of the different carrier transport regimes in twisted bilayer graphene. Electron–phonon scattering dominates except at very low temperature. At intermediate temperature ( TF, TBG) ≪ T < TVHS, resistivity is linear-in- T and at high temperature ( T ∼ εVHS), resistivity saturates as a function of temperature. Inset: electron-impurity scattering dominates the transport at low temperature. Resistivity increases for T ≪ TF and decreases for T ≳ TF. This nonmonotonicity arises from the charged impurity scattering 40, and when the system crosses over from the degenerate to the nondegenerate regime. ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 2 NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications conventional theory of phonon transport. We show here that both of these features are actually essential features of phonon- limited transport in tBG (see Fig. 1). While a microscopic theory showing Planckian dissipation has not yet been developed for tBG (it has for other systems, see e.g., ref. 45,46), we can assume a Planckian mechanism to make predictions for the transport. We show here that a Planckian mechanism also gives a saturation in resistivity at high temperature and linear-in-temperature beha- vior at a low temperature consistent with experimental observa- tions. In this work, we focus on the metallic regime that is far from the superconducting regime. A detailed analysis ultimately shows that the dominant scattering mechanism in this regime is not Planckian for several reasons including (a) the Planckian theory also predicts a strong carrier density dependence (absent in the experiment); (b) the experiment and the phonon mechanism both show the resistivity saturation at high tem- perature is set by the VHS energy, while for the Planckian theory this saturation is intrinsic (i.e., independent of bandstructure); (c) the twist angle dependence of Fermi velocity as extracted from experiment for phonon-limited scattering is consistent with the continuum theory 12–14, while it is orders of magnitude off for the Planckian theory; and most signi ficantly, (d) the extracted value of the scattering time from the experiment using the Planckian theory contradicts the assumptions of the Planckian theory. Our work shows that the phonon interpretation of ref. 42 is consistent with the theory we develop here. Results and discussion Comparison of resistivities . In Fig. 2 we compare data from ref. 42 (middle panels) with both a phonon-limited theory (left panel) and a Planckian theory (right panel). Similar to the experimental data, the phonon-mediated theory has weak density dependence. By contrast, the resistivity of the Planckian theory has strong density dependence (not seen in the experiment) that results from the density of states dependence of the Drude weight, which, unlike electron –phonon, remains uncompensated by the scattering time. We note that both the phonon-limited theory and the Planckian theory are linear-in- T at low temperature, and saturate at high temperature (qualitatively similar to what is seen experimentally). However, the origin of the saturation is very different. For phonon scattering, the saturation is set by the electronic bandwidth 2 εVHS, while for Planckian dissipation it is mostly independent of εVHS and set by Planckian strength C, which is expected to be somewhat universal and C ≤ 1. This illustrates that both the phonon-limited theory and the Planckian theory provide robust predictions that can be tested against experiment. Near magic angle, tBG has electron-hole asymmetry that arises from the second-nearest hoping in the effective tight- binding model for graphene. The superlattice potential renor- malizes the kinetic energy scales including the asymmetry to lower energies 12. We include this effect in the theory by fitting separately for the electron and hole side (see below). Details of how we fit the experimental data to electron –phonon and Planckian theory are provided in Methods. Linear-in-T resistivity. Since phonons appear to dominate the transport properties, we are motivated to carefully consider the role of electron –phonon scattering. In this paper, we inves- tigate the problem of normal state electronic transport in tBG focusing on the role of electron –phonon collision, which we find is the most important scattering mechanism in the relevant temperature and density regimes. The usual treatment of the electron–phonon does not take into account interband scattering, that we find below to be crucial near the magic angle. Secondly, dynamical screening of phonons is completely neglected because typically v F ≫ cph, which again no longer holds true near the magic angle. Thirdly, geometric enhancement of the gauge pho- non mode remains poorly addressed, which we find below to be the most dominant phonon mode in tBG. Lastly, the usual treatment of the electron –phonon problem has limited validity and fails beyond the linear regime (near VHS) that pushes us to go beyond the Dirac approximation and include nonlinear lattice effects. An accurate treatment of the electron –phonon scattering requires us to go well beyond the usual treatment, whereby including both interband and intraband processes, we show, for example, that the interband process allows for a linear-in- T behavior well below the Bloch –Gruneisen temperature and the transition between the two is accompanied by several orders of magnitude decrease in the resistivity at a critical angle θ cr, distinct from the magic angle θM. By considering the finite-temperature dynamical screening of the electron –phonon matrix element, we show, for example, that only the antisymmetric gauge phonon mode survives at a low twist angle; and by going beyond the linear Dirac dispersion, we show that the van Hove singularity causes saturation in resistivity as a function of temperature. In addition to explaining the observations in currently available experimental data, our theory makes concrete predictions that can be tested in ongoing experiments. To our knowledge, the Boltzmann transport theory for acoustic phonon scattering in monolayer graphene was first developed by Hwang and Das Sarma 47. In this work, we adopt the same formalism with several extensions appropriate for tBG. First, we consider both intraband and interband processes. As we show in the Supplementary Information, close to the magic angle, only interband scattering is operational which was not considered in ref. 47. Second, while the linear Dirac Hamiltonian was an appropriate model for monolayer graphene, the reduced energies in tBG requires us to use an effective two-band Hamiltonian first proposed by ref. 48 that captures the physics near the van Hove singularity. Finally, we do the full finite-temperature and finite frequency RPA screening of the electron –phonon matrix elements (which is necessary due to the diverging density of states close at the magic angle). This demonstrates that it is the off-diagonal (or so-called gauge phonon contributions 40,49,50)o f the acoustic phonon matrix element that dominates the transport properties. How to screen the electron –phonon matrix element in two dimensions has long remained controversial. The issue is that prior to the present work, calculating the full dynamical polarizability at finite temperature has been challenging. Since the phonon propagator couples at a particular frequency, without dynamical screening, it is unclear how the electrons screen the deformation potential. This led to speculation in the theoretical literature as to whether the deformation potential should be screened or left unscreened. For example, in ref. 51, Okuyuma and Tokuda argue that experimental data for GaAs 2DEGs is better fit using the unscreened theory, while in a later work, Kawamura and Das Sarma argue that once correctly done, static screening gives excellent agreement with experimental data 52. Even for monolayer graphene, ref. 50 argue that the deformation potential is completely screened, while the ref. 53 argue for no screening. In this work, we demonstrate conclusively that as anticipated by Kawamura and Das Sarma, the static screening approximation is closer to the correct dynamically screened result than the commonly used unscreened approximation. The Supplementary Information has details of the screening as well a discussion of the geometric enhancement due to the moiré that is present for gauge phonons, but not scalar phonons. Interplay between intraband and interband scattering . There are two qualitatively distinct regimes depending on whether the NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 ARTICLE NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications 3 phonon velocity cph is greater than or smaller than the Fermi velocity vF. The crossover from vF > cph to vF < cph is expected because the renormalized vF vanishes at magic angles θM. For the largest magic angle, this crossover occurs at θ = θcr ~1.15∘and separates the regimes of interband ( vF < cph) and intraband (vF > cph) scattering. Within the Dirac regime, the theory for intraband scattering is now well established 47,50,53–55. The resis- tivity shows a Bloch –Grüneisen behavior similar to metals and is given by ρe/C0ph ¼½ 16ζðθÞ2kF=ðe2μscphv2 FÞ/C138FðTBG=TÞ, where TBG = 2ℏcphkF is the Bloch –Grüneisen temperature, which is the char- acteristic crossover scale over which the temperature dependence of the resistivity due to electron –phonon scattering changes from T4 below TBG to T − linear above TBG. This change in the resistivity occurs due to the restricted scattering phase space of phonons at low temperatures when their quantum nature becomes important, compared to higher temperatures where the phonon distribution is quasiclassical. The graphene mass density is μs, and FðxÞ¼ R 1 0 dy½xy4 ffiffiffiffiffiffiffiffiffiffiffiffiffi 1 /C0y2p exy/C138= exy /C01ðÞ 2 (This form of the integral first appeared in ref. 53). For T ≫ TBG, the quanti- zation of the lattice phonon modes is irrelevant and the scattering is expected to be proportional to the amplitude of lattice vibrations and is linear-in-T. In fact, one can show exactly that the linearity persists up to temperatures as low as TBG/4. The interband phonon scattering rate (dominant close to magic angle) shares some similarities with the intraband scattering: it is density-independent and T − linear at high- T. Moreover, it vanishes when vF → 0 as the scattering phase space tends to zero. However, qualitatively the interband and intraband scattering are quite different (see Supplementary Information). For example, while the intraband scattering rate within the Dirac model shows a monotonic increase with energy, the interband scattering rate is non-monotonic highlighting the suppression of interband scattering for energies larger than kBT. Most impor- tantly, the temperature scale for the scattering rate to be T − linear is not set by the Bloch –Grüneisen temperature TBG, but rather by the Fermi temperature TF (which is the maximum phonon energy at the Fermi surface allowed by kinematic constraints). We note that TF and TBG are de fined in such a way so that at θcr, TF = TBG/2. Therefore, close to the magic angle when vF < cph, we have TF < TBG/2 and the electron –phonon scattering becomes T − linear for temperatures well below TBG The reversal of the temperature scales is due to the fact the TF depends on the twist angle dependence (through vF). On the -20 0 20 0 0.1 0.2 0.3 0.4 0.5 -20 0 20 -20 0 20 0 50 100 0 0.1 0.2 0.3 0.4 0.5 0 50 100 0 50 100 Fig. 2 Experimental observations are consistent with the electron-phonon scattering theory. The electron –phonon scattering theory (left panels) correctly captures the a carrier density and b temperature dependence of experimentally observed resistivity (middle panels), unlike the Planckian theory (right panel) that shows a stronger density dependence. Experimental data is taken from ref. 42 for θ = 1.11∘ (comparison for devices with other twist angles is shown in the Supplementary Information). Solid lines in the electron –phonon and Planckian theory are for a two-band effective model that includes the van Hove singularity, while the dashed lines are for the linear Dirac model. For electron –phonon scattering, the linear-in- T resistivity at low temperature is captured by the Dirac model, while the saturation at higher temperature requires the van Hove singularity. For the Planckian theory, the Dirac model a nd the two-band model are quantitatively similar and show much stronger density dependence compared to the experiment. In this case, the saturation at h igh temperature is set not by the van Hove singularity, but by a universal value ρðT !1 Þ¼ C=ð8ln 2Þh=e2 (the coefficient C ≤ 1 for Planckian dissipation). For most experimental data, including those shown here, C ≥ 1. Taken together with the weak density dependence seen experimentally, this suggests that phonon scattering rather than Planckian dissipation is the dominant scattering mechanism at play in twisted bilayer graphene. ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 4 NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications other hand, TBG depends on phonon velocity and density, both of which do not change with twist angle. Thus at small twist angles, TF can become signi ficantly smaller than TBG. For T ≫ TF,w e find ρe/C0ph inter ¼ h e2 2~β 2 Av2 FkBT _2μsc6 ph ; ð1Þ where ~βA is the twist angle-dependent enhanced gaugefield coupling constant. Close to magic angle, we expect ρe/C0ph inter ðT /C29 TFÞ/ v4 F (where the additional v2 F comes from ~β 2 A, see Supplementary Information), and vanishes at the magic angle due to the lack of scattering phase space. Figure 3a shows intraband resistivity for a chosen θ > θcr. The linear-in-T resistivity persists down to TBG/4 consistent with earlier results47. A comparison of the scales of TBG and TF is also done. Figure 3b shows the interband resistivity for a chosen θ < θcr, comparing the scales of TF and TBG. The linearity in T is observed to persist down to very low temperatures even when T ≪ TBG, which is very different from the known theory of electron–phonon scattering in a typical Fermi liquid. The empirical observation of T − linear resistivity well below TBG has been attributed to the strange metallicity of non-Fermi liquids41, however, we find that there is nothing mysterious about this feature, it is merely the qualitative change in the nature of electron–phonon scattering when vF < cph i.e., the lower of the two energy scales switches from lattice vibrational energy to electronic energy when vF crosses below cph. Lattice relaxation effects . For small twist angles, atoms on both layers will tend to move away from their nominal positions in order to minimize the total energy of the system 13. This relaxa- tion process will increase the fraction of atoms with AB stacking relative to regions with AA stacking. This rearrangement of atoms also changes the relative strength of moiré coupling between different sublattices across the twisted interface. Relaxed atomic positions are calculated either using a continuum elasticity theory56 for what we call “relaxed model A ” or using the mole- cular dynamic approach as implemented in LAMMPS 57,58 that we call “relaxed model B ”. Once we know the relaxed atomic positions, the moiré coupling parameters are obtained by Fourier 012 0 2 4 6 10-5 05 1 0 0 5 10 15 10-3 (a) (b) Fig. 3 The persistence of linear-T behaviour to low temperature is due to inter-band scattering. Electron–phonon resistivity for tBG within the Dirac model for a intraband scattering and b interband scattering. Dotted and dashed lines indicate TBG/4 and TF/2, respectively. Interband resistivity shows a transition to linear-in- T at ∼TF/2 compared to intraband resistivity which shows a transition around ∼TBG/4. In the interband regime, the striking persistence of linear-in- T behavior for T ≪ TBG is observed. 0.6 0.8 1 1.2 0 0.1 0.2 0.3 0.6 0.8 1 1.2 0 0.1 0.2 0.3 0.6 0.8 1 1.2 0 0.1 0.2 0.3 Fig. 4 Close to the magic angle, the electron –phonon resistivity in tBG is very sensitive to the twist angle exhibiting a variation of several orders of magnitude when vF = cph and unrelated to Mott insulation or superconductivity. The number of these sharp dips in resistivity and the angles at which they occur provide information about lattice relaxation 13. a The rigid lattice model of Bistritzer and MacDonald 12 predicts that the resistivity will have multiple dips with decreasing twist angle (corresponding to three dips per magic angle). b The relaxation model of ref. 64 gives only three sharp dips close to a single value of the magic angle. c The relaxation model of ref. 58 gives three sharp dips close to a single value of the magic angle and another dip when vF becomes quite close to cph but never goes below it. NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 ARTICLE NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications 5 transforming the matrix that couples the orbitals in layer 1 and layer 2. It is these “relaxed” hopping parameters that are then used in the continuum model Hamiltonian to calculate electronic bandstructure and corresponding renormalized Fermi velocity (see Fig. 4) that is used as the input for our Boltzmann transport calculation. While the two relaxation models have some quanti- tative differences such as the position of the magic angle, they both give the same qualitative description for the role of relaxa- tion on the moiré bandstructure. In Fig. 4 we plot the electron –phonon resistivity for both the rigid continuum model 12 of tBG as well as including the lattice relaxation effects using the two different models that we label as relaxed model A and relaxed model B. We first note that whenever vF = cph, there are sharp dips in the resistivity pro file, which can span a few orders of magnitude. Secondly, at each magic angle, there is another large dip. For the rigid lattice model, there is a broad window where vF < cph around the magic angle θM ~1.06∘, and the family of magic angles implies multiple crossings of vF = cph, and hence multiple dips in the resistivity as the twist angle is lowered. With the inclusion of lattice relaxation effects, there is only a single stable magic angle (that is shifted slightly from the original magic angle), and there are therefore only three dips in the resistivity near each magic angle. The window of the interband scattering regime ( v F < cph) is also narrower. We find that for relaxed model B, near 0.6 ∘, the Fermi velocity becomes quite close to cph but never goes below it. This results in another sharp dip near that angle. We emphasize that the following robust features survive with or without the inclusion of relaxation effects and are independent of the relaxation modeling: (i) existence of the critical angle θcr by construction, (ii) huge drop in resistivity at θcr and at θM resulting in a total of three sharp resistivity dips near each magic angle. Since close to magic angle, the resistivity is highly sensitive to twist angle exhibiting a variation of several orders of magnitude, one must be careful to experimentally distinguish this phonon effect from possible Mott insulation or superconductivity. Beyond the Dirac model . While the Dirac model captures the physics of tBG at low density, in order to make accurate pre- dictions for larger density and temperature, we need to extend the model to capture the VHS. We use an effective two-band Hamiltonian48 HðkÞ¼/C0 _vF jΔKj 0 k/C32 /C0ðΔK/C3=2Þ2 k2 /C0ðΔK=2Þ2 0 ! ; ð2Þ where vF is the Fermi velocity of tBG at the Dirac point, k = kx + iky, ΔK is wave vector separation between the two Dirac points which are located at K and Kθ, which magnitude is given by kθ /C17jΔKj¼ 2kD sinðθ=2Þ, with kD being the wave vector separation between the two Dirac points in monolayer graphene and θ is twist angle. The two-band model is valid when the bandwidth is much larger than ℏvF∣ΔK∣ which is a good approximation for small twist angles. The eigenenergies of this Hamiltonian are given by ε k;λ ¼ λð1=4Þ_vF ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi k2 θ þ 8ðk2 x /C0k2 y Þþ 16½ðk2 x þ k2 y Þ=kθ/C138 2 q ,a n d are anisotropic. Within an isotropic approximation for the momentum-transport cross section, we calculate the electron–phonon scattering rate for both intraband ( vF > cph)a n d interband ( vF < cph) with only intravalley electron –phonon scat- tering (see Supplementary Information for details). We find where zν = cν/vF, qðr; ϕ; r0; ϕ0Þ is momentum transferred by phonon, and f 0 p0 ;λðnq;ν Þ is the equilibrium Fermi (Bose –Einstein) distribution function describing the electron (phonon) popula- tion. Here, ξ = ±1 refers to phonon absorption ( +1) and emission (−1), respectively. The resistivity for this Hamiltonian is 1 ρij ¼ 8e2 ∑ λ¼ ±1 R 1 0 dr R π /C0π dϕ J ðr;ϕÞ ð2πÞ2 vðiÞ k;λvðjÞ k;λ ´ τe/C0phðr; ϕÞ/C0 ∂f 0 k;λ ∂εk;λ /C16/C17 ð5Þ where vðjÞ k;λ ¼ð 1=_Þð∂εk;λ=∂kjÞ is band velocity in j direction and λ = ±1 refers to conduction ( +1) and valence ( −1) band, respectively (see Supplementary Information for the explicit form of J ðr; ϕÞ). We plot the results of our calculation of hρi¼ ffiffiffiffiffiffiffiffiffiffiffiffiρxxρyy p for twist angle of θ = 1.1∘ in the left panel of Fig. 2. The electron–phonon resistivity is linear-in- T at low temperature but saturates at high temperature. The slope of resistivity with temperature in the low temperature regime is set by vF while the saturation of resistivity with temperature is set by the bandwidth. Within this effective Hamiltonian, the bandwidth 2 εVHS and vF are not independent and related by 2 εVHS = (1/2)ℏvFkθ. While this relation between εVHS and vF is speci fic to our model, we believe that the main conclusions are generic i.e., the low temperature linear-in- T behavior is set by the Fermi velocity, while the saturation is set by the VHS. The saturation of electron–phonon resistivity can be simpli fied deep in the intraband ( vF ≫ cph) and interband ( vF ≪ cph) regime as ρe/C0ph intra ðT !1 Þ/ ~β 2 A e2_μs 1 v2 F c2 ph εVHS ð6Þ ρe/C0ph inter ðT !1 Þ/ ~β 2 A e2_μs v2 F c6 ph εVHS ð7Þ Resistivity from the Planckian model . We take a simple phe- nomenological model with _τ/C01 Pl ¼ Ck BT. The resistivity is 1 τe/C0ph intra ðr;ϕ;λÞ ¼ ∑ ξ¼ ±1 ν¼TA;LA ~β 2 Akθ 8πμs _cν vF Z 1 0 dr0r0 Z π /C0π dϕ0 qðr; ϕ; r0; ϕ0Þsin2 ϕ0 /C0ϕ /C0/C1 ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1 þ r0 cos ϕ0Þ2 þð r0 sin ϕ0Þ2 q ð1 /C0ξÞ=2 þ ξf 0 p0 ;λ þ nq;ν /C16/C17 ´ δλ ðr0 /C0rÞ/C0 4ξz qðr;ϕ;r0 ;ϕ0 Þ kθ /C16/C17 ð3Þ 1 τe/C0ph inter ðr;ϕ;λÞ ¼ ∑ ν¼TA;LA ~β 2 Akθ 8πμs_cν vF Z 1 0 dr0 Z π /C0π dϕ0 qðr; ϕ; r0; ϕ0Þsin2 ϕ0 /C0ϕ /C0/C1 ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1 þ r0 cos ϕ0Þ2 þð r0 sin ϕ0Þ2 q ð1 þ λÞ=2 /C0λf 0 p0 ;/C0λ þ nq;ν /C16/C17 ´ δ r0 þ r /C04z qðr;ϕ;r0 ;ϕ0 Þ kθ /C16/C17 ; ð4Þ ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 6 NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications obtained from the Boltzmann equation by using the relaxation time approximation, with τPl as the relaxation time (see discus- sion in ref. 43), and calculating the appropriate thermal average with the density of states and Fermi velocity (see Eq. 19). We can do this both for the Dirac model and the two-band Hamiltonian. Surprisingly, both the Dirac Hamiltonian and the two-band model give very similar results indicating that the van Hove singularity is not important for the Planckian theory. This phe- nomenological model exhibits a linear-in- T resistivity at low temperature that saturates at higher temperature (qualitatively similar to what is seen experimentally). We find that the slope of resistivity with temperature in the low temperature regime is set by both vF and C, while the saturation of resistivity with tem- perature is set only by C. Within the Dirac approximation, the resistivity is 1 ρPl ¼ 1 C 4e2 h ∑ λ¼ ±1 ln 1 þ exp λ μ kBT /C18/C19/C20/C21 ; ð8Þ with low and high temperature asymptotes ρPl ¼ C h 4e2 T TF 1 þ π2 6 T TF /C16/C17 2/C20/C21 ; T /C28 TF C h e2 1 8ln 2 1 /C0 1 128 ln 2ðÞ 3 TF T /C16/C17 4/C20/C21 ; T /C29 TF 8

< : ð9Þ For the effective two-band model, the Planckian resistivity is given by an expression similar to Eq. 5. It can be simpli fied to 1 ρij Pl ¼ e2 h 1 C Kj n nVHS ; kBT εVHS /C18/C19 ; ð10Þ where the function Kj is computed numerically (see Supplementary Information). In the right panel of Fig. 2a, b, we show the Planckian resistivity as a function of density and temperature, respectively at twist angle of 1.11 ∘. The Planckian resistivity for both models saturate at ρPlðT !1 Þ¼ C=8ln 2, independent of the bandwidth εVHS, which is in sharp contrast to electron–phonon scattering. Figure 5 shows the results of fitting 1 1.2 1.4 1.6 1.8 2 0 5 10 15 1 1.2 1.4 1.6 1.8 2 0 2 4 6 8 10 1 1.2 1.4 1.6 1.8 2 0 2 4 6 8 10 ( ): Relaxed model A, electron (holes) ( ): Relaxed model B, electron (holes) 1 1.2 1.4 1.6 1.8 2 0 5 10 15 (a) (b) (d)(c) Fig. 5 The parameters obtained from fit to refs. 2, 42, 62 (see Methods for details). We find that the data agree much better with the expectations from the electron –phonon theory than the Planckian model. a Relative deviations of Fermi velocity as a function of twist angle obtained by fitting the experimental resistivity to the electron –phonon theory (b) obtained by fitting to the Planckian theory, c Relative deviations of the effective electron –phonon coupling constant ~βA, and d Planckian strength C. The error bars represent 95% con fidence interval of the parameters. The Planckian bound ( C ≤ 1) is violated for small twist angles rules out the Planckian model as the dominant transport mechanism in tBG. -20 0 20 0 0.1 0.2 0.3 -20 0 20 0 0.1 0.2 0.3(a) (b) Fig. 6 The weak density dependence at high temperature and strong density dependence at a low temperature of the resistivity close to charge neutrality seen in experiment 42 is well captured within the Dirac theory of electron –(gauge) phonon scattering 40. Experimental data (top squares) for resistivity at θ = 1.24∘ vs density at various temperatures ( a) compared to the Dirac theory (bottom solid lines) of electron –phonon scattering ( b). NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 ARTICLE NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications 7 eight data sets of varying twist angle, temperature, and carrier density to both the phonon and Planckian models (full fits are shown in Supplementary Information). We find that the phonon scattering theory (but not the Planckian model) gives fit parameters for both the Fermi velocity and model parameters that are consistent with theoretical expectations. The nature of electron –phonon scattering investigated here gives rise to several new features: Interband electron –phonon scattering, which should be kinematically forbidden in monolayer graphene, is shown to occur in tBG below a critical twist angle θcr, when vF < cph. The critical angle ( θcr) is a sweet spot where both the interband and intraband scattering phase space (and thus resistivity) both drop to zero, giving rise to multiple dips in the resistivity unrelated to insulating or superconducting states that could be investigated in the future experiments. We derive explicit analytical expressions for the interband scattering rate and show its qualitative dissimilarity from the intraband scattering. Importantly, we show that this explains the linear- in-T resistivity well below the Bloch –Gruneisen temperature TBG, a previously unexplained experimental puzzle. In this work, we also provide additional theoretical veri fication of our earlier claim 40 that the gauge phonon modes are enhanced by the moiré geometry and are not screened, while the scalar phonon modes have neither property. We also provide experi- mental verification (see Figs. 6, 7) of our previous predictions that charged impurity scattering takes over as the dominant scattering mechanism at very low temperatures (below 20 K) and low carrier densities (below ~10 11cm−2). Taken together with the present work on phonon scattering (that applies at high temperature and high carrier density), this now presents a complete theory for the carrier transport for tBG in the metallic regime. Methods Dynamical screening of electron –phonon coupling . The basic building block of RPA screening is the polarizability bubble ΠC(q, ω), which can be written in the most general form as ΠC ðq; iωmÞ¼ 1 β g A ∑ k;iωn Tr ½Gðk; iωnÞGðk þ q; iωn þ iωmÞ/C138; ð11Þ where the summation iωn is over the imaginary frequency, G(k, iω) is the Green ’s function, and the trace is over the sublattice degrees of freedom. Performing the trace and the Matsubara summation we obtain Π C ðq; iωÞ¼ lim η!0þ g A ∑ λ;λ0 ;k f 0 k;λ /C0f 0 kþq;λ0 _ω þ εk;λ /C0εkþq;λ0 þ iη Fλλ0 k;kþq ð12Þ where Fλλ0 k;kþq is the tBG chirality factor, which for Dirac model is given by ð1 þ λλ0 cos θk;kþqÞ=2 and f 0 k;λ ¼½ expfðεk;λ /C0μÞ=ðkBTÞg þ 1/C138/C01 is the equilibrium Fermi distribution function. The dielectric function ϵ(q, iωm) is given by the RPA sum- mation and is related to the basic pair bubble as ϵ(q, ω) = 1 − VqΠC(q, iω), where Vq = 2πe2/κq is the Fourier transform of the Coulomb potential, κ being the dielectric constant. We evaluate ΠC and ϵ(q, ω) for finite frequencies and temperatures, semi-analytically. To the best of our knowledge, this has not been done previously in the literature, at least in the context of phonons (see Supple- mentary Information). Boltzmann transport . The distribution function fk,λ is evaluated within the Boltzmann transport formalism, which reads /C0eE _ ∇k f k;λ ¼ St½f k;λ/C138 ð13Þ The collision integral can be written as St f k;λ hi ¼ ∑ k0 ;λ0 ;ν Pλ0 λ k0 k;ν f k0 ;λ0 ð1 /C0f k;λÞ /C0Pλλ0 kk0 ;ν f k;λð1 /C0f k0 ;λ0 Þ; ð14Þ where Pλ0 λ k0 k;ν is the scattering probability from state λ0; k0/C12/C12 /C11 to λ; kji within phonon branch ν (TA or LA), which is given by Pλλ0 k;kþq;ν ¼ 2π _ gλ;λ0 k;kþq;ν /C12/C12 /C12 /C12/C12 /C12 2 nq;ν δðεkþq;λ0 /C0εk;λ /C0_ωq;ν Þ h þð1 þ nq;ν Þδðεkþq;λ0 /C0εk;λ þ _ωq;ν Þ i ; ð15Þ where the two terms account for absorption and emission of phonons, nq,ν is the Bose–Einstein distribution function describing the phonon population, ℏωq,ν = ℏcνq is the phonon energy, and gλ;λ0 k;kþq;ν is the electron –phonon coupling, which can be expressed as gλ;λ0 k;kþq;ν ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi _ 2Aμs ωq;ν s Mλ;λ0 k;kþq ð16Þ where A is the area of the graphene layer, μs is the mass density, and Mλ;λ0 k;kþq is the matrix element for scattering between initial and final states, which is given by Mλ;λ0 k;kþq ¼ ζqF λλ0 k;kþq hi 1=2 ; ð17Þ where ζ is the effective deformation potential and Fλλ0 k;kþq is the tBG chirality factor. The effective deformation potential ζ stands for either the effective scalar potential ~DA or twice the effective gauge potential 2 ~βA (see Supplementary Information). We use ansatz f k;λ ¼ f 0 k;λ þ neEτk;λ cos θk vF ∂f 0 k;λ ∂εk;λ ð18Þ to obtain the transport scattering time τk,λ. Finally the resistivity ρe−ph is obtained from the scattering time τe/C0ph k;λ by 1 ρij e/C0ph ¼ e2 g/C3 A ∑ k;λ¼ ±1 vðiÞ k;λvðjÞ k;λτe/C0ph k;λ /C0 ∂f 0 k;λ ∂εk;λ ! ; ð19Þ where g* = g for Dirac model and g* = g/2 for two-band effective model ( g = 8i s degeneracy), vðjÞ k;λ ¼ð 1=_Þð∂εk;λ=∂kjÞ is band velocity in j direction, and λ = ±1 refers to conduction ( +1) and valence ( −1) band, respectively. Resistivity within the Dirac model .I nF i g . 6, we show the comparison of experimental resistivity vs density from ref. 42 to theoretical resistivity due to electron–phonon and electron-impurity scattering within Dirac model40.A tl o w temperatures, electron-impurities dominates over electron–phonon interaction40, therefore we only fit the experiment to the electron-impurity limited resistivity. These charged impurities also give rise to carrier density inhomogeneities59, which cures the 0 1 02 03 04 05 0 0 0.1 0.2 0.3 0 1 02 03 04 05 0 0 0.1 0.2 0.3 (a) (b) Fig. 7 The non-monotonic temperature dependence at low temperatures and linear-in- T dependence at intermediate temperatures seen in experiment is well captured within the Dirac theory of electron –(gauge) phonon scattering40. Resisitivity vs temperature from a ref. 42 and b Dirac theory40 at twist angle of θ = 1.24∘. The non-monotonicity at low temperatures is due to the crossover from impurity-dominated transport to phonon-dominated transport. ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 8 NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications otherwise divergent resistivity at the charge neutrality. We fixed the parameters from the effective medium theory (EMT) 60 fit to experimental resistivity at the lowest temperature, i.e., 2 K. In the fit, we have included a short-range scattering component of the conductivity σs61, charged impurity density nimp, and charge density fluctua- tions nrms.T h e first two parameters are used to calculate the Boltzmann-RPA con- ductivity σB[σs, nimp], which takes into account both the dominant scattering mechanism of screened Coulomb impurities and additional scattering mechanisms due to short-range scatterers, such as point defects and line defects. The parameter nrms enters through the EMT equations. We find that the 1.24 ∘ device of ref. 42 has charged impurity density nimp = 2.1 × 1011cm−2, short-ranged conductivity σs = 104e2/h and charge density fluctuations nrms = 3.6 × 1010cm−2. At a temperature of 10 K and density of ~10 11cm−2, electron–phonon interaction start to produce a noticeable influence on the resistivity. We observe the reversal of temperature dependence trend at 10 K that is followed by linear-in- T dependence at intermediate temperatures (see Fig.7) and weaker density dependence at higher temperatures in an experiment. These were all predictions we made in our previous work 40 that have been now shown experimentally to be correct. Fitting of experimental data . In Fig. 2, we compare the electron –phonon resis- tivity hρi¼ ffiffiffiffiffiffiffiffiffiffiffiffiρxx ρyy p from the two-band effective model (left panel) with the experiment from ref. 2,42,62 (middle panel) at a twist angle of θ = 1.1∘. The fitting parameters are Fermi velocity vF and the enhanced gauge field coupling constant ~βA , which are obtained from fitting the temperature dependence of resistivity at density n = 1011cm−2 for a fixed twist angle. The dependence of the relative deviations of vF and ~βA on twist angle are plotted in Fig. 5a, c, respectively. The same parameters are used to plot the temperature dependence of electron –phonon resistivity at higher density as well as its density dependence at several fixed temperatures (see the left panel of Fig. 2b, a, respectively). We fit the electron-side and hole side separately due to a slight asymmetry between them. We find that in the low temperature regime, our results coincide with the Dirac model as expected, since the Hamiltonian in Eq. 2 reduces to a Dirac Hamiltonian at low energies. However, in the high temperature regime, only two band effective model agrees with the experiment (see Fig. 2). We see a clear saturation of the resistivity at higher temperatures as well as a weak density dependence in the electron and the hole side, as seen in the experiment. Both of these features are not captured within the Dirac model. From the fitting of the electron–phonon theory to the experiment (see Fig. 5c), we obtain βA around 1 − 9 eV, which is in good agreement with the accepted values for the gauge field coupling constant in monolayer graphene 50,63.T h e fit values for the velocity ratio vF/v0 also agree well with theoretical estimates 12. Hence, we have developed a transport theory of tBG which explains all the salient features observed in the entire metallic regime (intermediate and high temperatures) of the experiment

Comparison with Planckian Resistivity . We also compared the Planckian resis- tivity to the same experiment in Fig. 2 (right panel and middle panel, respectively). The fitting parameters are Fermi velocity vF and Planckian strength C. The fitting parameters of the resistivity fit to both the electron –phonon and Planckian theory for all the data sets are shown in Table 1. These are then used to plot the tem- perature dependence of Planckian resistivity at higher densities, as well as its Table 1 Fitting parameters for the resistivity fit to the electron –phonon and Planckian theory. The parameters are Fermi velocity (vF) and electron –phonon gauge coupling ð~βAÞ for electron–phonon theory and Fermi velocity ( vF) and Planckian strength ( C) for Planckian theory. Twist angle Electron –phonon Planckian Ref veðhÞ F =v0 ~β eðhÞ A ðeVÞ veðhÞ F =v0 Ce(h) 1.06∘ 0.14 ± 0.01(0.17 ± 0.01) 47 ± 2(51 ± 3) 0.30 ± 0.03(0.35 ± 0.02) 2.6 ± 0.1(2.5 ± 0.1) 42 1.07∘ 0.15 ± 0.01(0.15 ± 0.01) 32 ± 1(32 ± 1) 0.30 ± 0.01(0.30 ± 0.01) 1.06 ± 0.03(1.06 ± 0.03) 62 1.08∘ 0.2 ± 0.1(0.2 ± 0.1) 110 ± 60(100 ± 60) 0.5 ± 0.2(0.4 ± 0.3) 9 ± 3(9 ± 3) 2 1.11∘ 0.13 ± 0.02(0.13 ± 0.01) 39 ± 3(39 ± 3) 0.3 ± 0.1(0.31 ± 0.05) 2.0 ± 0.1(2.0 ± 0.1) 42 1.24∘ 0.18 ± 0.03(0.17 ± 0.03) 41 ± 4(39 ± 4) 0.5 ± 0.1(0.5 ± 0.1) 1.7 ± 0.1(1.7 ± 0.1) 42 1.59∘ 0.426 ± 0.001(0.426 ± 0.001) 24.6 ± 0.1(24.6 ± 0.1) 5 ± 4(5 ± 4) 1 ± 1(1 ± 1) 42 2.02∘ 0.491 ± 0.002(0.491 ± 0.002) 26.6 ± 0.2(26.6 ± 0.2) 2 ± 1(2 ± 1) 0.5 ± 0.3(0.5 ± 0.3) 42 01 0 2 0 3 0 4 0 5 0 -0.6 -0.4 -0.2 0 0 1 02 03 04 0 -0.5 -0.4 -0.3 -0.2 -0.1 0 0 1 02 03 04 0 -0.4 -0.2 0 0 1 02 03 04 0 -0.5 -0.4 -0.3 -0.2 -0.1 0 (a) (b) (c) (d) Fig. 8 Alternative method for comparing the phonon theory (blue), Planckian theory (red), and experimental data (black) focusing only on the change in resistivity with carrier density. We take four representative experimental curves from the full data set used to make Fig. 5,i . e . ,a θ = 1.07∘ (ref. 62), b θ = 1.08∘ (ref. 2), c θ = 1.11∘ (ref. 42), and d θ = 1.24∘ (ref. 42). Theory curves are the same as before, except we use the best theoretical estimates for βA and C. Also included is the Planckian theory with effective medium theory to include the effect of disorder smearing (green). Similar to the conclusions in F ig. 5,t h i s method shows that the phonon theory has better agreement with the experiment. NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 ARTICLE NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications 9 dependence on density at several fixed temperatures (see the right panel of Fig. 2b, a, respectively). Although the Planckian theory exhibits similar linear-in- T beha- vior at low temperature and saturation at high temperature, its density dependence is much stronger than that of the experiment, especially at low temperatures (see Fig. 8). The Fermi velocity vF extracted from the fit is somewhat larger than the theoretical prediction (see Fig. 5). Moreover, we observe that the Planckian bound (C ≤ 1) is violated for a small twist angle (see Fig. 5c). This rules out the Planckian theory as the dominant mechanism of metallic transport in tBG. We illustrate this point in Fig. 5 where we compare the relative deviations between theory and experiment. In order to compare the various theories and experiment on equal footing, these deviations were determined by shifting the magic angle for the theoretical Fermi velocity such that they all have a common magic angle that minimized the relative deviations of the parameters from their theoretical expec- tations— the x-axis could be therefore be interpreted as twist angle with the magic angle taken to be θ M = 0.72∘. Data availability The data used to generate the figures are available at https://github.com/indrayudhistira/ tBG_metallic_phonon. Code availability The code used to generate the figures are also available at https://github.com/ indrayudhistira/tBG_metallic_phonon. Received: 24 June 2020; Accepted: 17 August 2021; References

  1. Cao, Y. et al. Unconventional superconductivity in magic-angle graphene superlattices. Nature 556, 43 (2018).
  2. Cao, Y. et al. Correlated insulator behaviour at half- filling in magic-angle graphene superlattices. Nature 556, 80 (2018).
  3. Yankowitz, M. et al. Tuning superconductivity in twisted bilayer graphene. Science 363, 1059 –1064 (2019).
  4. Gibney, E. How ‘magic angle ’ graphene is stirring up physics. Nature 565, 15–18 (2019).
  5. Lopes dos Santos, J. M. B., Peres, N. M. R. & Castro Neto, A. H. Graphene bilayer with a twist: electronic structure. Phys. Rev. Lett. 99, 256802 (2007).
  6. Tarruell, L., Greif, D., Uehlinger, T., Jotzu, G. & Esslinger, T. Creating, moving and merging dirac points with a fermi gas in a tunable honeycomb lattice. Nature 483, 302 –305 (2012).
  7. Tang, H.-K. et al. The role of electron-electron interactions in two- dimensional dirac fermions. Science 361, 570 –574 (2018).
  8. Trambly de Laissardiére, G., Mayou, D. & Magaud, L. Localization of dirac electrons in rotated graphene bilayers. Nano Lett. 10, 804 –808 (2010).
  9. Hicks, J. et al. Symmetry breaking in commensurate graphene rotational stacking: comparison of theory and experiment. Phys. Rev. B 83, 205403 (2011).
  10. Li, G. et al. Observation of van hove singularities in twisted graphene layers. Nat. Phys. 6, 109 –113 (2010).
  11. Wong, D. et al. Local spectroscopy of moiré-induced electronic structure in gate-tunable twisted bilayer graphene. Phys. Rev. B 92, 155409 (2015).
  12. Bistritzer, R. & MacDonald, A. H. Moiré bands in twisted double-layer graphene. Proc. Natl Acad. Sci. USA 108, 12233 –12237 (2011).
  13. Jung, J., DaSilva, A. M., MacDonald, A. H. & Adam, S. Origin of band gaps in graphene on hexagonal boron nitride. Nat. Commun. 6,1 –11 (2015).
  14. Carr, S., Fang, S., Zhu, Z. & Kaxiras, E. Exact continuum model for low-energy electronic states of twisted bilayer graphene. Phys. Rev. Res. 1, 013001 (2019).
  15. Jiang, Y. et al. Charge order and broken rotational symmetry in magic-angle twisted bilayer graphene. Nature 573,9 1 –95 (2019).
  16. Kerelsky, A. et al. Maximized electron interactions at the magic angle in twisted bilayer graphene. Nature 572,9 5 –100 (2019).
  17. Xie, Y. et al. Spectroscopic signatures of many-body correlations in magic- angle twisted bilayer graphene. Nature 572, 101 –105 (2019).
  18. Choi, Y. et al. Electronic correlations in twisted bilayer graphene near the magic angle. Nat. Phys. 15, 1174 –1180 (2019).
  19. Liu, X. et al. Spin-polarized correlated insulator and superconductor in twisted double bilayer graphene. Nature 583, 221 –225 https://www.nature.com/ articles/s41586-020-2458-7 (2020).
  20. Tang, Y. et al. Wse2/ws2 moir\ ’e superlattices: a new hubbard model simulator. Nature 579, 353 –358 https://www.nature.com/articles/s41586-020- 2085-3 (2020).
  21. Wang, L. et al. Magic continuum in twisted bilayer wse2. Nature Materials 19, 861 https://www.nature.com/articles/s41563-020-0708-6 (2020).
  22. Adak, P. C. et al. Tunable bandwidths and gaps in twisted double bilayer graphene system on the verge of correlations. Phys. Rev. B 101, 125428 https:// journals.aps.org/prb/abstract/10.1103/PhysRevB.101.125428(2020).
  23. He, M. et al. Symmetry breaking in twisted double bilayer graphene. Nat. Phys. 17,2 6 –30 (2021).
  24. Choi, Y. W. & Choi, H. J. Strong electron-phonon coupling, electron-hole asymmetry, and nonadiabaticity in magic-angle twisted bilayer graphene. Phys. Rev. B 98, 241412 (2018).
  25. Lian, B., Wang, Z. & Bernevig, B. A. Twisted bilayer graphene: a phonon- driven superconductor. Phys. Rev. Lett. 122, 257002 (2019).
  26. Sharma, G., Trushin, M., Sushkov, O. P., Vignale, G. & Adam, S. Superconductivity from collective excitations in magic-angle twisted bilayer graphene. Phys. Rev. Res. 2, 022040 (2020).
  27. Xu, C. & Balents, L. Topological superconductivity in twisted multilayer graphene. Phys. Rev. Lett. 121, 087001 (2018).
  28. Roy, B. & Juricic, V. Unconventional superconductivity in nearly flat bands in twisted bilayer graphene. Phys. Rev. B 99, 121407 https://journals.aps.org/prb/ abstract/10.1103/PhysRevB.99.121407 (2019). 2 9 . P o ,H .C . ,Z o u ,L . ,V i s h w a n a t h ,A .&S e n t h i l ,T .O r i g i no fm o t ti n s u l a t i n gb e h a v i o r and superconductivity in twisted bilayer graphene.Phys. Rev. X 8, 031089 (2018).
  29. Koshino, M. et al. Maximally localized wannier orbitals and the extended hubbard model for twisted bilayer graphene. Phys. Rev. X 8, 031087 (2018).
  30. Kang, J. & Vafek, O. Symmetry, maximally localized wannier states, and a low-energy model for twisted bilayer graphene narrow bands. Phys. Rev. X 8, 031088 (2018).
  31. Guo, H., Zhu, X., Feng, S. & Scalettar, R. T. Pairing symmetry of interacting fermions on a twisted bilayer graphene superlattice. Phys. Rev. B 97, 235453 (2018).
  32. Isobe, H., Yuan, N. F. Q. & Fu, L. Unconventional superconductivity and density waves in twisted bilayer graphene. Phys. Rev. X 8, 041041 (2018).
  33. Xie, M. & MacDonald, A. H. On the nature of the correlated insulator states in twisted bilayer graphene. Phys. Rev. Lett 124, 097601 https://journals.aps.org/ prl/abstract/10.1103/PhysRevLett.124.097601 (2020).
  34. Laksono, E. et al. Singlet superconductivity enhanced by charge order in nested twisted bilayer graphene fermi surfaces.Solid State Commun. 282,3 8–44 (2018).
  35. Padhi, B., Setty, C. & Phillips, P. W. Doped twisted bilayer graphene near magic angles: proximity to wigner crystallization, not mott insulation. Nano Lett. 18, 6175 –6180 (2018).
  36. Slater, J. Magnetic effects and the hartree-fock equation. Phys. Rev. 82, 538 (1951).
  37. Rademaker, L. & Mellado, P. Charge-transfer insulation in twisted bilayer graphene. Phys. Rev. B 98, 235158 (2018).
  38. Guinea, F. & Walet, N. R. Electrostatic effects, band distortions, and superconductivity in twisted graphene bilayers. Proc. Natl Acad. Sci. USA 115, 13174–13179 (2018).
  39. Yudhistira, I. et al. Gauge phonon dominated resistivity in twisted bilayer graphene near magic angle. Phys. Rev. B 99, 140302(R) (2019).
  40. Cao, Y. et al. Strange metal in magic-angle graphene with near planckian dissipation. Phys. Rev. Lett. 124, 076801 (2020).
  41. Polshyn, H. et al. Large linear-in-temperature resistivity in twisted bilayer graphene. Nat. Phys. 15, 1011 –1016 (2019).
  42. Bruin, J. A. N., Sakai, H., Perry, R. S. & Mackenzie, A. P. Similarity of scattering rates in metals showing t-linear resistivity. Science 339, 804 –807 (2013).
  43. Hartnoll, S., Lucas, A. & Sachdev, S. Holographic Quantum Matter (MIT Press, 2018).
  44. Tan, C. et al. Realization of a universal hydrodynamic semiconductor in ultra- clean dual-gated bilayer graphene. Preprint at arXiv:1908.10921 (2019).
  45. Patel, A. A. & Sachdev, S. Theory of a planckian metal. Phys. Rev. Lett. 123, 066601 (2019).
  46. Hwang, E. H. & Das Sarma, S. Acoustic phonon scattering limited carrier mobility in two-dimensional extrinsic graphene.P h y s .R e v .B77, 115449 (2008).
  47. de Gail, R., Goerbig, M. O., Guinea, F., Montambaux, G. & Castro Neto, A. H. Topologically protected zero modes in twisted bilayer graphene. Phys. Rev. B 84, 045436 (2011).
  48. Suzuura, H. & Ando, T. Phonons and electron-phonon scattering in carbon nanotubes. Phys. Rev. B 65, 235412 (2002).
  49. Sohier, T. et al. Phonon-limited resistivity of graphene by first-principles calculations: electron-phonon interactions, strain-induced gauge field, and boltzmann equation. Phys. Rev. B 90, 125414 (2014).
  50. Okuyama, Y. & Tokuda, N. Electron-phonon interactions in modulation- doped al xga1-x as/gaas heterojunctions. Phys. Rev. B 40, 9744 –9750 (1989).
  51. Kawamura, T. & Das Sarma, S. Phonon-scattering-limited electron mobilities in al xga1-x as/gaas heterojunctions. Phys. Rev. B 45, 3612 –3627 (1992).
  52. Efetov, D. K. & Kim, P. Controlling electron-phonon interactions in graphene at ultrahigh carrier densities. Phys. Rev. Lett. 105, 256805 (2010).
  53. Wu, F., Hwang, E. & Das Sarma, S. Phonon-induced giant linear-in- t resistivity in magic angle twisted bilayer graphene: ordinary strangeness and exotic superconductivity. Phys. Rev. B 99, 165112 (2019). ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 10 NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications
  54. Xie, H.-Y. & Foster, M. S. Transport coef ficients of graphene: interplay of impurity scattering, coulomb interaction, and optical phonons. Phys. Rev. B 93, 195103 (2016).
  55. Carr, S., Fang, S., Zhu, Z. & Kaxiras, E. Exact continuum model for low-energy electronic states of twisted bilayer graphene. Phys. Rev. Res. 1, 013001 (2019).
  56. Leconte, N., Javvaji, S., An, J. & Jung, J. Relaxation effects in twisted bilayer graphene: a multi-scale approach. Preprint at arXiv:1910.12805 (2019).
  57. Alezzi, M. & Adam, S. in Preperation (2021).
  58. Adam, S., Hwang, E. H., Galitski, V. M. & Das Sarma, S. A self-consistent theory for graphene transport. Proc. Natl Acad. Sci. USA 104, 18392 –18397 (2007).
  59. Rossi, E., Adam, S. & Das Sarma, S. Effective medium theory for disordered two-dimensional graphene. Phys. Rev. B 79, 245423 (2009).
  60. Jang, C. et al. Tuning the effective fine structure constant in graphene: opposing effects of dielectric screening on short- and long-range potential scattering. Phys. Rev. Lett. 101, 146805 (2008).
  61. Park, J. M., Cao, Y., Watanabe, K., Taniguchi, T. & Jarillo-Herrero, P. Flavour hund’s coupling, correlated chern gaps, and diffusivity in moiré flat bands. Nature 592,4 3 https://www.nature.com/articles/s41586-021-03366-w (2021).
  62. Kaasbjerg, K., Thygesen, K. S. & Jacobsen, K. W. Unraveling the acoustic electron-phonon interaction in graphene. Phys. Rev. B 85, 165440 (2012).
  63. Guinea, F. & Walet, N. R. Continuum models for twisted bilayer graphene: effect of lattice deformation and hopping parameters. Phys. Rev. B 99, 205134 (2019). Acknowledgements We acknowledge the Singapore Ministry of Education AcRF Tier 2 (Grant No. MOE2017- T2-2-140), the Singapore National Research Foundation Investigator Award (Grant No. NRF-NRFI06-2020-0003), and the use of the dedicated research computing resources at CA2DM. M.S.F. and S.A. acknowledge the support of the ARC through grants CE170100039 and DP200101345. It is a pleasure to thank C. Dean, E. Laksono, P. Jarillo- Herrero H. Mahalingam, N. Raghuvanshi, and M. Yankowitz for fruitful discussions. Author contributions G.S. and I.Y. contributed equally. S.A. supervised the project. G.S. calculated the dyna- mical screening of the electron –phonon coupling with assistance from D.Y.H.H. and G.V. I.Y. and N.C. calculated the electron –phonon resistivity beyond the Dirac model. I.Y. developed the Planckian model with assistance from M.S.F. G.S. and I.Y. calculated the electron–phonon resistivity with input from M.M.A.E. G.S., I.Y., N.C., and S.A. wrote the manuscript with input from all authors. Competing interests The authors declare no competing interests. Additional information Supplementary information The online version contains supplementary material available at https://doi.org/10.1038/s41467-021-25864-1. Correspondence and requests for materials should be addressed to Shaf fique Adam. Peer review information Nature Communications thanks Angelika Knothe and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Peer reviewer reports are available. Reprints and permission information is available at http://www.nature.com/reprints Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional af filiations. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article ’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/ licenses/by/4.0/. © The Author(s) 2021 NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-021-25864-1 ARTICLE NATURE COMMUNICATIONS | (2021) 12:5737 | https://doi.org/10.1038/s41467-021-25864-1 | www.nature.com/naturecommunications 11