# Fileset

[2024A00489G_Negative_compressibility_of_charge_islands_within_twisted_bilayer_graphene_main.pdf](https://mdr.nims.go.jp/filesets/1a3051f9-e906-48b3-9d1e-3bd080beb7dc/download)

## Creator

Robin J. Dolleman, Alexander Rothstein, Ammon Fischer, Lennart Klebl, Lutz Waldecker, [Kenji Watanabe](https://orcid.org/0000-0003-3701-8119), [Takashi Taniguchi](https://orcid.org/0000-0002-1467-3105), Dante M. Kennes, Florian Libisch, Bernd Beschoten, Christoph Stampfer

## Rights

©2025 American Physical Society.[In Copyright](http://rightsstatements.org/vocab/InC/1.0/)

## Other metadata

[Negative electronic compressibility in charge islands in twisted bilayer graphene](https://mdr.nims.go.jp/datasets/25caeb45-f3dc-44a6-a579-c84a1ec9446f)

## Fulltext

Negative electronic compressibility in charge islands in twisted bilayer grapheneRobin J. Dolleman,1, ∗ Alexander Rothstein,1, 2 Ammon Fischer,3 Lennart Klebl,4 Lutz Waldecker,1 Kenji Watanabe,5Takashi Taniguchi,6 Dante M. Kennes,3, 7 Florian Libisch,8 Bernd Beschoten,1 and Christoph Stampfer1, 2, †1JARA-FIT and 2nd Institute of Physics, RWTH Aachen University, 52074 Aachen, Germany, EU2Peter Grünberg Institute (PGI-9), Forschungszentrum Jülich, 52425 Jülich, Germany, EU3JARA-FIT and Institut für Theorie der Statistischen Physik,RWTH Aachen University, 52074 Aachen, Germany, EU4I. Institute of Theoretical Physics, University of Hamburg, Notkestrasse 9, 22607 Hamburg, Germany, EU5Research Center for Functional Materials, National Institute for Materials Science, 1-1 Namiki, Tsukuba 305-0044, Japan6International Center for Materials Nanoarchitectonics,National Institute for Materials Science, 1-1 Namiki, Tsukuba 305-0044, Japan7Max Planck Institute for the Structure and Dynamics of Matter,Center for Free Electron Laser Science, Hamburg, Germany, EU8Institute for Theoretical Physics, Vienna University of Technology, 1040 Vienna, Austria, EUWe report on the observation of negative electronic compressibility in twisted bilayer graphenefor Fermi energies close to insulating states. To observe this negative compressibility, we takeadvantage of naturally occurring twist angle domains that emerge during the fabrication of thesamples, leading to the formation of charge islands. We accurately measure their capacitance usingCoulomb oscillations, from which we infer the compressibility of the electron gas. Notably, we notonly observe the negative electronic compressibility near correlated insulating states at integer filling,but also prominently near the band insulating state at full filling, located at the edges of both theflat- and remote bands. Furthermore, the individual twist angle domains yield a well-defined carrierdensity, enabling us to quantify the strength of electronic interactions and verify the theoreticalprediction that the inverse negative capacitance contribution is proportional to the average distancebetween the charge carriers. A detailed analysis of our findings suggests that Wigner crystallizationis the most likely explanation for the observed negative electronic compressibility.I. INTRODUCTIONTwisting two layers of graphene creates a moiré su-percell with enlarged periodicity [1–4], leading to bandinsulating states at full filling (ν = ±4, where ν rep-resents the filling factor denoting the number of chargecarriers per supercell) [5] and flattening of the electronicbands [2, 4]. Near the so-called magic angle of ≈ 1.1◦,the ratio of Coulomb repulsion to kinetic energy becomesmaximized [6], giving rise to correlated insulating statesat fractional fillings and integer values of ν [7–10]. Whenthe carrier density is slightly tuned away from these in-teger filling factors, strongly interacting itinerant chargecarriers emerge. The interaction strength of charge carri-ers is characterized by a dimensionless parameter knownas the Wigner-Seitz radius rs [11, 12], given by:rs =1√πaaB=am∗e24π3/2ε0ε′rℏ2. (1)Here, a = 1/√n is the average distance between chargecarriers, n the charge carrier density, aB the (effective)Bohr radius, ε0 the vacuum permittivity, ε′r the effectiverelative dielectric constant of the material, ℏ the reducedPlanck constant, m∗ the charge carrier effective mass ande the elementary charge. Near an energy gap, the low∗ dolleman@physik.rwth-aachen.de† stampfer@physik.rwth-aachen.decarrier density (and large a) together with a large m∗of itinerant charge carriers results in a large value of rs.This high rs is likely crucial for the formation of stronglycorrelated electronic phases in twisted bilayer graphene(tBLG), such as the superconducting phases adjacent tothe correlated insulating states, the origin of which re-mains not fully understood [8, 13–18]. However, pre-cise measurements near these gaps pose challenges dueto variations in the twist angle across the sample, whichhave been identified as a significant source of disorderin tBLG [6, 9, 14, 19–21]. These variations give rise totwist-angle domains within the sample, characterized byrelatively uniform twist angles but abrupt transitions atthe boundaries [22, 23]. Since the twist angle dictatesthe position of moiré-induced energy gaps, their locationsvary within the sample, complicating efforts to maintaina uniform itinerant charge carrier density and rs acrossthe sample geometry.In this study, we capitalize on the twist angle variationsto isolate single twist angle domains, accurately quan-tify their charging energy, and extract the interactionstrength between itinerant charge carriers. To achievethis, we utilize tBLG heterostructures with varying sizesand geometries, as described in Section II. We demon-strate that the inherent twist angle variations induce elec-trostatic confinement near the insulating states of tBLG(Section III). These confined regions function as chargeislands, and we precisely measure their capacitance us-ing Coulomb oscillations. An in-depth analysis of thiscapacitance reveals a negative electronic compressibility2750 nm500 nm350 nm200 nm(b) (c) Device D1Filling factor ν0.11Temperature T (K)4 52 30-1 1Rxx (kΩ)SCSCBI BICI CI CI CI0.33(d)(e)-3-4-550250-2 -1Gate voltage Vg (V)0 1T = 0.03 KT = 5.05 K2 3 4 5Filling factor ν-3-4-5 -2 -1 0 1 2 3 4 5-2-3-5 -44 52 30-1 1Gate voltage Vg (V)-2-3-5 -4(a)λλAAAA AA AA AAAA AASiliconGraphitegateCr/Au contactsHexagonalboron nitrideAA siteTwisted bilayer grapheneSiO2θ 0.1 1 10 100Rxx (kΩ)5 µmFIG. 1. (a) Schematic cross-section of a tBLG device, including an illustration of the moiré lattice forming between the twographene layers (bottom left) and the resulting stacking order (bottom right). Here, θ denotes the twist angle between thelayers and λ the moiré wavelength. (b) Illustration of the constricted Hall-bar device D1, with the width of the constrictionslabeled. (c) Optical microscopy image of device D1. (d) 4-point resistance as a function of gate voltage Vg and bulk fillingfactor ν at two temperatures (T ), highlighting the fragile superconducting region (SC). (e) 4-point resistance as a function ofgate voltage Vg and temperature of the 750-nm-wide constriction, measured with the contacts highlighted by arrows in panel(c). Band insulators (BI), correlated insulators (CI) and a fragile superconducting dome (SC) are visible.near both the band insulating states in the remote andflat bands, as well as near the correlated insulating states(Section IV). By fitting a model to the data, we find thatthe negative capacitance contribution is proportional tothe square root of the charge carrier density, consistentwith expectations for correlated carriers. Furthermore,as the compressibility of the charge island remains unaf-fected by magnetic fields and is consistent among differ-ent bands, we propose that the observed negative com-pressibility is best explained by the formation of a Wignercrystal (Section V). Thus, our technique provides in-sights into the intriguing properties of itinerant chargecarriers in tBLG when the Fermi energy approaches themoiré-induced energy gaps.II. SAMPLES AND SETUPThis study includes a total of six tBLG devices withdifferent twist angles and geometries. The tBLG is cre-ated using either the ”tear-and-stack” [5] or the ”cut-and-stack” technique [24], generating a moiré pattern withperiodicity λ, as illustrated in Fig. 1(a). The tBLG isencapsulated in hexagonal boron nitride (hBN), whichserves as an atomically flat protective layer with elec-trical insulation [25]. To maintain a uniform electricfield and minimize electrostatic potential disorder, weutilize a graphite gate, promoting atomically flat inter-faces [26, 27]. For low-resistance one-dimensional con-tacts, we employ selective reactive ion etching and met-allization techniques [28–30]. The resulting device struc-ture is depicted in Fig. 1(a), and an example device isshown in Figs. 1(b)–1(c). Further details on the fabrica-tion process and devices are available in Appendix A4.In the main part of this work, we mainly show ex-emplary results from a specific device labeled D1, whichincorporates a Hall-bar structure with constrictions, asshown in Figs. 1(b)–1(c). These constrictions allow usto differentiate between edge and bulk confinement (seeSupplemental Material S2 [31]). Within this device, wehave implemented constrictions of varying widths: 750nm, 500 nm, 350 nm, and 200 nm, denoted as C1, C2,C3, and C4, respectively.To confirm the nature of the 750-nm-wide constriction3Conductance G (µS)-5 -4BI-3 -2 -1 0Bulk filling factor ν1 2 3 4 5050100 CNPBIGate voltage Vg (V)-5 -2 -1 0 1 2 4 52004060 G (µS) G (µS) G (µS)-3.52 -3.51 -3.521-4.0 -3.5 -3.0Vg (V)Vg (V) Vg (V)Vg (V)-4 -33.51 3.52 3.530.811.23.0 3.5 4.0 4.53FIG. 2. Two-point conductance G as a function of gate volt-age Vg over the entire gate voltage range. Zoom-ins close tothe band insulating states reveal high-frequency conductanceoscillations.C1 (Device D1) as tBLG, we analyze the temperature-dependence of the resistance shown in Figs. 1(d)–1(e). Inthis measurement, we can distinctly identify the band in-sulators (BI) positioned around the filling factor ν = ±4.By pinpointing their precise locations in gate voltageand utilizing the gate lever arm, we can approximate thetwist angle as θ ≈ 1.02◦ (see Appendix A4). Moreover,in the vicinity of integer fillings of the partially occu-pied flat bands, we observe resistance peaks correspond-ing to fractional fillings ν = −2, 1, 2, 3. These resistivefeatures align with expectations for correlated insulating(CI) states [8, 32]. Close to ν = −2.7, the resistanceexhibits a significant reduction with decreasing temper-ature, although it remains finite. This behavior indi-cates the presence of a fragile superconducting state [33]in the sample, as evidenced in the Supplemental Mate-rial S1 [31]. Furthermore, magnetotransport measure-ments (Supplemental Material S1 [31]) uncover the exis-tence of Chern insulators featuring identical topologicalinvariants to those reported in other studies of tBLG nearthe magic angle [24, 34–38]. Consequently, we can con-clude that our tBLG sample distinctly exhibits correlatedphases akin to those observed in prior investigations nearthe magic angle [8, 15–18, 32].Gate voltage Vg (V) Vb (mV)0.70.80.91.01.11.2 G (µS)-3.55-3.53-3.54-3.55-3.56-3.545-3.555-0.50.501.01.52.02.5G (µS)(a)(b)FIG. 3. (a) Exemplary 2-point conductance (G) trace closeto the band insulators or constriction C1. (b) Conductanceas a function of gate voltage and bias voltage Vb in a smallgate voltage range.III. CHARGE ISLANDS REVEALED BYCOULOMB OSCILLATIONSEach experiment which now follows, consists of two-point conductance measurements where we increment thegate voltage in small steps to reveal high-frequency con-ductance oscillations. In Fig. 2, we show a conductance(G) trace obtained with 0.1 mV gate voltage steps mea-sured across constriction C1. The charge neutrality point(CNP) at a filling factor of ν = 0, and band insulators(BI) near full filling at ν = ±4 are prominently visiblewithin the trace. A zoom-in of the traces reveals regu-lar high-frequency conductance oscillations at the transi-tion towards the insulating states of our tBLG samples.These oscillations are observed in all our tBLG samplesthat show moiré induced energy gaps (see SupplementalMaterial S2, S4, S5, S7 [31]).A. Dependence on bias voltageTo identify the origin of the conductance oscillations,we perform gate-dependent bias spectroscopy measure-ments, as illustrated in Fig. 3. The combined influenceof bias and gate voltage forms a diamond-shaped regionwhere the conductance is suppressed [Fig. 3(b)], a clearmanifestation of the Coulomb blockade effect [39]. In theCoulomb blockade regime, transport is impeded by elec-trostatic repulsion within a region that confines chargecarriers. If the bias potential eVb is larger than the charg-ing energy Ec = e2/CΣ, where CΣ is the total capacitanceof the charge island, the Coulomb blockade is lifted [40].Due to the similarity in magnitude between the AC lockinexcitation (100 µV root mean square) and the step size inbias and gate voltage, the sharp features of the Coulombdiamonds appear blurred. Nevertheless, we can extract4-3.550 -3.548 -3.546-3.552-3.554Gate voltage Vg (V)02468Magnetic field B (T)02468Magnetic field B (T)-2000200-d2Vxx/dVg (V2/V2)20 2 4 6 8 10Frequency fg (x103 V-1)0FFT of -d2Vxx/dVg  (arb. units)102030(a)(b) 2FIG. 4. Magnetic field dependence of the oscillations. (a)Gradient of the voltage drop (measured in the 4-point con-figuration) as a function of gate voltage and magnetic field.(b) Power spectral density of the oscillations in panel (a) asa function of magnetic field.the bias voltage (Vb) required to lift the blockade, whichis eVb ≈ ±0.31 meV, as indicated by the dashed lines inFig. 3(b).The spacing between two Coulomb resonances on thegate-voltage axis is determined by two energy scales: theelectrostatic charging energy and the quantum level spac-ing. The remarkable regularity of the measured Coulomboscillations (Fig. 2) suggests that the charging energydominates over the quantum level spacing. Therefore,the distance ∆Vg between the Coulomb resonances isgiven by ∆Vg = e/C where C represents the capaci-tance between the charge island and the graphite backgate [40, 41]. From Fig. 3, we extract a periodicity of theoscillations e∆Vg ≈ 0.61 meV. From this analysis, we canconclude that C ≈ CΣ/2, suggesting that half of the to-tal capacitance of the charge island arises from couplingto charge carriers in the source/drain leads.B. Magnetic field dependence of the oscillationsNext, we explore the behavior of the Coulomb oscilla-tions under a perpendicular magnetic field B. We con-duct measurements by sweeping the gate voltage in anarrow window while incrementally increasing the mag-netic field. We observe a significant shift in Coulombresonance positions [Fig. 4(a)], while the spacing be-tween them appears unaffected by the magnetic field. Afast Fourier transform (FFT) at each magnetic field re-veals two distinct primary frequency components, withhigher harmonics arising from the non-sinusoidal shapeof the Coulomb resonance [Fig. 4(b)]. The power spec-tra’s evolution with magnetic field displays a pronouncedchange in peak amplitude: the frequency component near1600 V−1 decreases, and the component near 1100 V−1becomes more prominent. As neither component fullydisappears (which is evidenced by the phase coherenceof both frequency components in the full magnetic fieldrange, see Appendix A1), they must originate from twodistinct confined regions within the sample. The constantspacing further confirms that the charging energy dom-inates over the quantum level spacing in the Coulombblockade effect. If the quantum level spacing were moresignificant, the single-particle energies would evolve dif-ferently due to the lifting of spin- and valley degeneracywith the magnetic field [42, 43].The position of the oscillations in Fig. 4(a) exhibitsperiodicity in 1/B. This behavior, recently been demon-strated in bilayer graphene quantum dots [44], can beexplained by a classical electrostatic shift induced by den-sity of states oscillations near the confined region. Whenthe Fermi energy matches a Landau level, the density ofstates in the surrounding area increases, exerting a clas-sical electrostatic force that shifts the energy levels insidethe confined regions. Since an extremum in Fig. 4(a) sig-nifies a constant electron density inside the confined re-gion, the change in charge density must occur outside theconfined region. Therefore, these measurements providevaluable information about the Fermi surface of the car-riers in the source/drain leads that couple to the chargeisland [44], further discussed in Appendix A1.C. Origin of the Coulomb oscillationsThe presence of Coulomb oscillations requires an en-ergy gap for charge confinement, such as a band gap,and spatial variations in the position of this gap rela-tive to the Fermi level to create a confining potentialfor charge carriers. In Supplemental Material S2 [31],we investigate the Coulomb oscillations as a function ofconstriction width on device D1, finding a clear trend ofCoulomb oscillations vanishing with decreasing samplewidth. This suggests that their origin lies in bulk char-acteristics rather than edge effects. Furthermore, the en-capsulation of tBLG with hBN and the use of graphitegating should strongly suppress the effects of charge im-purities and other potential disorder [27]. Therefore, wepropose that the confinements arise due to variations inthe twist angle across the sample geometry [6, 9, 14, 19–21], leading to the formation of twist angle domains overthe sample [22, 23], as illustrated in Fig. 5(a). Local5EnergyFermi levelGapHolesElectronsDistanceTwist angleθ (arb. units)DistanceDistanceEnergyEnergyConfined region(a)(b)(c)FIG. 5. Schematic illustration of twist angle variations andresulting bandstructure variations over the sample geometry.(a) Illustration of the formation of twist angle domains overthe sample geometry, the local variation of the bandstructureis drawn on the right hand side. (b) The Fermi level increases,and in this example the areas with a small twist angle willbecome insulating since the Fermi level is in their band gap,while other areas remain conducting. (c) Further increasingthe Fermi level results in the formation of confinements, wherethe energy levels are quantized.variations in the twist angle result in local changes in theenergy gap offset relative to the Fermi level [see the po-tential landscape on the right-hand side of Fig. 5]. There-fore, as the Fermi level approaches the band gap, certaindomains become insulating earlier than others [Fig. 5(b)].Close to the insulating state, individual domains may re-main conducting while the surrounding area has alreadybecome insulating [Fig. 5(c)]. This scenario leads to indi-vidual charge islands and the observation of conductanceoscillations as a function of gate voltage.However, this leaves the question of how these chargeislands are accessible in the transport experiment. Anin-depth analysis of the quantum oscillation frequencyin Appendix A1 reveals that the carriers tunneling inand out of the confinement exhibit a complete absence ofmoiré-induced energy gaps. This observation points to-ward an important role of the boundaries between twistangle domains. Within these boundaries, strong disordermay be present on the length scale of the moiré superlat-tice. Consequently, locally, the moiré-induced gap mayvanish, leaving charge carriers which tunnel in and outof the confinement. Since these boundaries are expectedto appear over the entire sample geometry, we proposethat they form a network that allows to probe the localcharge confinement in electronic transport experiments.IV. NEGATIVE ELECTRONICCOMPRESSIBILITYNext, we study the spacing of the Coulomb oscilla-tions ∆Vg to determine the capacitance of the chargeislands. Since ∆Vg = e/C, the oscillation frequencyfg can be expressed as fg = C/e, directly proportionalto the back gate capacitance of the charge island. Tofacilitate the analysis, we convert conductance tracesto the frequency domain by calculating power spectraPω = |F{dG/dVg}|2, where F represents the Fouriertransform (for more details, see Appendix A4). The re-sulting gate-voltage-dependent power spectra, presentedin Fig. 6(a) for device C1, reveal multiple distinct fre-quency components near the insulating states (labeled asC1.1 to C1.4), each corresponding to a single charge is-land. While the oscillation period remains regular withinsmall gate voltage ranges (< 10 mV), a continuous fre-quency tuning is evident on larger voltage scales. Acrossall tBLG devices in this study, we consistently observe anincrease in frequency (and capacitance) as the Fermi levelapproaches the moiré-induced energy gap, representing acentral result of our work (additional examples are pre-sented in Supplemental Material S2, S4, S5, S7 [31]).To elucidate the observed change in capacitance, wefirst consider the geometric capacitance Cg of the chargeisland within a parallel-plate capacitor model, given byCg = ε0εrA/d, where A is the area of the island, and dis the distance between the tBLG and the graphite gate.Based on the twist angle domain model outlined in Sec-tion III C and Fig. 5, we anticipate that A and Cg shouldremain constant as a function of gate voltage if the chargeisland consists of a single domain. However, if the chargeisland comprises multiple domains, we expect to observea step-wise reduction of A and the capacitance. Alterna-tively, if the confinement arises from other types of po-tential disorder [45], we anticipate a continuous reductionof the area A, as more of the surrounding area should be-come insulating when the Fermi level enters the bandgapin the surrounding bulk. Given that these scenarios areinconsistent with our observations, we can conclude thatthe observed increase in capacitance is indicative of anegative compressibility of the charge carriers.A. Models for the negative capacitancecontributionTo investigate the negative compressibility contribut-ing to the observed increase in capacitance, we considertwo scenarios where negative compressibility can arise.First, we consider the exchange interaction in an elec-tron gas, where the reduced likelihood of finding chargecarriers with the same spin at the same position cre-ates an ”exchange hole” due to the opposite background6Gate voltage Vg (V)Conductance G (µS)-5-5-4-4-3-3-2-2-1-100Bulk filling factor ν1122334455050100(a)010002000300040005000Frequency f g (V-1)0 0.5Pω (S2/V2)100.10.20.30.40.50.60.7Capacitance C (fF)-4.0 -3.5 -3.0Gate voltage Vg (V)C1.1V0,C1.1 V0,C1.2C1.20 0.5 1Pω (S2/V2)(b)zoom in (inset)x3Strenght of interaction S (c)00 1 2 31234CountsIncluded in fitExcluded from fitGaussian fitµ (S) = 0.84σ (S) = 0.26C1.1C1.2C1.3C1.4FIG. 6. Observation of negative compressibility and determining the interaction strength S. (a) Power spectrum Pω versuscapacitance C and gate voltage Vg on constriction C1. The conductance trace (right vertical axis) is included for easy com-parison. (b) Zoom-in of Pω near the band insulators at ν = −4. Equation (5) is fit to two frequency components. The inset isscaled vertically by a factor 5 to highlight component C1.2. (c) Histogram of all the values of S determined from fitting anddetermining the mean of S.charge [12, 46, 47]. The resulting negative interactionenergy Ei,X is given by:Ei,X = − e2An3/23√2π3/2ε0ε′r[(1 + ξ)3/2 + (1− ξ)3/2], (2)where ξ represents the polarization of magnetic moments,ranging from 0 (unpolarized) to 1 (fully polarized). Sec-ondly, we consider the case of a Wigner crystal, wherecharge carriers minimize their potential energy by form-ing a solid phase with a triangular lattice [48]. The result-ing negative interaction energy Ei,W is given by [48–53]:Ei,W = −ηT e2An3/28πε0ε′r, (3)where ηT = 3.92 is a numerical constant associated withthe triangular lattice.In both scenarios, the negative interaction energy Eileads to a negative (thermodynamic) electronic compress-ibility κ, expressed as κ−1 = dµ/dn = (1/A)(d2Ei/dn2),where µ is the electrochemical potential. The nega-tive compressibility contributes to a negative capacitancecontribution Ci, where C−1i = κ−1/(Ae2), increasing thetotal capacitance C beyond the geometric contributionCg according to the relation [52, 54–56]:C−1 = C−1g + C−1i = C−1g +1A2e2d2Eidn2. (4)Since the carrier density dependence of both interactionenergies is the same [Eqs. (2)–(3)], we can use Eq. 4 toobtain a general model for the capacitance, including theeffect of correlated itinerant charge carriers:1C=dε0εrA− Sε0ε′rA√|n|, (5)where the dimensionless parameter S characterizes thestrength of the correlation, which can be experimentallydetermined through fitting. The expressions and valuesfor S for the different theoretical scenarios can be foundin Table I. To fit Eq. 5 to our data, we express the carrierdensity as |n| = α|Vg − V0|, where Vg is the gate voltage,V0 is the gate voltage where the itinerant charge carrierdensity is zero, and α is the lever arm. We determine α7TABLE I. Expressions for the S in different correlated elec-tron models compared to the experimental value. The abbre-viation WC denotes the Wigner crystal scenario, while HFdenotes the Hartree-Fock model in the exchange scenario.S ValueWC3ηT32π0.12HF, unpolarized12√2π3/20.064HF, polarized1.4142√2π3/20.090Experiment - 0.84± 0.26from the Landau levels that emerge in magnetotransport(see Supplemental Material S1 [31]).B. Determining the strength of the interactionFigure 6(b) shows the power spectrum of constrictionC1 close to the band gap at ν = −4. We identify and de-termine the maxima of two frequency components and fitEq. (5) to extract S. The resulting excellent fit demon-strates that the n−1/2 dependence in our model pro-vides a good description of the gate-voltage-dependentfrequency observed in the experiments. We perform thisanalysis on 12 frequency components from five differentsamples that are close to the magic angle: constrictionsC1, C2 and C3 on Device D1; and Devices D2 and D3.The obtained values of S are summarized in Fig. 6(c) andthe complete fitting results are presented in the Supple-mental Material S3 [31]. The mean value of S is found tobe µ(S) = 0.84 with a standard deviation of σ(S) = 0.26.Three outliers in Fig. 6(c) were excluded from calculat-ing the mean value due to relatively large uncertaintybounds resulting from the fitting procedure.C. Analysis with a fixed strength of the interactionNext, we refine our model by fixing the parameter Sto the mean value of S = 0.84, i.e., we account for allour data from devices close to the magic angle (the con-strictions C1, C2, C3 and devices D2, D3) and band gaptransitions using one single, fixed value of S. This al-lows us to fit the subset of frequency components thatoccur within a narrow gate voltage range, expanding ourdataset to include 18 frequency components. Note thatin this refined model, both the vertical offset and steep-ness of the curve are solely determined by the size of theconfinement A, while V0 determines the horizontal posi-tion. The fitting results are presented in Figure 7, wherethe vertical axis represents the inverse capacitance 1/C,and the horizontal axis represents the average distancebetween the charge carriers a = α√|V − V0|. Each fre-quency component is labeled using the format ”xx.y”,where ”xx” is the constriction or device number, and”y” is the number of the frequency component. Accord-ing to Eq. (5), the resulting maxima should fall along thestraight line determined by the fit. Remarkably, using afixed value of S provides excellent fits to all the data,regardless of the filling factor or the band where the fre-quency component is observed. The remaining plots ofthe frequency components, including the fits and the fit-ted values of A and V0, are presented in the SupplementalMaterial S3 [31].It is noteworthy that we also find periodic Coulomboscillations that are well-described by Eq. (5) close tothe insulating states that form at fractional superlatticefillings [Fig. 7(c)]. These include components D2.4 foundnear ν = 2, D3.2 near ν = −2 and D3.3 near ν = 3. Thelatter findings clearly demonstrate that a single-particleband gap is not necessary for confined regions to arise inthe sample; instead, the charge gap arising from the cor-related insulating states is sufficient. The fact that theyare well-described by Eq. (5), assuming a zero itinerantcharge carrier density at V0 in the partially filled band,is also in-line with the opening of a correlation-inducedenergy gap near the integer fillings of the superlattice.D. Size and twist angle of the charge islandsWithin this subsection, we delve deeper into the re-sults of our modeling and fitting analysis to evaluate theiralignment with the scenario we proposed in Section IIIinvolving twist angle domains. By comparing V0 to thevoltage where we pinpoint the center of the band gap inthe bulk, we can estimate the deviation of the twist angle.We do this for constriction C1, and find that componentsC1.1, C1.2, C1.3, C1.4 deviate −0.023◦, 0.016◦, 0.011◦and −0.003◦ respectively from the bulk value of 1.022◦.These twist angle variations are in-line with observationsin the literature [23]. This shows that minute twist an-gle variations of only a few percent from the average bulkvalue suffice to produce the confinements observed in thiswork.Figure 8 displays the empirical cumulative distribu-tion function (CDF) of the sizes (√A) of the 18 chargeislands extracted from our fits on all frequency compo-nents. These sizes are consistent with the dimensions ofour samples and similar to the experimental twist anglemaps obtained by Uri et al. [23]. A notable differenceis that we find no substantially smaller (√A < 125 nm)or larger (√A > 425 nm) areas in our experiments. Wenote that small areas are usually embedded into a largertwist angle domain that exhibits more uniformity, caus-ing twist angle boundaries that can couple to the chargeisland to be absent [23]. On the other hand, large do-mains are less likely to form a detectable charge island,since the large circumference makes it less probable thatthe entire surrounding region is insulating at the same8Average distance a (nm)Average distance a (nm) Average distance a (nm)Average distance a (nm)Inverse capacitance 1/C (  1016 F-1) 1/C (  1016 F-1)1/C (  1016 F-1)1/C (  1016 F-1)Near ν = ±4, flat bands Near ν = ±4, remote bandsFractional fillings15 20 25 30 35 40 45 50 55 6000.511.522.533.544.55D2.1C1.2C1.1C3.2C2.3C1.3C1.4C2.5C2.1C3.1C2.2C2.4D2.3D2.2D3.1D3.3D2.4D3.220 25 30 350.30.40.50.6inset20 30 40 5000.511.522.510 20 30 40 5000.20.40.60.81(a) (b)(c)Ek kkEEFIG. 7. Maxima in the power spectrum used to fit Eq. (5), plot with 1/C on the vertical axes and the average distance betweencarriers 1/√n = a on the horizontal axes such that Eq. (5) becomes a straight line. The fits to each frequency componentare made with a fixed value of S = 0.84. The frequency components are divided into panels such that: (a) shows componentsfound near ν = ±4, on the side of the flat band (see the simplified bandstructure in the inset) (b) shows components foundnear ν = ±4, on the side of the remote band and (c) shows frequency components which are found near fractional fillings ofthe superlattice.0.10 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1Size of confinement (µm)00.20.40.60.81Cumulative probabilityOur workUri et al. [23]leadsFIG. 8. Empirical cumulative distribution function (CDF) ofthe size of the confinements√A compared to the sizes of twist-angle domains estimated from the results of Uri et al. [23].Insert shows a schematic of the confinements and the ”leads”along the twist angle boundaries where transport is thoughtto occur.time. Considering these factors, our size distributionaligns reasonably well with the literature, validating ourmodel and the obtained value of S ≈ 0.84.V. NATURE OF THE CORRELATED STATEOur work provides compelling evidence for a neg-ative capacitance contribution within confined regionsof tBLG. However, the nature of the correlated stateremains uncertain. The analysis has established the1/√|n| dependence of the inverse capacitance, consis-tent with both an exchange contribution or a Wignercrystal phase. However, the value of S ≈ 0.84 is signif-icantly larger than the expectation for both the Wignerand exchange scenarios (Table I). Nevertheless, severalexperimental observations may help clarify the nature ofthe correlated state:(1) The observation of regular Coulomb oscillationssuggests that the spin- or valley degree of freedom is notsignificant for the addition energy in the confined region.(2) Substantial changes in the polarization of magnetic9moments with magnetic field, including the disappear-ance of correlated insulating phases at high magneticfields (see Supplemental Material S1 [31]), indicate thatthe Zeeman energy exceeds the exchange-induced energygap. At a magnetic field of 9 T, the Zeeman energyin the order of 1 meV should dominate over other spincontributions, leading to a fully polarized electron gas.Consequently, if the electron gas is not polarized at zerofield, the change in polarization between 0 T and 9 Twould lead to a frequency increase. For example, wewould expect a frequency increase of 72% as ξ → 1 forthe 1100 V−1 component in Fig. 4(b). Instead, an indica-tion of a change in the Coulomb oscillation frequency isobserved in none of our measurements, additional exam-ples are included in the Supplemental Material S9 [31].(3) If the charge carriers are fully out-of-plane polar-ized at B = 0 T, no change in the Coulomb oscillationfrequency is anticipated in the exchange scenario. How-ever, no manifestation of such ferromagnetism such as ananomalous Hall effect [57] is observed, nor do we observehysteresis when sweeping the magnetic field from -9 to9 T and vice versa (see Supplemental Material S8 [31]).(4) The carrier density dependence of the polarizationof magnetic moments in the fluid phase is expected tobe strong [14], potentially leading to significant devia-tions from the observed 1/√n-dependence of the inversecapacitance in the exchange scenario.(5) Describing frequency components with the samevalue of S regardless of whether they are found near acorrelated or band insulator is unexpected in the case ofexchange contributions, according to theoretical simula-tions [58], because the spin/valley polarization varies foreach energy gap that opens in tBLG [14].(6) Theoretical models that take short-range interac-tions into account indicate that strong negative contribu-tions to the capacitance are not anticipated close to theband insulating states at ν = ±4 [58].The absence of clear signatures expected in the ex-change scenario and the consistent value of S across dif-ferent bands suggests that a Wigner crystal phase mayprovide a better explanation for the observed negativecompressibility. The exchange energy in a Wigner crys-tal follows a scaling law of the form EX ∝ exp−γ√rs [59,60], where γ is of the order 1. As a consequence, the in-teraction energy of a Wigner crystal does not rely onthe polarization ξ, providing an explanation for the sixobservations listed above. The higher value of S than ex-pected for a Wigner crystal may be related to the averagedistance a being in the same order to the moiré wave-length λ, causing a fraction of charge carriers to obtaineven lower energy states within the non-uniform poten-tial landscape provided by the moiré lattice. Importantly,as shown in the Appendix A2, this effect also leads tothe same 1/C ∝ 1/√n dependence, consistent with theexperimental findings. Further investigations, consider-ing additional factors and refinements to the model, mayprovide a more detailed understanding of the presentlyobserved correlated state in tBLG.VI. CONCLUSIONSThis study demonstrates a negative electronic com-pressibility in confined regions of tBLG. The observeddependence of the inverse capacitance on the charge car-rier density, characterized by 1/C ∝ 1/√n, aligns withthe presence of strongly correlated itinerant charge car-riers. The magnetic field dependence and the similarityof the correlation strength parameter S across differentbands suggest a limited role of exchange contributions,pointing towards Wigner crystallization as a likely expla-nation for the observed negative compressibility. Thesefindings highlight the role of naturally occurring electro-static confinements in tBLG, which enable precise inves-tigation of the ground state energy of correlated statesin this material.Our research reveals two distinct phenomena rootedin the moiré physics of tBLG. The first phenomenon in-volves electrostatic confinements arising from variationsin the twist angle. This effect hinges solely on the pres-ence of moiré-induced energy gaps, while the flatness ofthe band is not important. Therefore, this phenomenonis not exclusive to samples near the magic angle, as ev-idenced in tBLG samples far from this angle (Supple-mental Material S4, S6 [31]). The second phenomenonis the manifestation of negative compressibility withinthese confinements. In this case, the band flattening intBLG plays a crucial role, since rs is required to be large.This gains support from control experiments, particularlyone conducted on a sample where the twisted bilayer hasrelaxed. In this case, we observe the absence of moiréminigaps, while the charging spectrum of a confinementin this sample shows no negative compressibility (Sup-plemental Material S6 [31]). Furthermore, we extend ourinvestigation to samples with twist angles significantlydeviating from the magic angle (θ ∼ 0.65◦ in Supple-mental Material S5 and θ > 1.29◦ in Supplemental Ma-terial S7 [31]). While these samples do exhibit an increasein frequency as the Fermi level approaches an energy gap,their evolution with gate voltage deviates from the trendspecified in Eq. (5). In these instances, the parameterrs may not be sufficiently large to induce negative com-pressibility in the full gate voltage range where the con-finement is formed, a factor that can be attributed to thereduced band flattening compared to the magic angle.While it is commonly believed that correlated effectsare not significant in the remote bands due to the ab-sence of correlated insulating states [61], our observa-tion of negative compressibility in these bands [Fig. 7(b)]challenges this notion. An analysis of the band struc-ture, including a Hartree-Fock correction, suggests thatthe observed negative compressibility may be attributedto a significant correlation-induced flattening of the bandwhen the filling factor exceeds |ν| = 4 (see Appendix A3).The size of a typical charge island, as illustrated inFig. 8, is notably larger than the moiré wavelength (λ =13.8 nm for a 1.02◦ twist angle). This suggests that theperiodicity of the moiré potential remains a crucial factor10within the charge island, and the negative compressibilityis not solely a consequence of electrostatic confinement.Our experimental evidence strongly supports this idea,revealing a significant suppression of the kinetic energiesof carriers within the confinement compared to single-layer or Bernal-stacked bilayer graphene. In essence, wehave demonstrated how moiré-induced confinement canoffer quantitative insights into the interaction energies ofcorrelated effects in tBLG by accurately measuring thelocal charging energy. This is the principal achievementof our work and has significant implications for under-standing and manipulating correlated phases in moirématerials.The source data and Matlab code underlying thispaper are available at Ref. [62] upon publication of thismanuscript.ACKNOWLEDGMENTSThe authors thank C. Volk and J. Sonntag for exper-imental support; F. Volmer, T. Ouaj and M. Schmitzfor assistance during sample fabrication; S. Trellenkampand F. Lentz for electron-beam lithography, B. Shklovskiifor correcting Eq. (2); and L. Gaudreau and M. Mor-genstern for discussions. This work was supported bythe FLAG-ERA grants 436607160 2D-NEMS, 437214324TATTOOS and 471733165 PhotoTBG by the DeutscheForschungsgemeinschaft (DFG, German Research Foun-dation) and by the Deutsche Forschungsgemeinschaft(DFG, German Research Foundation) under Germany’sExcellence Strategy – Cluster of Excellence Matter andLight for Quantum Computing (ML4Q) EXC 2004/1– 390534769, by the Deutsche Forschungsgemeinschaft(DFG, German Research Foundation) within the PriorityProgram SPP 2244 (535377524) and from the EuropeanResearch Council (ERC) under the European Union’sHorizon 2020 research and innovation programme (grantagreement No. 820254). A.F, L.K and D.M.K acknowl-edge funding by the Deutsche Forschungsgemeinschaft(DFG, German Research Foundation) under RTG 1995and within the Priority Program SPP 2244 “2DMP”.F.L. acknowledges support by the Austrian Science Fund(FWF), project I-3827-N36 and by the COST associa-tion, COST action CA18234. This work was supportedby the Max Planck-New York City Center for Nonequilib-rium Quantum Phenomena. Fabrication of the sampleswas supported by the Helmholtz Nano Facility [63]. K.W.and T.T. acknowledge support from JSPS KAKENHI(Grant Numbers 19H05790, 20H00354 and 21H05233).R.J.D. and A.R. performed the experiments. A.R. fab-ricated the samples, L.W. built the laser cutting setup;and K.W. and T.T. grew the hBN crystals. R.J.D. per-formed the data analysis and the results were interpretedwith the input of A.R., A.F., L.K., L.W., D.M.K., F.L.,B.B. and C.S. R.J.D wrote the manuscript with the inputof A.R., A.F., L.K., L.W., D.M.K., F.L., B.B. and C.S.Appendix A1: Fermi surface area extracted from thequantum oscillationsThe electrostatically-induced shifts in the position ofthe Coulomb resonances provides important clues aboutthe transport characteristics between the metal leads andthe confined regions. We plot the phase of the Coulomboscillations (which is proportional to the Coulomb res-onance position) obtained from the FFT analysis ofthe two prominent frequency components [Fig. 4(b)] inFig. A1(a). This reveals a quantum oscillation that isperiodic in 1/B, and the quantum oscillation frequencycan be extracted from an additional FFT analysis of thephase [inset of Fig. A1(a)]. By monitoring the phase ofthe oscillations at different gate voltages, we can extractthe quantum oscillation frequency and plot these as redsquares in Figs. A1(b)–A1(c). Since these quantum oscil-lations have the same microscopic origin as the Subnikov-de Haas (SdH) oscillations observed in the sample’s bulkmagnetoresistance (see Supplemental Material S1 [31]),we also add the SdH frequency as a comparison [bluedots in Fig. A1(b)]. The quantum oscillation frequenciesobtained from the Coulomb oscillations form a straightline that intersects the charge neutrality point at zeromagnetic field in Fig. A1(c). In contrast, the SdH fre-quency intersects with fractional superlattice fillings atzero magnetic field due to the Dirac revival effect intBLG, which reconstructs the Fermi surface area at frac-tional fillings [14]. We observe the same disparity in asecond device, labeled D4, which is presented in the Sup-plemental Material S7 [31]. At higher gate voltages, wenote that the lever arm of the gate decreases, which ex-plains why the quantum oscillation frequency from theCoulomb oscillations exhibits a more gradual slope [reddashed lines in Figs. A1(b)–A1(c)] compared to the SdHfrequency component emerging from the charge neutral-ity point [blue lines in Fig. A1(b)].The Onsager relation establishes a direct proportion-ality between the frequency of quantum oscillations andthe extremal Fermi surface area in momentum space [64].In proximity to an energy gap, it logically follows thatboth the Fermi surface area should approach 0 nm−2and the frequency of quantum oscillations approaches0 T. Therefore, the fact that the quantum oscillation fre-quency forms a linear relationship emanating from thecharge neutrality point (CNP) is evidence for the absenceof moiré-induced energy gaps for the charge carriers thatare involved in the tunneling processes in- and out of thecharge islands.110Magnetic field B (T)-500050010001500200025003000Phase (o ) 1100 V-11600 V-12 4 6 8Frequency (T)0FFT of phase(arb. units)0 5017.6 T(a)051015202530-3-4-5 -2 -1 0 1 2 3 4 5SdH oscillationsCoulomb Gate voltage Vg (V) Gate voltage Vg (V)Frequency (T)Frequency (T)(b) (c)15202543 3.5-3.5 -3-4FIG. A1. (a) Phase of 2 prominent frequency components extracted from the fast Fourier transform in Fig. 4(b). (b) Quantumoscillation frequency of the magnetoresistance oscillations (see Supplemental Material S1 [31]) compared to the quantumoscillations in the phase of the Coulomb oscillations, as a function of gate voltage. (c) Zoom-in of the quantum oscillationfrequency in the case of the Coulomb oscillation.Appendix A2: Ground state energy of a Wignercrystal in a moiré superlatticeIn this Appendix we show that a periodic moiré poten-tial may lead to a lowering of the ground state energy ofa Wigner crystal. Due to the periodic moiré potential, afraction of charge carriers in the triangular Wigner lat-tice can obtain an even lower energy within the potential,lowering the ground state energy even further than com-pared to a uniform background. If we assume a potentialwith minima ∆E, we can write this additional interactionenergy contribution as:Emoiré = −An∆ER, (A1)where R represents the effective fraction of charge car-riers that are in the optimal position to profit fromthe moiré potential. For simplicity, we assume thatR = amoiré/a (where amoiré is the size of the moiré su-percell) such that R = 1 if amoiré = a. Since 1/a =√n,the moiré energy contribution becomes:Emoiré = −An3/2∆Eamoiré, (A2)this leads to an additional capacitance contribution:1Cmoiré= −3∆Eamoiré4e2A√n. (A3)with the same dependency on area and density as theinteraction energy for a uniform background.To estimate how much the ground state energy is low-ered, we compare to a tight-binding model of the moirésuperstructure [65]. For small twist angles, the local ro-tation between unit cells in the top and bottom layers canbe neglected in favor of only considering a rigid displace-ment vector d(r) between the unrotated top and bottomlayer. In the small angle approximation, d(r) can thenbe explicitly written as d(r) ≈ −θẑ × r. All displace-ments d(r) lie within the unit cell of the pristine lattice.Consequently, one can map each local configuration at apoint r of the moiré supercell to a rigid displacement d(r)in so-called configuration space mapped on the unit cellof the pristine lattice. Following the model we have out-lined in [65], we map out tight-binding parametrizationsin configuration space using a 10 × 10 grid in configu-ration space. We use a continuum elasticity model firstsuggested by Nam et al. [66] to calculate the effects oflattice reconstruction in tBLG. From our parametriza-tion we extract the variations in on-site potential at atwist angle of ≈ 1 degrees, yielding ∆E ≈ 14 meV andamoiré = 14 nm. From this analysis, we find that the newcorrelation strength parameter S = 0.16, or an increaseof 33% compared to the Wigner crystal with a uniformcharge background. This result gives an idea of the orderof magnitude of the correction, but is likely an underes-timation of the correction for two reasons. First, theinteraction with the positive background is not included,and will also result in an additional contribution due tothe non-uniform charges of the underlying lattice. Sec-ond, we do not include effects of elasticity in the Wignercrystal, the Wigner lattice will likely deform to ensure ahigher fraction is in the optimal position.Appendix A3: Suppression of kinetic energy inremote bandsThe atomic and electronic structure of tBLG is cap-tured within an atomistic modeling approach [3, 67, 68]that relies on commensurate moiré unit cells with twistanglecos θ =n2 + 4nm+m22 (n2 + nm+m2), (A1)where m,n ∈ N. For the simulations, we use (n,m) =(32, 33) corresponding to a twist angle of θ = 1.018◦ with12Energy E (eV)0.750.800.850.90Energy E (eV)0.700.750.800.850.90ν = 0ν = 4.0ν = 4.1 ν = 0ν = -4.0ν = -4.1 Κ Γ Μ Κ‘(a)(b)EF (ν = 4.1)EF (ν = -4.1)FIG. A2. Band structure of 1.018◦ tBLG in the presenceof long-ranged electron-electron interactions. Hartree cor-rections render the bandstructure of tBLG filling-dependent,which leads to a pinning of the van Hove singularitiy to theFermi energy EF and additional band flattening at the tip ofthe remote valence and conduction bands.N = 12′676 carbon atoms per moiré unit cell. The posi-tions of the carbon atoms are relaxed using classical forcefields as outlined in Ref. [69]. The electronic structure ismodeled by a Slater-Koster tight-binding model of thecarbon pz-orbitals using the parametrization adopted inRef. [67]. Near the magic-angle, long-ranged Coulombinteractions were shown to significantly renormalize thesingle-particle flat bands of tBLG [70–73] if the systemis filled with electrons (holes), which can be capturedwithin a self-consistent Hartree theory. Ref. [71] demon-strated that within atomistic modeling approaches, theHartree potential can effectively be parameterized by anon-site term of the formV H(r) = V0ν∑jcos(Gj · r), (A2)where ν denotes the electronic filling ν = −4 . . . 4 of theflat bands with respect to charge neutrality (ν = 0) andGj are the three non-equivalent moiré reciprocal latticevectors that differ by rotations around 120◦. The valueof the Hartree potential V0 was found to be V0 = 5 meVfor the unscreened Coulomb interaction [71, 74].The (filling-dependent) bandstructure of 1.018◦ tBLGalong the high symmetry path K −Γ−M −K ′ is shownin Fig. A2. At half-filling, the systems features a setof flat bands that are well separated from the remotevalence and conduction bands by a well-defined energygap. Filling the flat bands of tBLG with electrons (up-per panel) or holes (lower panel), shifts the energies atthe K,K ′ points to higher (lower) energies due to theHartree potential. Therefore, the flattest sections of thebands follow the Fermi energy EF , which leads to a pin-ning of the van Hove singularities [73]. Furthermore, theHartree potential affects the flatness of the tip of thevalence (conduction) band manifold as indicated by thedashed line at filling factor ν = ±4.1. This effect re-duces the kinetic energy of charge carriers at the edge ofthe remote band, which may account for the prominentobservation of a negative compressibility in this regime.Appendix A4: Methods1. SamplesExfoliation: The flakes used for fabrication were me-chanically exfoliated onto a silicon wafer with 90 nm thickthermally grown silicon dioxide [75]. Graphite flakes(”graphenium”) where obtained from NGS NaturgraphitGmbH.Stacking: Device D1 was fabricated using the stack-and-tear method [5]. For Device D1, we used polyvinylalcohol (PVA) and polydimethylsiloxane (PDMS), andthe stacking of the flakes was performed using the pa-rameters described in Ref. [21]. For samples D2, D3, D5and D6 we used a polybisphenol A carbonate (PC) stampon top of a PDMS stamp [28]. In addition, the single-layer graphene flakes of these devices were pre-cut usinga laser. Device D4 was produced using a poly bisphenola carbonate (PC) stamp on PDMS [28].For each sample, the thicknesses of the flakes are mea-sured in tapping mode atomic force microscopy and theresults are shown in Table A1.TABLE A1. Thicknesses of the flakes used in each sample.Each thickness was determined by measuring the step-heightat the edge of the flake after fabrication in an atomic forcemicroscope. The uncertainty on each thickness is 1 nm.Sample Bottom hBN (nm) Top hBN (nm) Graphite (nm)D1 29 32 8D2 32 24 5D3 32 26 4D4 35 37 10D5 31 22 2D6 22 22 3Laser cutting: Laser cutting was performed using afocused supercontinuum laser, coupled into an optical mi-croscope setup. The output of the laser was spectrally fil-tered to contain wavelengths of 400-550 nm, which max-imises the ratio of absorption in graphene compared tothe absorption in Si for the used oxide thickness of 90 nm.The pulse duration on the sample is estimated to be few1310s of picoseconds. The laser has a maximum repeti-tion rate of 20 kHz, but is typically operated at 4 kHz.The power can be varied using absorption filters, witha higher power increasing the width of the cut. Typi-cal pulse energies used for laser cutting are 4 nJ, focuseddown to a submicrometer spotsize, equivalent to a max-imum intensity of approximately 2× 1010 W/cm2.Fabrication: For fabrication, we used a 50K/950Kpolymethyl methacrylate (PMMA) double-layer as ourresist system (Allresist 631.09 and 679.04 both spin-coated at 4000 rpm and baked at 150◦ for 3 minutesper layer, with a mixture of 3 parts isopropylalcoholand 1 part water as a developer) and e-beam lithogra-phy (Vistec EBPG5200+, 100 keV, clearance dose of 500µC/cm2). The Hall-bar device was fabricated by firstpatterning holes into the top hBN layer. These holeswere etched in a reactive ion etcher by first using a shortoxygen plasma step (Oxford PL 100 / ICP at 20 W RFpower, 40 sccm, as low pressure as possible ∼ 8 µbarfor ∼5 s) to remove contamination, followed by a CF4plasma etch (10 W, 40 sccm, low pressure), which sig-nificantly slows down at the graphene layer [29]. Theduration of the CF4 plasma step was adjusted to the tophBN thickness to prevent over-etching into the bottomhBN layer. This was followed by another brief oxygenplasma step to remove the graphene, after which chromeand gold were deposited using e-beam evaparation, fol-lowed by a lift-off in warm acetone (without sonication).This results in a clean one-dimensional edge contact tothe tBLG [28]. Subsequently, the Cr/Au metal contactsand bond pads were patterned and evaporated using thesame lift-off process.For device D4, the fabrication was performed usinga two-step process, where first the device geometry wasstructured using the CF4 plasma, and the electrodes weredeposited in the second step. The contact resistances inthis process are significantly higher than in the other de-vices, and this fabrication approach was was not pursuedfurther.2. Twist angle determinationDevice D1 To determine the twist angle in the 750 nmconstriction (C1), we extract the superlattice filling ns =2.43 × 1012 ± 0.04 × 1012 cm−2 from the Landau fan(Supplemental Material S1 [31]), and use the equationns = 8θ2/√3a2l [8] (where al = 0.246 nm is the lat-tice constant of graphene) to find the twist angle inθ = 1.022◦. For the remaining constrictions, we usedthe position of the band insulating features to determinethe twist angle and find twist angles of 1.07◦, 0.97◦ and0.91◦ for constrictions C2, C3 and C4, respectively.Device D2 consists of a 1 µm-wide Hall bar geometrywith a twist angle estimated from the Landau fan to be0.97◦. Additionally, we found an additional alignmentbetween one of the graphene layers and the hBN, with atwist angle of ∼ 0.7◦. We found no effects from this addi-tional alignment on the charging spectra studied in thiswork. An optical image of this device, the conductancetraces and power spectra of this device are presented inthe Supplemental Material S4 [31].Device D3 is a similar 1 µm-wide Hall bar geome-try without constrictions. The device broke down duringthe Landau fan measurement, therefore we could only ex-tract the lever arm of the back gate, but not the positionof the additional Landau fan. From the position of theinsulating states, we are nevertheless able to estimate thetwist angle to be 1.14◦. An optical image of this device,the conductance traces and power spectra of this deviceare also presented in the Supplemental Material S4 [31].Device D4 due to high contact resistances, it was notpossible to extract the lever arm or the twist angle frommagnetotransport experiments. Therefore, the twist an-gle was estimated using a parallel-plate capacitor modeland the position of the insulating states. We estimatethe twist angle to vary between 1.29◦ and 1.45◦ on thissample. Sample imaging and measurement results arepresented in Supplemental Material S7 [31].Device D5 on this device, we did not obtain a clearadditional Landau fan, but the lever arm could be esti-mated from the charge neutrality point. The twist anglewas estimated from the position of high resistance fea-tures found in transport measurements, which occur atν = ±8 in the low-twist-angle regime of this device [8].We find a twist angle of 0.65◦, and the results on thisdevice are presented in Supplemental Material S5 [31].Device D6 This device incorporated an additionalWSe2 (HQ graphene) into the stack which was pickedup after picking up the tBLG. The tBLG graphene wasrelaxed back to near-Bernal stacking, since no evidenceof moiré induced satellite peaks was found. This deviceserves as a control sample, which is presented in the Sup-plemental Material S6 [31].3. Experimental setupElectrical characterization of all samples was per-formed in a 3He/4He wet dilution refrigerator (OxfordKelvinoxMX400) with a base temperature of 32.5 mK.All electrical wiring is directly connected to the samplewith only a 1 kOhm pre-resistor in the BNC connectorbox used to connect to the instruments. The home-buildamplifiers and IV converters are each placed in a shieldedbox outside the refrigerator. The total resistance be-tween the BNC connector box to the sample holder is1.24 kΩ for each line at room temperature. An out-of-plane magnetic field is applied with a superconductingmagnet mounted in the liquid helium bath. We locatethe sample inside the coil of the magnet while avoidingmagnetic materials in the insert to ensure a uniform mag-netic field.144. Data aquisition and analysisMeasurement of temperature-dependent re-sistance: The temperature-dependent resistance inFigs. 1(d)–1(e) was measured using a homebuild IV-converter with a gain of 10 × 106 connected to two sidecontacts. This IV converter also applied a symmetricAC bias voltage to these contacts with an amplitude of100 µV (RMS) through a 1/10,000 voltage divider. Theexcitation signal was supplied by a Stanford SR830 lockinamplifier operating at 69.6 Hz and this apparatus also de-tected the signal from the IV converter. On the oppositeside-contact pair, a differential amplifier with a gain of1000 is connected to measure the voltage drop and thissignal is measured with a second lockin amplifier. A volt-age source (Yokogawa 7651) was connected to the gatethrough a 1 MΩ resistor. The measurement was per-formed during the condensing of the mixture, circulatingof the mixture and cooldown to base temperature of thedilution refrigerator. During this procedure the gate volt-age was constantly swept and at each data point the tem-perature on the mixing chamber plate is recorded. Dur-ing the cooldown procedure the mixing chamber rapidlycools once the condensed mixture enters, therefore reli-able measurements of the temperature between ∼ 2.5 Kand ∼ 4.5 K were not possible. The resulting data isplotted on a meshed grid with interpolation using theMatlab ’pcolor’ function.2-point conductance measurements: The conduc-tance traces used to construct the power spectra are ob-tained by a 2-point measurement of the conductance us-ing the homebuild IV converter and an AC bias of 100µV RMS at a frequency of 69.6 Hz. The voltage sourcefor the gate (Yokogawa 7651) was set to a range of 10 V,giving a resolution of 100 µV, and this is supplied to thesystem through a twisted-shielded cable with a low-passfilter with a cut-off frequency of 1.6 kHz to reduce outputnoise. The resolution of 100 µV is also the step size atwhich the gate voltage is swept, resulting in a Nyquistfrequency of 5000 V−1. For the bias spectroscopy datain Fig. 3, an additional voltage source was used to applya symmetric DC bias.Determination of the power spectra: To calculatethe power spectra Pω, we first perform a windowed auto-correlation on the gradient of the conductance with theMatlab function corrgram [76]. We use a window of 200samples, with an maximum lag of 200 points and overlapbetween the windows of 90%. Note, that with this win-dow we can only resolve frequency components fg > 50V−1. We then calculate the Fourier transform of the au-tocorrelation function using the Matlab implementation’ezfft’ [77]. By the convolution theorem this approach isequivalent to calculating |F{dG/dVg}|2.Fitting to the power spectra: To fit Eq. (5) to ourdata, we find the maxima in the power spectra at eachgate voltage. Due to the limited frequency resolution,this can lead to multiple gate voltage values having thesame capacitance in Fig. 7, but this is not a problem forthe least squares fitting procedure. We then load thisinto the curve fitting tool ’cftool’ in Matlab. We thenmanually select the correct maxima that belong to a sin-gle frequency component, and fit Eq. (5). cftool auto-matically calculates the 95% confidence interval, whichdefines the error in the fitting results presented in thiswork. We find the hBN thickness d from AFM scanningthe stack before fabrication, their values are shown in Ta-ble A1 under the column: ”bottom hBN”. The lever armα is extracted during the extraction of the twist angle,as described above.Magnetic field dependence: For Figs. 4(b) the fastFourier transform was calculated at each magnetic fieldusing the standard FFT implementation in Matlab. Thisfast Fourier transform gives a complex number, of whichthe angle can be calculated to obtain the phase. Theresulting phase was unwrapped to obtain the continuoussignal shown in Fig. A1(a). To find the frequency of thephase in 1/B as shown in the insert of Fig. A1(a), thephase signal was interpolated using a spline interpolationon a grid that is equally spaced in 1/B. After this, theFFT can be calculated to find the quantum oscillationfrequency. The resulting frequencies are shown as redsquares in Figs. A1(b)–A1(c). The same interpolationapproach is taken on the magnetoresistance data (pre-sented in Supplemental Material S1 [31]). After this, wefind the maxima at each B-field and manually select thecorrect frequency components, which are shown as bluedots in Fig. A1(b).[1] J. M. B. Lopes dos Santos, N. M. R. Peres, and A. H.Castro Neto, Graphene Bilayer with a Twist: ElectronicStructure, Phys. Rev. Lett. 99, 256802 (2007).[2] E. Suárez Morell, J. D. Correa, P. Vargas, M. Pacheco,and Z. Barticevic, Flat Bands in Slightly Twisted BilayerGraphene: Tight-Binding Calculations, Phys. Rev. B 82,121407(R) (2010).[3] G. Trambly de Laissardière, D. Mayou, and L. Magaud,Localization of Dirac Electrons in Rotated Graphene Bi-layers, Nano Lett. 10, 804 (2010).[4] R. Bistritzer and A. H. MacDonald, Moiré Bands inTwisted Double-Layer Graphene, Proc. Natl. Acad. Sci.U.S.A. 108, 12233 (2011).[5] Y. Cao, J. Y. Luo, V. Fatemi, S. Fang, J. D. Sanchez-Yamagishi, K. Watanabe, T. Taniguchi, E. Kaxiras,and P. Jarillo-Herrero, Superlattice-Induced InsulatingStates and Valley-Protected Orbits in Twisted BilayerGraphene, Phys. Rev. Lett. 117, 116804 (2016).[6] A. Kerelsky, L. J. McGilly, D. M. Kennes, L. Xian,M. Yankowitz, S. Chen, K. Watanabe, T. Taniguchi,J. Hone, C. Dean, A. Rubio, and A. N. Pasupathy,Maximized Electron Interactions at the Magic Angle in15Twisted Bilayer Graphene, Nature 572, 95 (2019).[7] Y. Xie, B. Lian, B. Jäck, X. Liu, C.-L. Chiu, K. Watan-abe, T. Taniguchi, B. A. Bernevig, and A. Yazdani,Spectroscopic Signatures of Many-Body Correlations inMagic-Angle Twisted Bilayer Graphene, Nature 572, 101(2019).[8] Y. Cao, V. Fatemi, A. Demir, S. Fang, S. L. Tomarken,J. Y. Luo, J. D. Sanchez-Yamagishi, K. Watanabe,T. Taniguchi, E. Kaxiras, R. C. Ashoori, and P. Jarillo-Herrero, Correlated Insulator Behaviour at Half-Fillingin Magic-Angle Graphene Superlattices, Nature 556, 80(2018).[9] Y. Choi, J. Kemmer, Y. Peng, A. Thomson, H. Arora,R. Polski, Y. Zhang, H. Ren, J. Alicea, G. Refael, F. vonOppen, K. Watanabe, T. Taniguchi, and S. Nadj-Perge,Electronic Correlations in Twisted Bilayer Graphenenear the Magic Angle, Nat. Phys. 15, 1174 (2019).[10] Y. Jiang, X. Lai, K. Watanabe, T. Taniguchi, K. Haule,J. Mao, and E. Y. Andrei, Charge Order and BrokenRotational Symmetry in Magic-Angle Twisted BilayerGraphene, Nature 573, 91 (2019).[11] J. Falson, I. Sodemann, B. Skinner, D. Tabrea,Y. Kozuka, A. Tsukazaki, M. Kawasaki, K. von Klitzing,and J. H. Smet, Competing Correlated States Around theZero-Field Wigner crystallization Transition of Electronsin Two Dimensions, Nat. Mater. 21, 311 (2022).[12] B. Tanatar and D. M. Ceperley, Ground State of theTwo-Dimensional Electron Gas, Phys. Rev. B 39, 5005(1989).[13] Y. Cao, V. Fatemi, S. Fang, K. Watanabe, T. Taniguchi,E. Kaxiras, and P. Jarillo-Herrero, Unconventional Su-perconductivity in Magic-Angle Graphene Superlattices,Nature 556, 43 (2018).[14] U. Zondiner, A. Rozen, D. Rodan-Legrain, Y. Cao,R. Queiroz, T. Taniguchi, K. Watanabe, Y. Oreg, F. vonOppen, A. Stern, E. Berg, P. Jarillo-Herrero, and S. Ilani,Cascade of Phase Transitions and Dirac Revivals inMagic-Angle Graphene, Nature 582, 203 (2020).[15] M. Yankowitz, S. Chen, H. Polshyn, Y. Zhang, K. Watan-abe, T. Taniguchi, D. Graf, A. F. Young, and C. R. Dean,Tuning Superconductivity in Twisted Bilayer Graphene,Science 363, 1059 (2019).[16] X. Lu, P. Stepanov, W. Yang, M. Xie, M. A. Aamir,I. Das, C. Urgell, K. Watanabe, T. Taniguchi, G. Zhang,A. Bachtold, A. H. MacDonald, and D. K. Efetov, Su-perconductors, Orbital Magnets and Correlated States inMagic-Angle Bilayer Graphene, Nature 574, 653 (2019).[17] A. L. Sharpe, E. J. Fox, A. W. Barnard, J. Finney,K. Watanabe, T. Taniguchi, M. A. Kastner, andD. Goldhaber-Gordon, Emergent Ferromagnetism nearThree-Quarters Filling in Twisted Bilayer Graphene, Sci-ence 365, 605 (2019).[18] M. Serlin, C. L. Tschirhart, H. Polshyn, Y. Zhang, J. Zhu,K. Watanabe, T. Taniguchi, L. Balents, and A. F. Young,Intrinsic Quantized Anomalous Hall Effect in a MoiréHeterostructure, Science 367, 900 (2020).[19] N. Tilak, X. Lai, S. Wu, Z. Zhang, M. Xu,R. de Almeida Ribeiro, P. C. Canfield, and E. Y. Andrei,Flat Band Carrier Confinement in Magic-Angle TwistedBilayer Graphene, Nat. Commun. 12, 4180 (2021).[20] N. P. Kazmierczak, M. Van Winkle, C. Ophus, K. C.Bustillo, S. Carr, H. G. Brown, J. Ciston, T. Taniguchi,K. Watanabe, and D. K. Bediako, Strain Fields inTwisted Bilayer Graphene, Nat. Mater. 20, 956 (2021).[21] A. Schäpers, J. Sonntag, L. Valerius, B. Pestka, J. Stras-das, K. Watanabe, T. Taniguchi, L. Wirtz, M. Morgen-stern, B. Beschoten, R. J. Dolleman, and C. Stampfer,Raman Imaging of Twist Angle Variations in TwistedBilayer Graphene at Intermediate Angles, 2D Mater. 9,045009 (2022).[22] C. N. Lau, M. W. Bockrath, K. F. Mak, and F. Zhang,Reproducibility in the Fabrication and Physics of MoiréMaterials, Nature 602, 41 (2022).[23] A. Uri, S. Grover, Y. Cao, J. A. Crosse, K. Bagani,D. Rodan-Legrain, Y. Myasoedov, K. Watanabe,T. Taniguchi, P. Moon, M. Koshino, P. Jarillo-Herrero,and E. Zeldov, Mapping the Twist-Angle Disorder andLandau Levels in Magic-Angle Graphene, Nature 581,47 (2020).[24] J. M. Park, Y. Cao, K. Watanabe, T. Taniguchi, andP. Jarillo-Herrero, Flavour Hund’s Coupling, Chern Gapsand Charge Diffusivity in Moiré Graphene, Nature 592,43 (2021).[25] C. R. Dean, A. F. Young, I. Meric, C. Lee, L. Wang,S. Sorgenfrei, K. Watanabe, T. Taniguchi, P. Kim, K. L.Shepard, and J. Hone, Boron nitride substrates for high-quality graphene electronics, Nat. Nanotechnol. 5, 722(2010).[26] R. Ribeiro-Palau, S. Chen, Y. Zeng, K. Watan-abe, T. Taniguchi, J. Hone, and C. R. Dean, High-Quality Electrostatically Defined Hall Bars in MonolayerGraphene, Nano Lett. 19, 2583 (2019).[27] E. Icking, L. Banszerus, F. Wörtche, F. Volmer,P. Schmidt, C. Steiner, S. Engels, J. Hesselmann,M. Goldsche, K. Watanabe, T. Taniguchi, C. Volk,B. Beschoten, and C. Stampfer, Transport Spectroscopyof Ultraclean Tunable Band Gaps in Bilayer Graphene,Adv. Electron. Mater. 8, 2200510 (2022).[28] L. Wang, I. Meric, P. Y. Huang, Q. Gao, Y. Gao, H. Tran,T. Taniguchi, K. Watanabe, L. M. Campos, D. A. Muller,J. Guo, P. Kim, J. Hone, K. L. Shepard, and C. R.Dean, One-Dimensional Electrical Contact to a Two-Dimensional Material, Science 342, 614 (2013).[29] T. Uwanno, T. Taniguchi, K. Watanabe, and K. Na-gashio, Electrically Inert h-BN/Bilayer Graphene Inter-face in All-Two-Dimensional Heterostructure Field Ef-fect Transistors, ACS Appl. Mater. Interfaces 10, 28780(2018).[30] M. Schmitz, S. Engels, L. Banszerus, K. Watanabe,T. Taniguchi, C. Stampfer, and B. Beschoten, High Mo-bility Dry-Transferred CVD Bilayer Graphene, Appl.Phys. Lett. 110, 263110 (2017).[31] See Supplemental Material at xxx for additional bulktransport data on the 750-nm-wide constriction (S1),dependence of the Coulomb oscillations on the samplewidth (S2), the fitting results to determine the value ofS (S3), the remaining power spectra and analysis with afixed value of S (S4), additional results on a 0.65◦-twistedsample D5 (S5), control experiments in the absence ofmoiré on sample D6 (S6), observation of Coulomb oscil-lations in sample D4 (S7), results on the Coulomb oscilla-tions with magnetic field reversal on sample D4 (S8); andadditional results on the negative electronic compressibil-ity in a magnetic field (S9).[32] P. Stepanov, I. Das, X. Lu, A. Fahimniya, K. Watanabe,T. Taniguchi, F. H. L. Koppens, J. Lischner, L. Levitov,and D. K. Efetov, Untying the Insulating and Supercon-ducting Orders in Magic-Angle Graphene, Nature 583,16375 (2020).[33] J. M. Park, Y. Cao, L.-Q. Xia, S. Sun, K. Watanabe,T. Taniguchi, and P. Jarillo-Herrero, Robust Supercon-ductivity in Magic-Angle Multilayer Graphene Family,Nat. Mater. 21, 877 (2022).[34] K. P. Nuckolls, M. Oh, D. Wong, B. Lian, K. Watanabe,T. Taniguchi, B. A. Bernevig, and A. Yazdani, StronglyCorrelated Chern Insulators in Magic-Angle Twisted Bi-layer Graphene, Nature 588, 610 (2020).[35] S. L. Tomarken, Y. Cao, A. Demir, K. Watanabe,T. Taniguchi, P. Jarillo-Herrero, and R. C. Ashoori, Elec-tronic Compressibility of Magic-Angle Graphene Super-lattices, Phys. Rev. Lett. 123, 046601 (2019).[36] P. Stepanov, M. Xie, T. Taniguchi, K. Watanabe, X. Lu,A. H. MacDonald, B. A. Bernevig, and D. K. Efetov,Competing Zero-Field Chern Insulators in Superconduct-ing Twisted Bilayer Graphene, Phys. Rev. Lett. 127,197701 (2021).[37] Y. Saito, J. Ge, L. Rademaker, K. Watanabe,T. Taniguchi, D. A. Abanin, and A. F. Young, HofstadterSubband ferromagnetism and Symmetry-Broken Cherninsulators in Twisted Bilayer Graphene, Nat. Phys. 17,478 (2021).[38] S. Wu, Z. Zhang, K. Watanabe, T. Taniguchi, and E. Y.Andrei, Chern Insulators, van Hove Singularities andTopological Flat Bands in Magic-Angle twisted BilayerGraphene, Nat. Mater. 20, 488 (2021).[39] C. Stampfer, E. Schurtenberger, F. Molitor, J. Güttinger,T. Ihn, and K. Ensslin, Tunable Graphene Single Elec-tron Transistor, Nano Lett. 8, 2378 (2008).[40] A. A. M. Staring, H. van Houten, C. W. J. Beenakker,and C. T. Foxon, Coulomb-Blockade Oscillations in Dis-ordered Quantum Wires, Phys. Rev. B 45, 9222 (1992).[41] C. W. J. Beenakker, Theory of Coulomb-Blockade Os-cillations in the Conductance of a Quantum Dot, Phys.Rev. B 44, 1646 (1991).[42] J. Güttinger, C. Stampfer, F. Libisch, T. Frey,J. Burgdörfer, T. Ihn, and K. Ensslin, Electron-HoleCrossover in Graphene Quantum Dots, Phys. Rev. Lett.103, 046810 (2009).[43] M. Eich, F. Herman, R. Pisoni, H. Overweg, A. Kurz-mann, Y. Lee, P. Rickhaus, K. Watanabe, T. Taniguchi,M. Sigrist, T. Ihn, and K. Ensslin, Spin and Valley Statesin Gate-Defined Bilayer Graphene Quantum Dots, Phys.Rev. X 8, 031023 (2018).[44] L. Banszerus, T. Fabian, S. Möller, E. Icking, H. Heim-ing, S. Trellenkamp, F. Lentz, D. Neumaier, M. Otto,K. Watanabe, T. Taniguchi, F. Libisch, C. Volk, andC. Stampfer, Electrostatic Detection of Shubnikov–deHaas Oscillations in Bilayer Graphene by Coulomb Res-onances in Gate-Defined Quantum Dots, Phys. StatusSolidi B 257, 2000333 (2020).[45] C. Gold, A. Kurzmann, K. Watanabe, T. Taniguchi,K. Ensslin, and T. Ihn, Scanning Gate Microscopy of Lo-calized States in a Gate-Defined Bilayer Graphene Chan-nel, Phys. Rev. Res. 2, 043380 (2020).[46] J. P. Eisenstein, L. N. Pfeiffer, and K. W. West, NegativeCompressibility of Interacting Two-Dimensional Electronand Quasiparticle Gases, Phys. Rev. Lett. 68, 674 (1992).[47] K. Steffen, Compressibility and Capacitance of QuantumSystems, Ph.D. thesis, Universität Augsburg (2017).[48] G. Meissner, H. Namaizawa, and M. Voss, Stability andImage-Potential-Induced Screening of Electron Vibra-tional Excitations in a Three-Layer Structure, Phys. Rev.B 13, 1370 (1976).[49] L. Bonsall and A. A. Maradudin, Some Static and Dy-namical Properties of a Two-Dimensional Wigner Crys-tal, Phys. Rev. B 15, 1959 (1977).[50] M. S. Bello, E. I. Levin, B. I. Shklovskii, and A. L. Efros,Density of Localized States in the Surface Impurity Bandof a Metal-Insulator-Semiconductor Structure, J. Exp.Theor. Phys. 53, 822 (1981).[51] I. M. Ruzin, S. Marianer, and B. I. Shklovskii, Pinning ofa Two-Dimensional Wigner Crystal by Charged Impuri-ties, Phys. Rev. B 46, 3999 (1992).[52] H. Fu, B. I. Shklovskii, and B. Skinner, Correlation Ef-fects in the Capacitance of a Gated Carbon Nanotube,Phys. Rev. B 91, 155118 (2015).[53] S. Joy and B. Skinner, Wigner Crystallization at LargeFine Structure Constant, Phys. Rev. B 106, L041402(2022).[54] B. Skinner and B. I. Shklovskii, Anomalously Large Ca-pacitance of a Plane Capacitor with a Two-DimensionalElectron Gas, Phys. Rev. B 82, 155111 (2010).[55] L. Wang, Y. Wang, X. Chen, W. Zhu, C. Zhu, Z. Wu,Y. Han, M. Zhang, W. Li, Y. He, et al., Negative quan-tum capacitance induced by midgap states in single-layergraphene, Scientific reports 3, 1 (2013).[56] S. Ilani, L. A. K. Donev, M. Kindermann, and P. L.McEuen, Measurement of the Quantum Capacitance ofInteracting Electrons in Carbon Nanotubes, Nat. Phys.2, 687 (2006).[57] N. Nagaosa, J. Sinova, S. Onoda, A. H. MacDonald, andN. P. Ong, Anomalous Hall effect, Rev. Mod. Phys. 82,1539 (2010).[58] G. Rai, L. Crippa, D. Călugăru, H. Hu, L. d. Medici,A. Georges, B. A. Bernevig, R. Valent́ı, G. Sangio-vanni, and T. Wehling, Dynamical correlations andorder in magic-angle twisted bilayer graphene, arXiv10.48550/arXiv.2309.08529 (2023), 2309.08529.[59] B. Spivak and S. A. Kivelson, Phases Intermediate Be-tween a Two-Dimensional Electron Liquid and WignerCrystal, Phys. Rev. B 70, 155114 (2004).[60] M. Roger, Multiple Exchange in 3He and in the WignerSolid, Phys. Rev. B 30, 6432 (1984).[61] I. Das, X. Lu, J. Herzog-Arbeitman, Z.-D. Song,K. Watanabe, T. Taniguchi, B. A. Bernevig, andD. K. Efetov, Symmetry-broken Chern insulators andRashba-like Landau-level crossings in magic-angle bilayergraphene, Nat. Phys. 17, 710 (2021).[62] Dataset: Negative electronic compressibility in charge is-lands in twisted bilayer graphene, Zenodo 10.5281/zen-odo.7409456 (2024).[63] W. Albrecht, J. Moers, and B. Hermanns, HNF -Helmholtz Nano Facility, JLSRF 3, A112 (2017).[64] L. Onsager, Interpretation of the de Haas-van Alphen Ef-fect, The London, Edinburgh, and Dublin PhilosophicalMagazine and Journal of Science 43, 1006 (1952).[65] T. Fabian, M. Kausel, L. Linhart, J. Burgdörfer, andF. Libisch, Half-Integer Wannier Diagram and Brown-Zak Fermions of Graphene on Hexagonal Boron Nitride,Phys. Rev. B 106, 165412 (2022).[66] N. N. T. Nam and M. Koshino, Lattice Relaxation andEnergy Band Modulation in Twisted Bilayer Graphene,Phys. Rev. B 96, 075311 (2017).[67] M. Koshino, N. F. Q. Yuan, T. Koretsune, M. Ochi,K. Kuroki, and L. Fu, Maximally localized wannier or-bitals and the extended hubbard model for twisted bi-17layer graphene, Phys. Rev. X 8, 031087 (2018).[68] A. Fischer, L. Klebl, C. Honerkamp, and D. M.Kennes, Spin-fluctuation-induced Pairing in Twisted Bi-layer Graphene, Phys. Rev. B 103, L041103 (2021).[69] L. Klebl, Z. A. H. Goodwin, A. A. Mostofi, D. M. Kennes,and J. Lischner, Importance of Long-ranged Electron-electron Interactions for the Magnetic Phase Diagram ofTwisted Bilayer Graphene, Phys. Rev. B 103, 195127(2021).[70] F. Guinea and N. R. Walet, Electrostatic Effects, BandDistortions, and Superconductivity in Twisted GrapheneBilayers, Proc. Natl. Acad. Sci. U.S.A. 115, 13174 (2018).[71] Z. A. H. Goodwin, V. Vitale, X. Liang, A. A. Mostofi,and J. Lischner, Hartree Theory Calculations of Quasi-particle Properties in Twisted Bilayer Graphene, Elec-tron. Struct. 2, 034001 (2020).[72] L. Rademaker, D. A. Abanin, and P. Mellado, ChargeSmoothening and Band Flattening Due to Hartree Cor-rections in Twisted Bilayer Graphene, Phys. Rev. B 100,205114 (2019).[73] T. Cea, N. R. Walet, and F. Guinea, Electronic BandStructure and Pinning of Fermi Energy to Van Hove Sin-gularities in Twisted Bilayer Graphene: A Self-consistentApproach, Phys. Rev. B 100, 205113 (2019).[74] Z. A. H. Goodwin, V. Vitale, F. Corsetti, D. K. Efe-tov, A. A. Mostofi, and J. Lischner, Critical Role of De-vice Geometry for the Phase Diagram of Twisted BilayerGraphene, Phys. Rev. B 101, 165110 (2020).[75] K. S. Novoselov, D. Jiang, F. Schedin, T. J. Booth,V. V. Khotkevich, S. V. Morozov, and A. K. Geim,Two-Dimensional Atomic Crystals, Proc. Natl. Acad. Sci.U.S.A. 102, 10451 (2005).[76] N. Marwan, Windowed Cross Correlation (corrgram),https://www.mathworks.com/matlabcentral/fileexchange/15299-windowed-cross-correlation-corrgram (accessed:11 Jan 2022).[77] F. Moisy, EZFFT: An Easy to Use PowerSpectrum (FFT), https://www.mathworks.com/matlabcentral/fileexchange/22163-ezfft-an-easy-to-use-power-spectrum-fft (accessed: 18 Dec 2019).