Momentum shift current at terahertz frequencies in twisted bilayer graphene
MMomentum shift current at terahertz frequenciesin twisted bilayer graphene
Daniel Kaplan, Tobias Holder, and Binghai Yan ∗ Department of Condensed Matter Physics, Weizmann Institute of Science, Rehovot 7610001, Israel
The detection of terahertz (THz) radiation promises intriguing applications in biology, telecommu-nication, and astronomy but remains a challenging task so far. For example, semiconductor infrareddetectors (e.g., HgCdTe) utilize photo-excited electrons across the bandgap and hardly reach thefar-infrared terahertz regime because they are vulnerable to thermally excited carriers. In this work,we propose the THz sensing by the bulk photovoltaic effect (BPVE) in the twisted bilayer graphene(TBG). The BPVE converts light into a coherent DC at zero bias or an open-circuit voltage. As aquantum response from the wave function’s geometry (different from the p - n junction), BPVE ismore robust against temperature excitation and disorders. We predict that the TBG (bandgap ofseveral meV) exhibits a sizeable BPVE response in a range of 0.2 – 1 THz. Beyond the ordinaryshift current scenario, BPVE in TBG comes from a momentum-space shift of flat bands. Our workprovides a pathway to design twisted photonics for resonant terahertz detection. INTRODUCTION
The bulk photovaltaic effect (BPVE) refers to the gen-eration of a dc current from a homogeneous solid uponirradiation with light. The so-called shift current [1–4]is one of the most important mechanisms for the BPVEwhich has been appreciated for a long time for its po-tential in photovoltaic and photonic applications. Sofar, the BVPE has been studied predominantly for near-infrared and optical frequencies, but it may prove usefulto overcome the so-called ’terahertz gap’ [5] in the fieldof terahertz photonics [6]. This gap is the perceivedlack of robust, tunable and broad-frequency materials forterahertz detection. For instance, usual semiconductorinfrared detectors [7] harvest the photon-excited conduc-tion electrons, which are competing with the thermalexcitation across the energy gap, and thus cannot work inthe THz range. Recent developments attempt to detectplasmons [8–10] stimulated by the THz radiation.We point out that the shift current is quite robustagainst thermally excited carriers, because it is a coherentbulk effect driven by the quantum geometry [11–18] andworks even in metallic systems such as Weyl semimet-als [11, 19–26]. Most importantly, the created currentdoes not require any DC bias or an interface to becomeaccessible to measurement. Therefore, the thermal occu-pation of the conduction band at finite temperatures doesnot present a source of background current, meaning thatthe shift current does not compete with thermal currents,save for thermal noise. While the shift current is indeeddiminished by occupying the conduction band, this effectis only moderate in the terahertz regime. For example,even at room temperature the thermal occupation acrossa 4meV gap will lead to a shift current signal with a mag-nitude of still 5% of its zero temperature value. Thus,the shift current can be useful for higher temperatureterahertz sensing based on noncentrosymmetric materialsthat exhibit resonant BPVE response in the meV window.The shift current may substantially improve detectioncapabilities for terahertz radiation, increase the photo- voltaic efficiency for energy harvesting, and have furtherapplications in medical imaging, single-photon circuitsand novel electronic devices [27, 28].Twisted bilayer graphene (TBG) exhibits flat bandsand narrow band gaps of several meV, and thus pro-vides an ideal platform to investigate the THz BPVE.It has attracted an immense amount of attention fol-lowing the recent discovery of correlated insulating andsuperconducting phases at small twist angle (magic angle) θ ∼ ◦ [29, 30]. While the root cause for this exciting be-havior is widely accepted to be tied to a highly quenchedband structure [31–33], opinions differ about the mecha-nisms associated with the various correlated phases whichhave since been documented in the system [34–46]. Animportant stepping stone towards the understanding ofthese phases is the characterization of the quantum geo-metric and dispersive features of the flat bands. Becauseof inversion symmetry breaking, TBG and similar twistedbilayers exhibit intriguing nonlinear phenomena [47–53]such as the bulk photovaltaic effect (BPVE) [1, 54] andnonlinear anomalous Hall effect [55, 56]. These nonlinearprobes are ideally suited for the investigation of quasipar-ticle properties in flat bands as they are not suppressedby a vanishing Fermi velocity [18].With regards to semiconductors, the shift current hadoriginally been interpreted as a real-space shift of wavefunction centers upon excitation. It has gained renewedinterest in recent studies of Weyl semimetals [57, 58], be-cause the shift is believed to get enhanced by the Berrycurvature near Weyl points [11, 19–21, 23–26]. Very re-cently, we proposed to instead view the shift current as aresult of the anomalous quasiparticle acceleration, whichis determined by the quantum geometric properties [18].In the latter interpretation, the small effective mass in aWeyl semimetal leads to strong acceleration in the field oflight and hence generates large photocurrent, consistentwith the shift current interpretation [11]. Flat bands, asfor example in TBG, represent a kind of opposite limit(large quasiparticle mass) compared to Weyl bands (nearlymassless quasiparticles). Using the language of a quasipar- a r X i v : . [ c ond - m a t . m e s - h a ll ] J a n Figure 1. A schematic comparison of the wavefunction inreal-space and momentum-space of electronic states at theFermi surface for three types of gapless materials. For aWeyl semimetal, states are well localized at few points in themomentum space but quite extended in the real space. Incontrast, for a flat band all momenta are equally occupied,leading to sharply localized peaks of the wavefunction in realspace. The ordinary metal with a large Fermi surface repre-sents the intermediate region between the Weyl semimetal andthe flat-band system. If the shift current is viewed as a resultof the anomalous acceleration, the Weyl semimetal exhibitsa large real-space shift while the flat band displays a majormomentum-space shift. ticle shift, it is therefore less obvious to deduce the shiftcurrent for a flat-band system. In contrast, in the semi-classical acceleration picture, a state can be displaced (i.e.accelerated) in either real-space and momentum-space,as illustrated in Fig. 1. Although the real-space shift isindeed small in real-space due to the large quasiparticlemass, for flat bands the acceleration (i.e. the displacementin momentum space) can still be very large. Therefore,we are motivated to investigate the shift current of TBGfrom the acceleration point of view.A previous study has reported a non-linear conductiv-ity above 30meV due to the orbital magnetization at 3/4filling [47]. In this work, we study shift current generationbelow 10meV in magic angle TBG across the single parti-cle gaps. In a numerical analysis, find a large photocurrentresponse for light in the low terahertz regime. This isconsistent with the acceleration picture outlined above.To further support such a connection, we formally expressthe response in terms of the real-space shift current andmomentum-space shift current. Due to the flat-band dis-persion, the dominant contribution is constituted by themomentum-space shift, which we explicitly show to bedissimilar to an ordinary dispersive material (e.g., MoS ). Distinct from the nonlinear anomalous Hall effect, theshift current of TBG originates from the properties ofthe quantum geometry of the band structure that is un-related to the Berry curvature dipole. The magnitudeof the obtained current, tunability of the Hall response,the broadness of the resonances, and the terahertz regimein which they are all observed hold promise of the uti-lization of TBG devices in terahertz technologies, suchphotodetection and photovoltaic DC current generation. RESULTSA. Shift current theory
In the clean limit, the Bloch states produce an intrinsicdc-current response to linearly polarized light [2]. Thisnonlinear conductivity σ aa ; c ( s ) is commonly formulated [2, 3]as σ aa ; c ( s ) (0; ω, − ω ) = πe (cid:126) (cid:90) k (cid:88) mn f mn S cmn | r amn | δ ( ω ± ε mn ) , (1)where m, n are the band indices, ε mn = ε m − ε n theband energy difference, f mn = f ( ε m ) − f ( ε n ) the Fermidistribution function difference, and r amn ≡ (cid:104) m | r a | n (cid:105) thedipole transition matrix element. Here, a, c representsthe light field and current directions, respectively. S mn is the so-called shift vector, S cmn = ( r cmm − r cnn ) + ∂ k c arg r cmn , (2)where ( r cmm − r cnn ) is the shift of the wave function centers,which is gauge dependent.The second term in Eq. 2 is the phase derivative of r cmn , which ensures gauge-invariance for the shift vector. S mn has then been interpreted as the real-space shift ofthe mass center of a quasiparticle upon excitation fromband m to n [2]. As we will show now, the phase term in(2) gives rise to a shift in momentum space. Therefore,the shift current is actually the result of a shift of thequasiparticle excitation in both real space and momentumspace, which is naturally captured using the language ofthe anomalous quasiparticle acceleration [18].To this end, we express the nonlinear conductivity inthe form of geometric properties of the band structure asfollows, σ aa ; c ( s ) (0; ω, − ω ) = − πe (cid:126) (cid:90) k R aac shift + K aac shift (3) - - - + [ a . u . ] (a) (b)(c) Figure 2. Band structure, Berry curvature and nonlinear conductivity in the mini-Brillouin zone of TBG. (a) Band structureof TBG with inversion-breaking at twist angle θ = 1 . ◦ . Solid lines show the unstrained case, the dashed lines are at strain (cid:15) = 0 . K . The colorbar is in logarithmicscale. The Berry curvature in each valley is C z symmetric, which leads to a vanishing Berry curvature dipole in each valleyseparately. (c) Modulus of the nonlinear conductivity as given by σ yy ; x ( s ) [Eq. (8)], in logarithmic scale. At frequency ω = 1 . R aac shift = (cid:88) mn f mn ( r cmm − r cnn ) | r amn | δ ( ω ± ε mn ) (4) K aac shift = (cid:88) mn f mn (cid:104) ir amn λ acnm − ir cmn λ aanm (cid:105) δ ( ω ± ε mn ) − ∂ k a Ω ac | ε mn = ± ω (5) ∂ k a Ω ac | ε mn = ± ω = (cid:88) mn f mn (cid:104) i r amn Ω acnm − ir cmn λ aanm + ir amn λ acnm (cid:105) δ ( ω ± ε mn ) (6)where we use the symmetrized derivative λ abnm = (cid:0) ∂ k a r bnm + ∂ k b r anm (cid:1) . The two pieces R aac shift and K aac shift correspond to the contributions from the first termand the second term, respectively, of the shift vectorin Eq. (2). The second term in K aac shift is the Berrycurvature dipole, and the first term is related to theother geometric quantity, the quantum metric, whichcharacterizes the quantum distance between two states.Thus, K aac shift originates directly from the quantum ge-ometry of the wave function in k -space, and it repre-sents the anomalous acceleration in momentum space [18].The quantity λ abnm encodes the skewness of the accel-eration: The Berry curvature dipole involves deriva-tives of the type ∂ c ( r amn r bnm − r bmn r anm ), which doesnot cover the entire motion in momentum space. Theremaining terms involving λ ab in Eq. (5) can be con-nected to skew-symmetric derivatives of the structure( ∂ c r amn ) r bnm − r amn ( ∂ c r bnm ) + ( ∂ c r bmn ) r anm − r bmn ( ∂ c r anm ).Note that in the special case of a two-band Dirac dis- persion, it can be shown that the general expressions forthe shift current present here can be connected by a pull-back mapping to the Christoffel symbols of the geometricconnection on the generalized Bloch sphere [17]. Whilethis reinforces the association of the shift current with asemiclassical acceleration, for the general multi-band casediscussed here, the more convenient expression is the onepresented in Eq. (3) in terms of λ ab . B. Shift current of TBG
We model TBG using a modified form of the Bistrizer-Macdonald continuum model [33, 48]. We attach twomonolayer graphene sheets, with first rotating them withrespect to one another with an angle θ = 1 . o , and thenintroducing inter-layer coupling. By construction, thismodel is endowed with C z symmetry, and since bothmonolayers are inversion symmetric, the model as a whole - - (a) (b) Figure 3. (a) σ yy ; x of TBG with ε = 0%, with C symmetry intact. The conductivity is plotted for 3 values of µ = − , , R yy ; x shift and K yy ; x shift pieces, respectively. Solid lines depict TBG, while thedashed lines show MoS , for comparison. The shift current is essentially zero for x <
1, has a large resonance around x ≈ x ≈ µ = − . ω = 1 . σ yy ; x = +198 µ AV − nm, and for µ = 6 . − µ AV − nm. has inversion symmetry. Inversion symmetry breakingis introduced by a staggered potential ∆, as it typicallyarises from encapsulation of the bilayer in hBN. Onthe other hand, C z is not necessarily broken because itrequires some uniaxial strain of magnitude ε in the bottomlayer. A detailed description of the model is given in thesupplementary material. The resulting band structure inthe mini-Brillouin zone of valley K is shown in Fig. 2a,both for the unstrained and strained TBG. The flat bandsnear half filling ( µ = 0 meV) are gapped at Γ and K dueto inversion symmetry breaking. Fig. 2b shows the Berrycurvature in the mini-Brillouin zone in unstrained TBG,which is C -symmetric. We note that the Berry chargesnear the K and K (cid:48) -points in the original dispersion haveopposite signs and cancel.While the shift vector S in Eq. (2) itself is not expectedto be further decomposable into gauge invariant pieces,the expression that enters the shift current can indeed bedecomposed. Writing σ ( s ) = σ ( s + σ ( s , with σ aa ; c ( s = πe (cid:126) (cid:90) k ∂ k a Ω ac | ε mn = ω (7) σ aa ; c ( s = − πe (cid:126) (cid:90) k ( R aac shift + K aac shift ) − σ aa ; c ( s . (8)Here, Eq. (7) contains the contributions from the Berrycurvature dipole ∂ a Ω ac , which vanishes in the presenceof C z symmetry. This symmetry renders only two com-ponents of σ ab ; c independent. For simplicity we focusin the main text on σ yy ; x , which encodes the transversenonlinear conductivity response. The other independentcomponent is σ xx ; y (Supplementary Figs. 2 and 3); thesymmetry analysis can be found in the SI. In Fig. 3awe present the total shift current in the unstrained case,for three values of the chemical potential which reside within the three gaps that are opened by the staggeredpotential. Within our numerical accuracy, the currentis indeed found to have no contributions from the Berrycurvature dipole, i.e. σ yy ; x ( s = σ yy ; x ( s ) . Irrespective of this,for frequencies which cross the single-particle gap theshift current reaches giant values nearly 200 µ AV − nm,far exceeding predicted values for other materials [13, 59].As we pointed out, from the quasiparticle shift it isnot immediately obvious where such a giant non-linearresponse could originate from. For this reason, we exam-ine the different contributions in Eq. (8). Fig. 2c showsthe momentum-space structure of the integrand in Eq.(8) at µ = 0 , ω = 1 . K s , both of which are fairlybroad in momentum space. Fig. 3b illustrates the relativecontribution of the pieces K aa ; c shift and R aa ; c shift in Eq. (8), asa function of the distance from the gap at 0. Over alarge range of frequencies, it holds that K aa ; c shift (cid:29) R aa ; c shift ,meaning that the shift current is produced mostly by thephase of the Berry connection and not the shift of wave-function centers. For comparison, we performed the samebreakdown for the much more dispersive, gapped materialMoS , which has very similar symmetry properties [Fig.3b]. In this latter case, K aa ; c shift (cid:28) R aa ; c shift near the bandedge, which supports the notion that the shift currentis resulting from a displacement in both real space andmomentum space, with the latter being greatly enhancedfor a flatband dispersion. As we mentioned in the be-ginning, these statements might seem questionable uponregauging, a possible shortcoming on which we commentin the discussion. - - - (a) (b) (c) Figure 4. Shift current response of strained TBG, σ yy ; x , for 3 chemical potentials within the gaps created by the staggeredpotential, and uniaxial strain ε = 0 . σ yy ; x at T = 0K (solid lines) and for T = 300K (dashed lines, shown x10 for clarity). For chemical potential µ = 0 the shiftpeaks at σ yy ; x = 315 µ AnmV − ( T = 0K) and 15 µ AnmV − ( T = 300K), respectively. Inset: relative size of R yy ; x shift and K yy ; x shift contributing to the total conductivity as a function of the distance from the flat-band band gap, xE g (here µ = 0), cf. Fig. 3b. C. Effect of strain
In the presence of uniaxial strain, the symmetry groupof TBG is reduced to C . All components of the conduc-tivity are now independent, but for clarity we continue toexamine only the transverse contribution σ yy ; x . Since theBerry curvature dipole contribution as given by Eq (7) isno longer zero, we present both its contribution (Fig. 4a)and the remaining terms in Eq. 8 for σ aa ; c ( s (Fig. 4b). Thetotal shift current is depicted in Fig. 4c, with the relativesizes of R aa ; c shift and K aa ; c shift shown in the inset.At moderate uniaxial strain of size (cid:15) = 0 . σ ( s due to the Berry curvature dipole becomescomparable in size to σ ( s , but it has consistently theopposite sign. Thus, while the total conductivity unsupris-ingly increases due to the reduced symmetry in the system,it has less accentuated resonances for the transitions at theband edges, with values up to 300 µ AV − nm across theflatband gap. More unexpectedly, the shift current doesnot seem to profit from imbalanced Berry charges, as theircontribution either subtracts from the remaining current,or is almost negligible for transitions across the flatbands.Our results establish that a large anomalous accelerationin the quasiparticle motion can arise even if ∂ k a Ω ac = 0,i.e. its existence does not rely on the presence of Berrycharges in the system. This is an important distinctionbetween the bulk photovoltaic effect and the anomalousHall effect in TBG. We further observe that in all casesthe shift current at frequencies ω (cid:29) T = 300K (dashed lines). While the transitions across the gaps of the dispersive bands are strongly suppressed,the large density of states in the flat bands supports astrong signal for (cid:126) ω = 4meV, with an amplitude of still5% of its value at T = 0. DISCUSSION
The shift vector S has previously been connected to thereal-space shift of the center-of-mass coordinate betweentwo eigenstates upon excitation from the conduction band m into the valence band n . To give some intuition, weexpand the real-space representation of the periodic eigen-functions | u n k (cid:105) in terms of local Wannier orbitals | w n R (cid:105) with center coordinate R [60] (cid:104) r | u n k (cid:105) = (cid:88) R e − i k ( R − r ) (cid:104) r | w n R (cid:105) . (9)Then, in momentum space the Berry connection is givenby r amn = (cid:88) RR (cid:48) e i k ( R − R (cid:48) ) (cid:90) cell dV (cid:104) w n R (cid:48) | r a | w n R (cid:105) . (10)Evaluating the shift vector S cmn = r cmm − r cnn + ∂ k c arg r cmn based on this representation yields r amm − r ann = R amm − R ann for the direct difference. This is supplemented bythe phase derivative ∂ k a arg r mn , whose integral over theBrillouin zone is a multiple of 2 π . If only two Wannierorbitals have a significant overlap, the modulus | S | isclearly bounded by | R amm − R ann | < a , with lattice constant a . If several orbitals overlap, the phase factors in thesum Eq. (10) become important, with slope of growthin momentum space being at most a . Then, the phasederivative is expected to contribute similarly at O ( a ) tothe shift vector. This already indicates that interpretationof the shift current as a result of the wavefunction shift isnarrow to some extent. For ease of illustration, imaginea set of Landau levels in symmetric gauge. Their center-of-mass coordinate R can be moved around freely inexchange for acquiring an additional phase factor. Indeed,for generic flat bands it is to be expected that there is agauge choice which makes the shift S only depend on thephase, because the center of the Wannier functions can berepositioned with an appropriate gauge transformation.This is inconsistent with the interpretation of the shiftcurrent as the real space shift of the wavefunction centerupon absorption of a photon, as the wavefunction onlysuffers a phase shift.If the shift current is instead viewed as the anoma-lous acceleration that a quasiparticle undergoes due tothe interaction with the electric field at second order,the photogalvanic response follows as a straightforwardgeneralization of the linear response formalism involvingthe anomalous velocity [18, 61], thus removing the directinference of a current from a real-space displacement. In-stead, both R shift and K shift appear as the result of thesame acceleration that changes both the position and thewavevector of the quasiparticle. As shown in the lastsection, this is consistent with our numerical findings for R shift and K shift using Bloch wavefunctions for a mostlysmooth gauge choice within each band [cf. Fig. 3b].While one might object against inspecting gauge-dependent quantities, we emphasize that R shift and K shift can still contain valuable information about the quasipar-ticle dynamics in the sense that while they are not unique,this does not at all imply that they are arbitrary. For this,recall that all types of topological bands (including flatbands) have no uniquely defined center-of-mass coordinatewithin the unit cell, because the topological nature of thebands prevents such an assignment. However, from thisit does not follow that the momentum space integral of R shift can take arbitrary values, because it contains muchmore specific information about the relative positionaldifference between two bands, summed for all momenta.Indeed, it was already pointed out a long time ago [3, 62]that for an arbitrary band structure one cannot generallyexpect to find a gauge such that r amm − r ann consistentlyvanishes for all momenta and all bands.Summarizing these observations, we therefore pose thata useful indicator for band flatness is that the integratedpositional difference r amm − r ann between Bloch wavefunc-tions can be made substantially smaller than the inte-grated phase contribution ∂ k a arg r mn . We further conjec-ture that for highly dispersive bands a similar statementshould hold about the smallness of the integrated phasecontribution. A paradigmatic example in the latter case is a two-band semimetal with one band crossing. There,a mostly smooth gauge is at the same time periodic (i.e.without phase jump at the Brillouin zone boundaries),thus completely eliminating the phase contribution fromthe integrated shift vector. This is the expected resultfor a quasiparticle with vanishing effective mass which ischanges position in an applied electric field.Before closing, we remark that the difficulties in separat-ing real-space and momentum-space effects of the acceler-ation into gauge invariant pieces are intrinsic to the morecomplicated semiclassical motion arising at second orderin the applied field. In particular, the analogous splittingof the quasiparticle velocity into the regular (dispersive)and anomalous velocity has the important distinctionthat these two components of the velocity are orthogonalto each other, making them linearly independent. Sucha decomposition is not straightforward for the acceler-ation, because it describes the changes to both regularand anomalous velocity components in both normal andperpendicular direction, thus mixing them. However, webelieve that the conclusions outlined above can be mademore rigorous by deriving a lower bound for R shift and K shift , which will be the subject of a future work.In summary, we report a giant photogalvanic currentfor TBG which is irradiated by linear polarized light inthe terahertz range. The resonance profile we observein both the strained and unstrained cases suggest thatTBG is a promising candidate for THz detection andcircuits, even at room temperature. As the root causebehind the large response we identified the anomalousacceleration due to the skew symmetric properties of thequantum geometry of the band structure as encoded in λ abmn , which always appear in the shift current, but aregreatly amplified in TBG due to the flat-band dispersion.The latter also turned out to be particularly important forretaining a large shift conductivity at higher temperatures.We expect the transverse dc-current reported here tobe accessible using current samples and measurementtechniques [63]. The line of reasoning developed in thiswork can potentially shed light on the quantum geometryof the band structure in similarly twisted van-der-Waalsmaterials with nearly flat bands, for example MoTe orWSe [64]. Acknowledgements
We thank J. S. Hofmann,R. M. Ribeiro, and R. Queiroz for useful discussions. B.Y.acknowledges the financial support by the Willner FamilyLeadership Institute for the Weizmann Institute of Sci-ence, the Benoziyo Endowment Fund for the Advancementof Science, Ruth and Herman Albert Scholars Programfor New Scientists, the European Research Council (ERCConsolidator Grant No. 815869, “NonlinearTopo”). ∗ [email protected][1] V. I. Belinicher and B. I. Sturman, “The photogalvaniceffect in media lacking a center of symmetry,” Sov. Phys.Usp. , 199–223 (1980). [2] R. von Baltz and W. Kraut, “Theory of the bulk photo-voltaic effect in pure crystals,” Phys. Rev. B , 5590–5596 (1981).[3] J. E. Sipe and A. I. Shkrebtii, “Second-order optical response in semiconductors,” Phys. Rev. B , 5337–5352 (2000).[4] Steve M. Young and Andrew M. Rappe, “First PrinciplesCalculation of the Shift Current Photovoltaic Effect inFerroelectrics,” Phys. Rev. Lett. , 116601 (2012),arXiv:1202.3168 [cond-mat.mtrl-sci].[5] Wojciech Knap, Mikhail Dyakonov, Dominique Coquil-lat, Frederic Teppe, Nina Dyakonova, Jerzy (cid:32)Lusakowski,Krzysztof Karpierz, MacIej Sakowicz, Gintaras Valusis,Dalius Seliuta, Irmantas Kasalynas, Abdelouahad El Fa-timy, Y. M. Meziani, and Taiichi Otsuji, “Field effecttransistors for terahertz detection: Physics and first imag-ing applications,” in Journal of Infrared, Millimeter, andTerahertz Waves , Vol. 30 (Springer, 2009) pp. 1319–1337.[6] Masayoshi Tonouchi, “Cutting-edge terahertz technology,”Nat. Photonics , 97–105 (2007).[7] Michael A Kinch, “Fundamental physics of infrared detec-tor materials,” Journal of Electronic Materials , 809–817(2000).[8] Kamalodin Arik, Sajjad AbdollahRamezani, and AminKhavasi, “Polarization Insensitive and Broadband Tera-hertz Absorber Using Graphene Disks,” Plasmonics ,393–398 (2017).[9] V. M. Muravev and I. V. Kukushkin, “Plasmonic de-tector/spectrometer of subterahertz radiation based ontwo-dimensional electron system with embedded defect,”Appl. Phys. Lett. , 082102 (2012).[10] Denis A. Bandurin, Dmitry Svintsov, Igor Gayduchenko,Shuigang G. Xu, Alessandro Principi, Maxim Moskotin,Ivan Tretyakov, Denis Yagodkin, Sergey Zhukov, TakashiTaniguchi, Kenji Watanabe, Irina V. Grigorieva, MarcoPolini, Gregory N. Goltsman, Andre K. Geim, andGeorgy Fedorov, “Resonant terahertz detection usinggraphene plasmons,” Nat. Commun. , 5392 (2018),arXiv:1807.04703 [cond-mat.mes-hall].[11] T. Morimoto and N. Nagaosa, “Topological nature ofnonlinear optical effects in solids,” Sci. Adv. , e1501524–e1501524 (2016), arXiv:1510.08112 [cond-mat.mes-hall].[12] Liang Z. Tan, Fan Zheng, Steve M. Young, FenggongWang, Shi Liu, and Andrew M. Rappe, “Shift currentbulk photovoltaic effect in polar materials—hybrid andoxide perovskites and beyond,” npj Comput. Mater. ,16026 (2016).[13] A. M. Cook, B. M. Fregoso, F. de Juan, S. Coh, and J. E.Moore, “Design principles for shift current photovoltaics,”Nat. Commun. , 14176 (2017).[14] G. B. Ventura, D. J. Passos, J. M. B. Lopes dos Santos,J. M. Viana Parente Lopes, and N. M. R. Peres, “Gaugecovariances and nonlinear optical responses,” Phys. Rev.B , 035431 (2017), arXiv:1703.07796 [cond-mat.other].[15] D. J. Passos, G. B. Ventura, J. M. Viana ParenteLopes, J. M. B. Lopes dos Santos, and N. M. R. Peres,“Nonlinear optical responses of crystalline systems: Resultsfrom a velocity gauge analysis,” Phys. Rev. B , 235446(2018), arXiv:1712.04924 [cond-mat.other].[16] Daniel E. Parker, Takahiro Morimoto, Joseph Orenstein,and Joel E. Moore, “Diagrammatic approach to nonlinearoptical response with application to Weyl semimetals,”Phys. Rev. B , 045121 (2019), arXiv:1807.09285 [cond-mat.str-el].[17] Junyeong Ahn, Guang-Yu Guo, and Naoto Nagaosa,“Low-frequency divergence and quantum geometry ofthe bulk photovoltaic effect in topological semimetals,”arXiv , arXiv:2006.06709 (2020), arXiv:2006.06709 [cond- mat.mes-hall].[18] Tobias Holder, Daniel Kaplan, and Binghai Yan, “Conse-quences of time-reversal-symmetry breaking in the light-matter interaction: Berry curvature, quantum metric, anddiabatic motion,” Phys. Rev. Research , 033100 (2020),arXiv:1911.05667 [cond-mat.mes-hall].[19] Liang Wu, S. Patankar, T. Morimoto, N. L. Nair, E. The-walt, A. Little, J. G. Analytis, J. E. Moore, andJ. Orenstein, “Giant anisotropic nonlinear optical re-sponse in transition metal monopnictide Weyl semimetals,”Nat. Phys. , 350–355 (2017), arXiv:1609.04894 [cond-mat.mtrl-sci].[20] Qiong Ma, Su-Yang Xu, Ching-Kit Chan, Cheng-LongZhang, Guoqing Chang, Yuxuan Lin, Weiwei Xie, Tom´asPalacios, Hsin Lin, Shuang Jia, Patrick A. Lee, PabloJarillo-Herrero, and Nuh Gedik, “Direct optical detec-tion of Weyl fermion chirality in a topological semimetal,”Nat. Phys. , 842–847 (2017), arXiv:1705.00590 [cond-mat.mtrl-sci].[21] Gavin B. Osterhoudt, Laura K. Diebel, Mason J. Gray,Xu Yang, John Stanco, Xiangwei Huang, Bing Shen, Ni Ni,Philip J. W. Moll, Ying Ran, and Kenneth S. Burch,“Colossal mid-infrared bulk photovoltaic effect in a type-I Weyl semimetal,” Nat. Mater. , 471–475 (2019),arXiv:1712.04951.[22] P. Hosur and X. Qi, “Recent developments in transportphenomena in Weyl semimetals,” C. R. Phys. , 857–870 (2013), arXiv:1309.4464 [cond-mat.str-el].[23] Katsuhisa Taguchi, Tatsushi Imaeda, Masatoshi Sato,and Yukio Tanaka, “Photovoltaic chiral magnetic effectin Weyl semimetals,” Phys. Rev. B , 201202 (2016).[24] Ching-Kit Chan, Netanel H. Lindner, Gil Refael, andPatrick A. Lee, “Photocurrents in Weyl semimetals,”Phys. Rev. B , 041104 (2017), arXiv:1607.07839 [cond-mat.mes-hall].[25] Fernando de Juan, Adolfo G. Grushin, Takahiro Mori-moto, and Joel E. Moore, “Quantized circular photo-galvanic effect in Weyl semimetals,” Nat. Commun. ,15995 (2017), arXiv:1611.05887 [cond-mat.str-el].[26] Yang Zhang, Hiroaki Ishizuka, Jeroen van den Brink,Claudia Felser, Binghai Yan, and Naoto Nagaosa, “Pho-togalvanic effect in Weyl semimetals from first principles,”Phys. Rev. B , 241118 (2018), arXiv:1803.00562 [cond-mat.str-el].[27] R. A. Lewis, “A review of terahertz detectors,” J. Phys.D , 433001 (2019).[28] J. P. Guillet, B. Recur, L. Frederique, B. Bousquet,L. Canioni, I. Manek-H¨onninger, P. Desbarats, andP. Mounaix, “Review of terahertz tomography techniques,”Journal of Infrared, Millimeter, and Terahertz Waves ,382–411 (2014).[29] Yuan Cao, Valla Fatemi, Shiang Fang, Kenji Watanabe,Takashi Taniguchi, Efthimios Kaxiras, and Pablo Jarillo-Herrero, “Unconventional superconductivity in magic-angle graphene superlattices,” Nature , 43–50 (2018),arXiv:1803.02342 [cond-mat.mes-hall].[30] Yuan Cao, Valla Fatemi, Ahmet Demir, Shiang Fang,Spencer L. Tomarken, Jason Y. Luo, Javier D.Sanchez-Yamagishi, Kenji Watanabe, Takashi Taniguchi,Efthimios Kaxiras, Ray C. Ashoori, and Pablo Jarillo-Herrero, “Correlated insulator behaviour at half-filling inmagic-angle graphene superlattices,” Nature , 80–84(2018), arXiv:1802.00553 [cond-mat.mes-hall].[31] J. M. B. L. Santos, N. M. R. Peres, and A. H. Castro Neto, “Graphene Bilayer with a Twist: Elec-tronic Structure,” Phys. Rev. Lett. , 256802 (2007),arXiv:0704.2128 [cond-mat.mtrl-sci].[32] J. M. B. L. Santos, N. M. R. Peres, and A. H. CastroNeto, “Continuum model of the twisted graphene bilayer,”Phys. Rev. B , 155449 (2012), arXiv:1202.1088 [cond-mat.mtrl-sci].[33] Rafi Bistritzer and Allan H. MacDonald, “Moir´e bands intwisted double-layer graphene,” PNAS , 12233–12237(2011), arXiv:1009.4203 [cond-mat.mes-hall].[34] Aaron L. Sharpe, Eli J. Fox, Arthur W. Barnard,Joe Finney, Kenji Watanabe, Takashi Taniguchi,M. A. Kastner, and David Goldhaber-Gordon,“Emergent ferromagnetism near three-quarters filling intwisted bilayer graphene,” Science , 605–608 (2019),arXiv:1901.03520 [cond-mat.mes-hall].[35] M. Serlin, C. L. Tschirhart, H. Polshyn, Y. Zhang,J. Zhu, K. Watanabe, T. Taniguchi, L. Balents, andA. F. Young, “Intrinsic quantized anomalous Hall effectin a moir´e heterostructure,” Science , 900–903 (2020),arXiv:1907.00261 [cond-mat.str-el].[36] Xiaobo Lu, Petr Stepanov, Wei Yang, Ming Xie, Mo-hammed Ali Aamir, Ipsita Das, Carles Urgell, KenjiWatanabe, Takashi Taniguchi, Guangyu Zhang, AdrianBachtold, Allan H. MacDonald, and Dmitri K. Efetov,“Superconductors, orbital magnets and correlated statesin magic-angle bilayer graphene,” Nature , 653–657(2019), arXiv:1903.06513 [cond-mat.str-el].[37] Guorui Chen, Lili Jiang, Shuang Wu, Bosai Lyu,Hongyuan Li, Bheema Lingam Chittari, Kenji Watan-abe, Takashi Taniguchi, Zhiwen Shi, Jeil Jung, YuanboZhang, and Feng Wang, “Evidence of a gate-tunableMott insulator in a trilayer graphene moir´e superlattice,”Nat. Phys. , 237–241 (2019), arXiv:1803.01985 [cond-mat.mes-hall].[38] Yuhang Jiang, Xinyuan Lai, Kenji Watanabe, TakashiTaniguchi, Kristjan Haule, Jinhai Mao, and Eva Y. An-drei, “Charge order and broken rotational symmetry inmagic-angle twisted bilayer graphene,” Nature , 91–95(2019), arXiv:1904.10153 [cond-mat.mes-hall].[39] Youngjoon Choi, Jeannette Kemmer, Yang Peng, AlexThomson, Harpreet Arora, Robert Polski, Yiran Zhang,Hechen Ren, Jason Alicea, Gil Refael, Felix von Oppen,Kenji Watanabe, Takashi Taniguchi, and Stevan Nadj-Perge, “Electronic correlations in twisted bilayer graphenenear the magic angle,” Nat. Phys. , 1174–1180 (2019),arXiv:1901.02997 [cond-mat.mes-hall].[40] Alexander Kerelsky, Leo J. McGilly, Dante M.Kennes, Lede Xian, Matthew Yankowitz, Shaowen Chen,K. Watanabe, T. Taniguchi, James Hone, Cory Dean,Angel Rubio, and Abhay N. Pasupathy, “Maximizedelectron interactions at the magic angle in twisted bilayergraphene,” Nature , 95–100 (2019).[41] Yonglong Xie, Biao Lian, Berthold J¨ack, Xiaomeng Liu,Cheng-Li Chiu, Kenji Watanabe, Takashi Taniguchi,B. Andrei Bernevig, and Ali Yazdani, “Spectroscopicsignatures of many-body correlations in magic-angletwisted bilayer graphene,” Nature , 101–105 (2019),arXiv:1906.09274 [cond-mat.mes-hall].[42] Guorui Chen, Aaron L. Sharpe, Eli J. Fox, Ya-Hui Zhang,Shaoxin Wang, Lili Jiang, Bosai Lyu, Hongyuan Li, KenjiWatanabe, Takashi Taniguchi, Zhiwen Shi, T. Senthil,David Goldhaber-Gordon, Yuanbo Zhang, and FengWang, “Tunable correlated Chern insulator and ferro- magnetism in a moir´e superlattice,” Nature , 56–61(2020), arXiv:1905.06535 [cond-mat.mes-hall].[43] U. Zondiner, A. Rozen, D. Rodan-Legrain, Y. Cao,R. Queiroz, T. Taniguchi, K. Watanabe, Y. Oreg, F. vonOppen, Ady Stern, E. Berg, P. Jarillo-Herrero, andS. Ilani, “Cascade of phase transitions and Dirac revivalsin magic-angle graphene,” Nature , 203–208 (2020),arXiv:1912.06150 [cond-mat.mes-hall].[44] Dillon Wong, Kevin P. Nuckolls, Myungchul Oh, BiaoLian, Yonglong Xie, Sangjun Jeon, Kenji Watanabe,Takashi Taniguchi, B. Andrei Bernevig, and Ali Yaz-dani, “Cascade of electronic transitions in magic-angletwisted bilayer graphene,” Nature , 198–202 (2020),arXiv:1912.06145 [cond-mat.mes-hall].[45] Matthew Yankowitz, Shaowen Chen, Hryhoriy Polshyn,Yuxuan Zhang, K. Watanabe, T. Taniguchi, David Graf,Andrea F. Young, and Cory R. Dean, “Tuning super-conductivity in twisted bilayer graphene,” Science ,1059–1064 (2019), arXiv:1808.07865 [cond-mat.mes-hall].[46] S. L. Tomarken, Y. Cao, A. Demir, K. Watanabe,T. Taniguchi, P. Jarillo-Herrero, and R. C. Ashoori,“Electronic Compressibility of Magic-Angle Graphene Su-perlattices,” Phys. Rev. Lett. , 046601 (2019),arXiv:1903.10492 [cond-mat.mes-hall].[47] Jianpeng Liu and Xi Dai, “Anomalous Hall effect,magneto-optical properties, and nonlinear optical proper-ties of twisted graphene systems,” npj Comput. Mater. , 57 (2020), arXiv:1907.08932 [cond-mat.mes-hall].[48] Wen-Yu He, David Goldhaber-Gordon, and K. T. Law,“Giant orbital magnetoelectric effect and current-inducedmagnetization switching in twisted bilayer graphene,” Nat.Commun. , 1650 (2020).[49] Jin-Xin Hu, Cheng-Ping Zhang, Ying-Ming Xie, andK. T. Law, “Nonlinear Hall Effects in StrainedTwisted Bilayer WSe ,” arXiv , arXiv:2004.14140 (2020),arXiv:2004.14140 [cond-mat.mes-hall].[50] Meizhen Huang, Zefei Wu, Jinxin Hu, Xiangbin Cai,En Li, Liheng An, Xuemeng Feng, Ziqing Ye, Nian Lin,Kam Tuen Law, and Ning Wang, “Giant nonlinear Halleffect in twisted WSe ,” arXiv , arXiv:2006.05615 (2020),arXiv:2006.05615 [cond-mat.mes-hall].[51] Maximilian Otteneder, Stefan Hubmann, Xiaobo Lu,Dmitry A Kozlov, Leonid E Golub, Kenji Watanabe,Takashi Taniguchi, Dmitri K Efetov, and Sergey DGanichev, “Terahertz photogalvanics in twisted bilayergraphene close to the second magic angle,” Nano Lett. , 7152–7158 (2020).[52] Yizhou Liu, Tobias Holder, and Binghai Yan, “Chiral-ity induced Giant Unidirectional Magnetoresistance inthe Twisted Bilayers,” arXiv , arXiv:2010.08385 (2020),arXiv:2010.08385 [cond-mat.mes-hall].[53] Cheng-Ping Zhang, Jiewen Xiao, Benjamin T. Zhou,Jin-Xin Hu, Ying-Ming Xie, Binghai Yan, and K. T.Law, “Giant nonlinear Hall effect in strained twisted bi-layer graphene,” arXiv e-prints , arXiv:2010.08333 (2020),arXiv:2010.08333 [cond-mat.mes-hall].[54] Robert W Boyd, Nonlinear optics (Elsevier, San Diego,United States, 2003).[55] J. E. Moore and J. Orenstein, “Confinement-InducedBerry Phase and Helicity-Dependent Photocurrents,”Phys. Rev. Lett. , 026805 (2010), arXiv:0911.3630[cond-mat.mes-hall].[56] Inti Sodemann and Liang Fu, “Quantum Nonlinear HallEffect Induced by Berry Curvature Dipole in Time-
Reversal Invariant Materials,” Phys. Rev. Lett. ,216806 (2015), arXiv:1508.00571 [cond-mat.mes-hall].[57] B. Yan and C. Felser, “Topological Materials: WeylSemimetals,” Annu. Rev. Condens. Matter Phys. ,337–354 (2017), arXiv:1611.04182 [cond-mat.mtrl-sci].[58] N. P. Armitage, E. J. Mele, and A. Vishwanath, “Weyland Dirac semimetals in three-dimensional solids,” Rev.Mod. Phys. , 015001 (2018), arXiv:1705.01111 [cond-mat.str-el].[59] Tonatiuh Rangel, Benjamin M. Fregoso, Bernardo S.Mendoza, Takahiro Morimoto, Joel E. Moore, and Jef-frey B. Neaton, “Large bulk photovoltaic effect and spon-taneous polarization of single-layer monochalcogenides,”Phys. Rev. Lett. , 067402 (2017).[60] David Vanderbilt, Berry Phases in Electronic StructureTheory: Electric Polarization, Orbital Magnetization andTopological Insulators (Cambridge University Press, 2018).[61] D. Xiao, M.-C. Chang, and Q. Niu, “Berry phase effectson electronic properties,” Rev. Mod. Phys. , 1959–2007 (2010), arXiv:0907.2021 [cond-mat.mes-hall].[62] J. E. Sipe and J. Zak, “Geometric phase for electricpolarization along ’rational’ directions in crystals,” Phys.Lett. A , 406–412 (1999).[63] Aaron M. Burger, Radhe Agarwal, Alexey Aprelev, Ed-ward Schruba, Alejandro Gutierrez-Perez, Vladimir M.Fridkin, and Jonathan E. Spanier, “Direct observationof shift and ballistic photovoltaic currents,” Sci. Adv. ,eaau5588 (2019). [64] Fengcheng Wu, Timothy Lovorn, Emanuel Tutuc, andA. H. MacDonald, “Hubbard model physics in transitionmetal dichalcogenide moir´e bands,” Phys. Rev. Lett. , 026402 (2018).[65] A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S.Novoselov, and A. K. Geim, “The electronic propertiesof graphene,” Rev. Mod. Phys. , 109–162 (2009),arXiv:0709.1163.[66] Vitor M. Pereira, A. H. Castro Neto, and N. M. R. Peres,“Tight-binding approach to uniaxial strain in graphene,”Phys. Rev. B , 045401 (2009).[67] F. Guinea, M. I. Katsnelson, and A. K. Geim, “Energygaps and a zero-field quantum Hall effect in grapheneby strain engineering,” Nat. Phys. , 30–33 (2010),arXiv:0909.1787 [cond-mat.mes-hall].[68] Nguyen N. T. Nam and Mikito Koshino, “Lattice re-laxation and energy band modulation in twisted bilayergraphene,” Phys. Rev. B , 075311 (2017), 1706.03908.[69] M. Monteverde, C. Ojeda-Aristizabal, R. Weil, K. Ben-naceur, M. Ferrier, S. Gu´eron, C. Glattli, H. Bouchiat,J. N. Fuchs, and D. L. Maslov, “Transport and ElasticScattering Times as Probes of the Nature of ImpurityScattering in Single-Layer and Bilayer Graphene,” Phys.Rev. Lett. , 126801 (2010).[70] M. Tinkham, Group Theory and Quantum Mechanics (Dover Publications, 2003).
Supplemental Material
I. CONTINUUM MODEL CONSTRUCTION
The continuum model for TBG [33, 48] is constructed by joining together two monolayer graphene layers at zeroeffective separation between them. We choose the following real space unit vectors, for each graphene layer, a = √ d (cid:32) , √ (cid:33) , a = √ d (cid:32) , − √ (cid:33) (S1)In this description, the A and B sublattices are located, respectively, at v A = (0 , v B = d (0 , b = 4 π d (cid:32) √ , (cid:33) , b = 4 π d (cid:32) − √ , (cid:33) , (S2)and the Brillouin zone corners hosting the low energy states are at K u = ± π d (cid:16) √ , (cid:17) . Within the BM model, thebilayer system is symmetric under C by construction, and inversion and time-reversal (when both K u valleys ofthe original graphene monolayers are included). In order to break inversion symmetry, we introduce a coupling toa substrate (for example, hBN), which lifts inversion symmetry but leaves C symmetry intact. This allows for afinite shift current which is entirely independent of the Berry curvature dipole, because the latter is set to zero by C symmetry. The Hamiltonian of the bilayer system is therefore given by, H = H t + H b + H tb , (S3)0where t,b,int denote the top, bottom layer, and interlayer hopping respectively. The top layer has the followingcontinuum Hamiltonian, for a given momentum q , H t,u ( q ) = (cid:126) v f (cid:88) s ua † t,s,u ( q ) R + q · σ a t,s,u ( q ) , (S4)where s, u designate the spin and valley degrees of freedom; i.e., s = ↑ , ↓ , u = ±
1, for the
K, K (cid:48) valleys of theoriginal graphene monolayers. a t/b ,s,u is the annihilation operator for an electron with spin s, and valley u, onthe A/B sub-lattices of the top/bottom layers. R ± is the rotation matrix for the top/bottom layers, given by: R ± = R (cid:0) ± θ (cid:1) = cos( θ ) ∓ iσ y sin( θ ), acting on the sub-lattice space, with θ ≈ . ◦ denoting the twist angle. TheFermi velocity is (cid:126) v f = 0 . σ = ( σ x , σ y , σ z ) are Pauli matrices. The momentum q is measuredrelative to the valley cenetered at K u . Inversion symmetry breaking and uniaxial strain are introduced in the bottomlayer. Throughout this work, the strain is applied along the zigzag direction of the bottom graphene sheet. TheHamiltonian of the bottom layer then takes the form, H b,u ( q ) = (cid:126) v f (cid:88) s a † t,s,u ( q ) (cid:0) R − (1 + (cid:15) )( q + u A ) · σ + ∆ σ z (cid:1) a t,s,u ( q ) , (S5)Here, (cid:15) is the uniaxial strain matrix, which has the form (cid:15) = (cid:15) (cid:0) − ν (cid:1) . A is the pseudo-gauge field resultant fromthe application of strain [66, 67]. ∆ = 17 meV, is the staggered potential generated by alignment with an hBN layer. (cid:15) = 0 , .
1% representing the strain-less and strained cases, respectively. With strain H tb is changed accordingly, asdiscussed in Methods. II. METHODS
The Bistritzer-MacDonald continuum model contains an inter-layer coupling term, which couples two momenta q , q (cid:48) ,if q − q (cid:48) = { q ,u , q ,u , q ,u } , where q ,u = | q | (0 , − q ,u = | q | (cid:16) √ , (cid:17) , q ,u = | q | (cid:16) − √ , (cid:17) , are the Moir´e latticevectors with q = π sin( θ )3 √ d , and d = 1 . H tb term in Eq. S3 as, H tb = (cid:88) q , q (cid:48) ,s,u a † t,s,u ( q )( T ,u ( q , q (cid:48) )+ (S6) T ,u ( q , q (cid:48) ) + T ,u ( q , q (cid:48) )) a b,s,u ( q (cid:48) ) . In the presence of the strain defined in the main text, the coupling matrices (acting on the valley index) become, T ,u = t (cid:18) (cid:19) δ q − q (cid:48) , q ,u (S7) T ,u = t (cid:32) e − iu π ( − (cid:15) ν ) e iu π ( − (cid:15) ν ) 1 (cid:33) δ q − q (cid:48) , q ,u (S8) T ,u = t (cid:32) e iu π ( − (cid:15) ν ) e − iu π ( − (cid:15) ν ) 1 (cid:33) δ q − q (cid:48) , q ,u . (S9)Throughout, we take t = 0 .
33 eV. Accordingly, the lattice vectors of the Moir´e superlattice are deformed in thepresence of strain. These have the form, q ,u = u π √ d (cid:0) (cid:15) cos θ , (2 + (cid:15) ) sin θ (cid:1) (S10) q ,u = u π d (cid:16) √ (cid:15) cos θ − − (cid:15)ν ) sin θ , (S11)3 (cid:15)ν cos θ − √ (cid:15) ) sin θ (cid:17) q ,u = − u π d (cid:16) θ (2 − ν(cid:15) ) − √ (cid:15) cos θ , (S12)3 ν(cid:15) cos θ + √ (cid:15) ) sin θ (cid:17) Figure S1. The Brillouin zone of the top (red) and bottom (blue) graphene sheets, with and without strain. Left: Twoundistorted ( ε = 0) grahpene layers are rotated one with respect to the other, forming a folded mini Brillouin zone (mBZ),when rotated by an angle ± θ/
2. The K point of each layer shifts to K t and K b for the top and bottom layers respectively.The separation between them is denoted by q = R + K u − R − K u . When the mBZ is refolded onto the center of the originalBrillouin zone, K t/b = | q | (cid:16) √ , ± (cid:17) . Middle: Upon introduction of uniaxial strain on the bottom layer, the Brillouin zonedeformes by expanding one, and contracting in the other direction. The K points transform according to Eq. S13. Right: whenthe strained Brillouin zone is rotated with respect to an unstrained one, a deformed mBZ is formed, as shown here. Consequently,tunneling vectors which depend on the positions of K t/b in the mBZ are modified, as shown in Methods. The angle formedbetween the unstrained q and the strained vector q (cid:48) is given by θ = cos − (cid:16) | q | | q (cid:48) | (cid:17) ≈ . o . The Dirac points transform under strain as, ¯ K u = (1 − ε T ) K u − u A (S13)We introduce a pseudo-gauge field which stems from the underlying two-center approximation for the tunneling matrix[68]. It is given by, A = − β(cid:15)d (1 + ν,
0) (S14)with ν = 0 . β = 1 .
57 obtained from monolayer graphene. Finally, the diagonalization of the Hamiltonian H in Eq.S3 is accomplished by recasting it in the form, H = (cid:88) q ,s,u A † s,u ( q ) h u ( q ) A s,u , (S15)where now A s,u ( q ) = [ a b,s,u ( q ) , a t,s,u ( q + q ,u ) , a t,s,u ( q + q ,u ) , a t,s,u ( q + q ,u )] T is the infinite-component operatorvector satisfying the constraints on momentum transfer. h u has the following truncated structure, after applying themomentum transfer relations ensuring non-zero T tunnelling, h u ( q ) = H b,u ( q ) T ,u T ,u T ,u T † ,u H (1) t,u T † ,u H (2) t,u T † ,u H (3) t,u , (S16)with H ( i ) t,u = H t,u ( q + q i,u ). In this work, we used 81 sites in reciprocal space for the construction of the hamiltonian,which results in a Hamiltonian which is 324 × ×
600 in the ( k x , k y ) plane of the mini Brillouin zone. 10 frequency point samplings in ω are carried out uniformly. Convergence is checked against the case with C symmetry, where verification is done bycomparing σ xx ; x with − σ yy ; x ; and σ yy ; y with − σ xx ; y . All equalities were verified to within 5%. The delta functions ofEqs. (7),(8) are broadened with a width Γ = 0 . T = 0K and Γ = 0 . T = 300K. This corresponds totransport lifetimes observed in bilayer graphene in the clean limit [69].2 - - Figure S2. Conductivity σ xx ; y as a function of frequency, with ε = 0 for 3 values of the chemical potential. Compare with σ yy ; x in Fig. 3. Here, σ xx ; y is negative for both µ = 6 . , − . σ xx ; y = − , − µ AnmV − ,respectively. For µ = 0, σ xx ; y = 28 µ AnmV − . A sizeable conductivity, | σ xx ; y | > µ AnmV − is obtained for a wide range offrequencies between ω = 1 . − . III. SYMMETRY PROPERTIES OF THE RESPONSE
In the main text, we observe that the response tensor σ ab ; c has only two independent components, and that theBerry curvature dipole, ∂ a Ω bc vanishes. We proceed to prove this. The generator of C z is given by [70], M = − − √ √ −
00 0 1 . (S17)The Berry curvature dipole (BCD), D abc is a gauge invariant material property of the system, and has the form D abc = ∂ k a Ω bc , making it a rank-3 pseudotensor. Firstly, we observe that this quantity is anti-symmetric in ( b, c ). ApplyingNeumann’s principle, we enforce D abc = (cid:80) αβγ M aα M bβ M cγ D αβγ . Using the anti-symmetry of D abc = − D acb , we focusonly on non-trivial components, D axy = (cid:80) αβγ M aα M xβ M yγ D αβγ = (cid:80) aα M aα (cid:0) D αxy − D αyx (cid:1) = (cid:80) aα M aα D αxy .Since a = x, y , we obtain the following set of equations, D xxy = − D xxy − √ D yxy , D yxy = √ D xxy − D yxy . Thisset admits only the solution D xxy = D yxy = 0, as required, demonstrating that the BCD is zero, under C z . For thegeneral rank-3 symmetric conductivity tensor σ ab ; c , we derive analogous symmetry constraints under C z . This resultsin two independent components overall, σ xx ; x = − σ yy ; x = − σ xy ; y = − σ yx ; y (S18) σ yy ; y = − σ xx ; y = − σ yx ; x = − σ xy ; x . (S19)We further note that under linear-polarized light, the conductivity tensor exhibits a special permutation symmetry, σ ab ; c = σ ba ; c . IV. ADDITIONAL DATA FOR THE TRANSVERSE COMPONENTS
For completeness, we provide the remaining transverse component σ xx ; y of the shift current, which agrees qualitativelyand in part quantitatively with the component σ yy ; x discussed in the main text. We recall that in general, with C z symmetry the conductivity tensor has only 2 independent components, therefore the longitudinal components can bededuced straightforwardly from the data presented here.With ε = 0, the remaining independent component is σ xx ; y , also a transverse component. This is presented in Fig.S2, for 3 chemical potential values. Note that the contribution of Eq. (7) is zero, due to C z symmetry.With finite strain, all terms in Eq. (3) of the main text contribute to the current. For consistency, we again presentthe transverse component, σ xx ; y , for 3 chemical potential values, in the same way as in Fig. 3 of the main text. Withthe introduction of strain, the Berry curvature dipole contribution is no longer zero (Fig. S3a). Although the netconductivity is not substantially enhanced by the introduction of strain (cf. Fig. S3c), a broad resonance in σ xx ; y - - - (a) (b) (c) Figure S3. Contributions to the conductivity σ xx ; y , with strain ε = 0 . R shift + K shift , as in Eq. (8). (c) Total conductivity. While the introduction of strainproduces a giant Berry curvature dipole, the magnitude of Eq. (8) also increases albeit with opposite sign. Consequently, thetotal conductivity remains comparable to the unstrained case. For µ = − . , . , . σ xx ; y = − , , − µ AnmV − , respectively. We note that for µ = 6 . ω = 2 − - - (a) (b) (c) Figure S4. Contributions to the conductivity σ xx ; y with ε = − . K shift + R shift term. (c) Total conductivity. . appears at chemical potential µ = 6 . σ yy ; x for strain with negative (i.e., compressive) magnitude, (cid:15) = − ..