Optical conductivity, Drude weight and plasmons in twisted graphene bilayers

We numerically calculate the optical conductivity of twisted graphene bilayers within the continuum model. To obtain the imaginary part, we employ the regularized Kramers-Kronig relation allowing us to discuss arbitrary twist angles, chemical potential and temperature. We find that the Drude weight $D$ as function of the chemical potential $\mu$ closely follows the shell structure of twisted bilayer displayed by the density of states. For certain angles, this results in a transport gap D=0 at finite $\mu$. We also discuss the loss function which, for low doping, is characterized by acoustic interband"plasmons"and transitions close to the van Hove singularities. For larger doping, the plasmon mode of decoupled graphene bilayer is recovered that is damped especially for small wave numbers.


Introduction
Since the realization that neutral graphene monolayers display a universal optical absorption of πα ≈ 2.3% in the infrared and visible spectrum (with α ≈ 1/137 the fine structure constant), [1,2,3,4] interest in the interactions between light and graphene has grown rapidly. [5,6,7] A remarkable aspect of this absorption, other than its suggestive relation to the universal constant α, is its high value, given the ultimate thinness of a graphene monolayer. In fact, it corresponds to an attenuation factor that is hundreds of times larger than that of good metals, such as gold and can be even enhanced up to 10% for photon energies around 4.5eV due to transitions close to the van Hove singularity. [8,9] The unique feature of neutral graphene responsible for its strong optical response is the existence of low energy particle-hole interband excitations with zero momentum, by virtue of the gapless Dirac spectrum, that act as an effective sink for photons of arbitrarily low energy. By the same token, all thermally excited plasmons in undoped graphene are damped by this excitation density, but nevertheless remain welldefined. [10] The unique bipolar tunability of graphene further allows one to adjust these optical properties quite simply, by shifting the Fermi energy µ away from the Dirac point, be it electrostatically or through chemical doping. This opens a gap of energy 2µ in the particle-hole excitations, making graphene transparent below this energy, useful, e.g., for broadband optical modulation. [11] It also allows for the emergence of undamped plasmons. [12,13,14,15,16] The spectral features involved in the unusual optical properties outlined above are also present in graphene bilayers. [17] These are often classified according to their stacking, which in the absence of strain have three possible (meta) stable configurations: AB'-type (or BA', also known as Bernal stacking), AA'-type, and twisted. The former two are minimal stackings (unit cell with only four carbon atoms) in which half or all sites in opposite layers are vertically aligned, respectively. The twisted bilayer is the most generic type, and is defined by an relative rotation of angle θ between the two layers. Each of the three stacking types exhibits a different low-energy electronic structure, although all three react strongly to light like the monolayer. The twisted bilayer, however, and in particular the case of low twist angles (around and below 1 o ), has the richest optical response, a consequence of its complex low energy spectral properties. The response, moreover, is highly tuneable through the angle θ, an additional control parameter not available in the other graphenes. As a striking example, the strong frequency-dependence of twisted bilayer conductivity and thus reflectivity leads to a visible coloration under white light, that varies with angle θ, of flakes on 100 nm thick SiO 2 , as recently observed experimentally. [18] In this work we study the optical response of twisted graphene bilayers at different twist angles and temperature responsible for this effect. We analyze its optical conductivity at and away from neutrality, the Drude weight and plasmon spectrum. Some of these aspects, specifically the real part of the conductivity, have been studied in recent works for large-angle twisted bilayers. [19,20] We will make connection with these, and extend them into the low angle regime. The paper is organized as follows. In Section II, we introduce the continuum model of twisted graphene bilayer and discuss the electronic structure. In Section III, we introduce the Kubo formalism and also comment on the regularization of the Kramers-Kronig relation for the conductivity, necessary in the case of unbounded Dirac Fermions. In Section IV and V, we present our results for the optical conductivity and the related Drude weight for different angles, chemical potential and temperature. In Sec. VI, we discuss the plasmonic spectrum based on the local conductivity by plotting the energy loss function and close with conclusions.

Continuum model of a twisted graphene bilayer
Twisted graphene bilayers are a peculiar structure from a crystallographic point of view, since for a generic (non-commensurate) angle the lattice is not exactly periodic, although crystalline order within each layer is extraordinarily robust. This is a consequence of the highly anisotropic binding forces in these systems, also characteristic of graphite. In principle only commensurate twist angles θ = θ mn , where [21] cos(θ mn ) = 3m 2 + 3mn + n 2 /2 3m 2 + 3mn + n 2 for integer m, n, are amenable to band structure calculations by application of the Bloch theorem. It turns out, however, that the bandstructure at low energies depends continuously on the value of the angle θ n,m itself, not on the specific n and m (note that the set of possible commensurate angles θ mn is dense). This is a consequence of the fact that the physical interlayer distance d > 0.3 nm is greater than the intralayer carboncarbon distance a cc = 0.14 nm. In this regime, all properties associated to precise commensurability are efficiently averaged out using microscopic tight-binding models. Consequently, a simple, low energy continuum theory can be formulated that may be applied to an arbitrary twist angle θ, and obviates all crystallographic complications to a very good approximation. [22,21] The first ingredient in the continuum model is the momentum shift between Dirac cones in different layers due to the rotation, ∆K = 2|K| sin(θ/2), where |K| = 4π/3a 0 and a 0 = √ 3a cc is the monolayer Bravais period. The second non-trivial ingredient is the smooth spatial modulation in the interlayer coupling. The twist θ produces a triangular Moiré pattern in the local stacking, of period L M = a 0 /[2 sin(θ/2)]. Hence, the coupling between layers becomes position dependent, alternating periodically in space between AB', BA' and AA'. The modulation is smooth on the scale a 0 , so the two valleys in each layer remain decoupled in the twisted bilayer. For each valley, the continuum model Hamiltonian reads [23,21,24] This matrix is written in the local basis A, B, A , B , comprised of the two sites in the unit cell of each of the two layers. Here, the Fermi velocity is v F = ( √ 3/2)a 0 γ 0 / ≈ 9 · 10 5 m/s, and Π ± = [k x + i(k y ∓ ∆K/2)] e i±θ/2 . The function V (r) = 1 3 γ ⊥ 3 i=1 e ig i ·r describes the periodic spatial variation of the interlayer coupling, with γ ⊥ ≈ 0.33eV. The Moiré lattice vectors are denoted by a 1,2 (so L M = |a 1,2 |), and conjugate momenta are g 1,2 , such that g i ·a j = 2πδ ij (we also define g 3 = 0). The AA' sublattice is centered at the origin, with the AB'/BA' sublattices offset given by ±δr = (a 1 − a 2 )/3.
We first review the bandstructure phenomenology of the twisted bilayer. The continuum model has two distinct regimes, denoted as large and small angles, each with very different phenomenology. The large angle regime, relevant for θ greater than a few degrees, is characterised by two low energy Dirac cones per valley at the superstructure K and K , with a reduced Fermi velocity [25,26,27] where α = γ ⊥ /6 M is an expansion parameter that represents the dimensionless interlayer coupling, and M (θ) = v F ∆K/2 = (2π/3) v F /L M is the energy scale associated to the Moiré period L M . Also, a low energy van-Hove singularity at energy vH ≈ M − γ ⊥ /3 is created [23,28,29,30,31]. At an angle θ m 1 ≈ 1.05 o , known as the "first magic angle", v * F becomes zero, the van-Hove singularity degenerates into a narrow miniband around zero energy. This is the onset of the small angle regime. From this point on, the electronic structure becomes considerably complex. The Fermi velocity becomes finite again, as the miniband moves away from zero energy. Higher energy van-Hove singularities appear and approach zero at subsequent magic angles θ m 2,3,4... , where the Fermi energy vanishes once more.
In Fig. 1, the density of states (DOS) is shown for various angles labeled by the index i obtained from Eq. (1) with m = n + 1 = i. The curves are normalized by ρ 0 = 4γ 0 such that the DOS of two uncoupled graphene layers is just a straight line with slope 1 (dashed line). Note that the DOS is almost but not exactly symmetric for positive and negative energies as a consequence of the underlying (not electron-hole symmetric) Hamiltonian.
In the upper panel, we show the DOS for "large" angles θ with i = 10 (blue) and i = 15 (red). Around the neutrality point, the spectrum is characterized by a van Hove singularity followed by a step singularity, forming a closed shell. For small energies the system is thus approximately described by a hexagonal tight-binding model with renormalized band-width. This shell structure is somehow repeated for larger energies until it fades away and the DOS for decoupled bilayer graphene is recovered. Note that the weight of the van Hove singularity is slightly decreasing with increasing For small angles around the first magic angle θ m 1 = 1.05 o (i = 31, lower panel in Fig. 1), more distinct van Hove singularities appear and the height of the first van Hove singularity is strongly enhanced with increasing i. A close-up of this energy regime in shown on the left hand side of Fig. 7.

Linear response theory
The above phenomenology in the electronic bandstructure has an associated signature in the optical response of the system as the angle decreases. The basic quantity governing the linear optical response of any electronic system is the optical conductivity σ αβ (ω) (where α, β = x, y), which gives the electric current arising in response to a weak and uniform, time-dependent electric field of frequency ω, The optical conductivity may be written, using linear response theory, in terms of the current response χ jαj β (q, ω) at zero-momentum, where, as is conventional, δ denotes a positive infinitesimal. The current response function χ jαj β consists, in general, of two terms: the diamagnetic and the paramagnetic contribution. [32] The general formalism for obtaining the diamagnetic contribution in a tight-binding model is outlined in Ref. [33], however, in the continuum model for the twisted bilayer, as in the monolayer, the diamagnetic contribution vanishes due to the fact that ∂ 2 H/∂k 2 α = 0. Hence, we have simply the paramagnetic contribution, given by the Kubo formula, Here, g s = g v = 2 are the spin and valley degeneracies. States |m, k are eigenstates of H in subband m and of momentum k in the first Brillouin zone of the superstructure. Their eigenenergies are mk and n F is the Fermi function. The paramagnetic current operator is j α = −∂H/∂k α . In the continuum model with local interlayer coupling, this current does not mix layers and thus independent of the coupling γ ⊥ . For a system with time-reversal symmetry, the conductivity tensor is diagonal and due to the threefold rotational symmetry of (twisted) bilayer graphene proportional to the unity tensor. We will, therefore, drop the tensor indices in what follows. The real part of the conductivity σ(ω) represents the power dissipation of the system upon incidence of light at frequency ω, while its imaginary part represents inductive retardation of the current with respect to the driving field. From Eq. (5), it is obvious that the conductivity is ill-defined at ω = 0 and δ → 0 if χ jj (q = 0, ω = 0) = 0. It is thus customary to split up the conductivity into a regular part and a term containing the delta-singularity: with the Drude weight (or charge stiffness) D defined as Re χ jj (q = 0, ω) = lim ω→0 ω Im σ(ω) .
Above we used the fact that the imaginary part of the current response is an odd function due to causality. We thus have Im σ = Im σ reg and will drop the subindex when confusion is unlikely.

Graphene monolayer
In the case of a graphene monolayer, described within a continuum Dirac model, we have the famous zero temperature result (µ is the chemical potential with respect to the neutrality point) Re σ reg (ω < 2µ) = 0.
This universal value σ 0 for the Dirac gas includes valley and spin degeneracy factors, and is valid for energies ω ∼ 1eV well below the γ 0 ∼ 3eV van-Hove singularity of the monolayer. As ω approaches this scale, a discrete model must be employed to compute the monolayer σ 0 (see Ref. [4] for general results). The imaginary part of the conductivity, in the same limit of low frequencies, reduces to [4,3] Im σ(ω where µ is the Fermi energy, measured from the Dirac point (charge neutrality). Using Eq. (9), we have a Drude weight which vanishes at neutrality. Note that the above expression only scales with the chemical potential and no other material constant (Fermi velocity) is involved. The Drude weight for twisted bilayer as function of µ will thus be independent of the twist angle in the regime where the electronic spectrum can be described by a (renormalized) Dirac model. This breaks down close to the first magic twist angle with the emergence of a flat band.

Twisted bilayer
In the more general case of a twisted bilayer, the computation of σ(ω) becomes a numerical task and it is usually simpler to compute Re σ reg than Im σ, since, by virtue of Eq. (6), the former reduces to an integral constrained to eigenstates within ω of the Fermi energy The real part of the conductivity differs from that of two decoupled monolayers only at frequencies smaller than or comparable to the interlayer coupling ω < γ ⊥ . Hence we have Re σ reg (ω γ ⊥ ) ≈ 2σ 0 . This result will be useful in the following discussion to obtain a regularized expression for the imaginary part of the conductivity.
To compute the imaginary part, and from that, the Drude weight D through Eq. (9), one could in principle employ Eq. (6). This is problematic, since it involves a summation, even for small ω, over the whole band, which in the continuum model extends to infinity. We will, therefore, make use of the the Kramers-Kronig (KK) relation, applicable to any response function. In terms of the conductivity, this yields where P denotes the Cauchy principal value. There is apparently no advantage to using Eq. (13) as compared to Eq. (6), since the integral above extends to infinity, too. Moreover, the integral is ill defined, since at high frequencies the continuum model yields a constant Re σ reg (ω γ ⊥ ) ≈ 2σ 0 . We will thus define the imaginary part relative to this divergent background following the usual procedure of ultraviolet field theories. In fact, this is one way of obtaining Eq. (14) and was discussed more generally in the context of the f-sum rule in Ref. [34]. The low-energy properties are not modified by this regularization and we can simply subtract the ground-state expectation value of the two independent graphene layers at zero doping. Since the imaginary part for the neutral system is zero for the continuum model, the final regularized definition of Im σ simply reads Numerically, Eq. (3.2) can be evaluated by introducing a finite cutoff Λ for which Re σ reg (Λ) ≈ 2σ 0 . Integrating the second term involving the "ground-state value" 2σ 0 then yields the final formula: Note that this expression is independent of Λ as it was obtained from Eq. (3.2). The extra terms compared to the original KK-relation (13) can be interpreted as the diamagnetic contribution that naturally arises in the Dirac model with finite band cutoff.
Using this expression, the Drude weight D in Eq. (9) may be written as In terms of D, Eq. (14) is then finally rewritten as From Eq. (15), we can also deduce the f-sum rule which is independent of the angle θ, chemical potential µ or temperature T .

Optical conductivity
The conductivity at neutrality µ = 0 is shown in Fig. 2 for two "large" twist angles. Much like the monolayer, the real part (black) exhibits a universal asymptotic value Re σ(ω → 0) = 2σ 0 . It also exhibits signatures of the van-Hove singularities, in form of a dip-peak structure as discussed in detail in Ref. [19]. The imaginary part (red) exhibits similar features, but vanishes at ω → 0, i.e., there is no Drude weight at µ = 0. This changes for finite temperature where the real part of σ is suppressed at ω = 0, [35] leading to a finite Drude weight necessary to fulfill the f-sum rule Eq. (17). The Drude weight for finite temperature can be deduced from Eq. (12) with the substitution µ → 2 ln 2k B T . [36] Temperature can further provoke a shift of spectral weight close to the van Hove singularity (see the circled region). Also in the case of finite doping µ = 0, see Fig. 3, the imaginary part of the conductivity develops a ∼ D(µ)/ω low frequency divergence, while the real part becomes suppressed below a threshold frequency 2µ for small enough µ (see black and red curves). Unlike for the monolayer, however, this is only true for µ below the van-Hove singularity, since it relies on the linearity of the bandstructure (and density of states) around µ = 0. For chemical potentials close to the van Hove singularity, i.e., µ ≈ M , the optical gap is again filled by spectral weight and the Drude weight reduces consistent with the f-sum rule (blue dashed curves). At finite temperature, the optical gap may persist for small twist angle and the Drude weight can thus even increase, see Fig. 6. For large chemical potentials µ M , the electronic spectrum (DOS) of the continuum model is just that of two uncoupled graphene layers and one could thus expect to recover the conventional graphene conductivity. But the Moiré coupling between layers induces additional damping, since it opens up the possibility of more interband dissipation channels, i.e., any particle hole excitation with momentum equal to a multiple of the momenta of the Moiré lattice. In order to understand the basic features of the conductivity on the basis of the uncoupled system, we can thus attribute a finite in-plane wavenumber q ∼ n∆K to the incoming photon (n > 0), see Fig. 4. The two-particle spectrum of interband and intraband excitations then predicts dissipation at low energies (black) and a broadened transition (blue) from the optical gap to the universal value 2σ 0 of the conductivity.
In Fig. 5, the optical conductivity at the first magic angle i = 31 is shown for various chemical potentials (solid lines) and compared to Re σ reg for i = 10 (dashed blue  Fig. 2(a). Note that Im σ(ω) has been multiplied by ω to reveal its ω → 0 asymptotics.  Fig. 2. For comparison, also Re σ for (i = 10) is shown (dashed blue line). Dotted lines resemble the limit for uncoupled layers. Note that Im σ(ω) has been multiplied by ω to reveal its ω → 0 asymptotics. line). For low chemical potential, the optical response differs considerably for small and large angles, but for increasing µ the dependence on the twist angle becomes less and Re σ is characterized by a large Moire-broadening of the optical gap around ω = 2µ. Still, there is the pronounced spectral weight at ω = vH almost independent of µ, which is especially striking for larger angles. This feature cannot be explained from Fig. 4 and is twist angle dependent.

Drude weight
Considering only intraband transitions of one unbounded conduction band, the Drude weight or charge stiffness is proportional to the chemical potential, independent of its dispersion relation. For parabolic band electrons, D is thus simply related to the charge density, D = e 2 n/m, whereas for Dirac Fermions one gets D ∼ √ n, see Eq. (12). In systems with a finite band-width, we expect a non-monotonous behavior expressing the the fact that the Fermi energy is first related to electrons and later to holes when sweeping the chemical potential through the band. For the one-dimensional tightbinding model, one finds, e.g., D ∼ 1 − (µ/2γ 0 ) 2 with 2γ 0 the band width. For the two-dimensional hexagonal tight-binding model, the system is characterized by a van Hove singularity at µ = γ 0 around which the Drude weight shows a maximum. [35] For twisted bilayer graphene described by the continuum model, we expect a combination of both features, because of the unbounded linear dispersion of Dirac Fermions and the superlattice structure of the Moiré-pattern. From Fig. 1 we see that the density of states is characterized by a van Hove singularity followed by a "real" steplike singularity which is repeated to higher doping levels. This pseudo-gap structure divides the spectrum in several bands and will be followed by the Drude weight. Nevertheless, these modulations are superposed onto a constant linear background as suggested by the unbounded continuum model.
The µ dependence of the Drude weight D(µ) in units of 2D 0 is presented in Fig.  6 for two large twist angles. As is apparent, for small chemical potential D(µ) follows twice the monolayer value 2D 0 independent of the twist angle, as expected for Dirac systems. For larger µ, the Drude weight drops and remains less than 2D 0 until it converges to the uncoupled bilayer system for µ γ 0 . As outlined above, D(µ) exhibits dips that accurately correlate with dips in the density of states (DOS), also shown in Fig. 6 (red). The resulting picture suggests the above mentioned approximate shell structure, albeit (in general) without gaps between shells, just as the suppressed DOS. The first shell can be fairly well approximated be the hexagonal tight-binding model.
At finite temperature T = 300K (and also for disorder), the Drude weight is smeared out for large twist angle i = 10, thus reaffirming the interpretation of the repeated shell-structure. Nevertheless, the shell structure becomes less pronounced with smaller θ which is also reflected by its temperature response.
For small angles, the Drude weight provides a sensitive fingerprint for the breakdown of Dirac physics that occurs around the magical angles due to the emergence of flat bands. This can be seen on the right hand side of Fig. 7 for small angles above (i = 26) and below (i = 36) the critical angle. The Drude weight becomes larger than the Drude weight of two monolayers before it drops indicating the band width of the flat band. At the critical angle (i = 31), this band width becomes zero, and interestingly, the scaling of D 0 is recovered. Note that for i = 26, we find a transport gap D = 0 at finite chemical potential µ = 0.06γ 0 .
The increase of the band width for twist angles larger than θ m 1 can be nicely appreciated from the Drude weight. This becomes considerably more difficult from the DOS, displayed on the left hand side of Fig. 7. The reason for this is that the Drude weight is obtained from integrating over the whole spectral region and not only from the low energy band structure. Nevertheless, the weight of the peak at zero energy is maximal at the magic angle (red curve), indicating its singularity.

Plasmons
A complementary aspect to the optical absorption is the plasmon excitation spectrum of the electron gas. Conceptually, a plasmon is an oscillating charge density mode (hence accompanied by an oscillating electric potential), which is sustained by the Coulomb interactions between electrons. The dispersion relation of these modes is given by the zeroes of the dielectric function (q, ω) = 0, that relates the screened electric potential V sc (r, t) to an externally applied potential V ext (r, t), A plasmon is, therefore, a finite solution V sc (q, ω) that requires no external driving, i.e., it is self-sustained. A Landau-damped plasmon corresponds to a peak in 1/ (q, ω), and decays due to dissipation into particle-hole excitations. The dielectric function (q, ω) may be computed within the Random Phase Approximation (RPA), which in two dimensions yields where χ ρρ is the charge susceptibility, and is the dielectric constant of the environment surrounding the two-dimensional electron gas. From the continuity equation between charge and current, it follows that [32] χ ρρ (q, ω) = q 2 ω 2 χ L jj (q, ω) where χ L jj is the longitudinal part of the current-current correlator in Eq. (6). In the long-wavelength limit, we may thus relate RPA directly to σ(ω) using Eq. (5), where α g = e 2 /(4π 0 v F ) ≈ 2.2 is graphene's fine structure constant. To compute ω p for small frequencies, we may expand Im σ in Eq. (16) by the Drude term, This results in the dispersion of a classical 2D plasmon For two decoupled monolayers, one just replaces D = 8µ/π σ 0 above, which is shown as dotted curved line in Fig. 8.

Plasmons in twisted Bilayer
In twisted graphene bilayer there may exist four types of plasmonic modes or resonances: i) undamped (conventional) graphene plasmons for chemical potentials with µ M , ii) interband "plasmons" present also for µ = 0, iii) damped plasmons for large chemical potential µ M , iv) transverse plasmons. Below, we will briefly outline the main properties of these different types: 6.1.1. Undamped plasmons At low doping, there exists a low energy window with Re σ = 0 that gives rise to undamped plasmons governed by the Drude weight D. In this regime, the twisted bilayer can be described by linear Dirac Fermions and D(µ) follows the Drude weight of the uncoupled system 2D 0 , as shown in Fig. 6. This energy window becomes smaller for decreasing twist angle since the van Hove singularity moves closer to the neutrality point and is only relevant for large angles.
6.1.2. Interband plasmons In monolayer graphene, a peak in the loss function associated to π → π * transitions around the van Hove singularity and with linear dispersion was first predicted by DFT-studies, [37] and later experimentally observed in suspended graphene by electron energy-loss spectroscopy (EELS). [38] Within the hexagonal tight-binding model and RPA, i.e., without including correlation or renormalization effects, no zero of the dielectric function is obtained. [39] The absorption peak would thus be merely due to interband transitions enhanced by a band-structure effect; nevertheless, in bi-or multilayer, RP A (q, ω) becomes zero around the Mpoint. [40] In the twisted graphene bilayer, there exist several van Hove singularities such that interband plasmons (or absorption enhancement) are to be expected.

Damped plasmons
In the regime of large chemical potential, µ M , the conventional intraband plasmon is recovered, albeit Moire-damped. This can be deduced from Fig. 5 which will be used to discuss the plasmon dispersion in terms of the loss function. The dispersion only depends slightly on the twist angle and becomes well-defined for large µ, extending into the Landau damped region just as for the monolayer.

Transverse Plasmons
There is also the possibility for transverse or transverse electromagnetic (TE) plasmons modes given by the condition [41] 1 − ωiσ(ω) For TE modes, the plasmon dispersion is closely pinned to the light cone which justifies our local approximation (q → 0) and there are (damped) solutions to Eq. (20) if Im σ < 0. Since they can only exist in twisted bilayer surrounded by a homogeneous dielectric background, [42] we will not discuss them any further.

Loss Function
The loss function, defined as S(q, ω) = −Im −1 (q, ω) is a measure of the spectral density of plasmon excitations. Thus, a sharp peak in S(q, ω) reveals a long-lived (transverse magnetic, TM) plasmon of that frequency and momentum. A proper undamped plasmon, defined by (q, ω) = 0, corresponds to a delta peak in S(q, ω), since −1 (q, ω) satisfies Kramers-Kronig relations. In Fig. 8, we present S(q, ω) within the long-wavelength RPA approximation of Eq. (18) at the first magic angle θ m 1 = 1.05 o (i = 31) for different chemical potential with = 1. Similar results are obtained for larger angles, as can be deduced from Fig.  (5). In all cases, the Dirac dispersion ω = v F q and the ω = 2E F − v F q from Fig. 4 is shown as dashed lines, indicating the onset of intra-and interband transitions. For zero or low chemical potential (Fig. 8a), the loss function indicates an acoustic interband "plasmon" with sound velocity v s ≈ πα g v F (solid line), ‡ and only for low energies signatures of the van Hove singularities are seen which are more pronounced for larger angles. This linear (interband) mode crosses over to the conventional √ q (intraband) mode, Fig. 8b) and c), and becomes well-defined for large chemical potential, Fig. 8d, although there are slight variation compared to Eq. (19) shown as solid curved line.

Conclusions
We have analyzed the density of states (DOS), optical conductivity, Drude weight, and plasmon spectrum in a twisted bilayer for various angles, doping levels and temperature, and compared them to the uncoupled system. We have found a universal low-frequency conductivity and peaks at finite frequency that correlate to van-Hove singularities in the DOS. These latter features, however, become scrambled at decreasing angles. We have also shown that the Drude weight in twisted bilayers is usually below that of two decoupled monolayers, 2D 0 , and exhibits a shell-like structure, following again van-Hove singularities in the DOS. Only for small angles close to the neutrality point, the Drude weight becomes larger than 2D 0 , signaling the emergence of a flat band. The Drude weight resembles a sensitive quantity to numerically distinguish the different electronic structure at critical angles since it is obtained from integrating over the whole spectrum and is directly related to transport experiments. We have finally discussed how plasmons may become damped by the interlayer Moiré coupling which opens up new dissipation channels. Significant spectral structure survives, however, including undamped plasmons and multiple long-lifetime resonances (interband plasmons) within the extended interband particle-hole continuum.