Fractional Chern Insulator States in Twisted Bilayer Graphene: An Analytical Approach

Recent experiments on bilayer graphene twisted near the magic angle have observed spontaneous integer quantum Hall states in the presence of an aligned hexagonal boron nitride (hBN) substrate. These states arise from valley ferromagnetism, and the complete filling of Chern bands. A natural question is whether fractional filling of the same bands would lead to fractional quantum Hall states, i.e. to fractional Chern insulators (FCIs). Here, we argue that the magic angle graphene bands have favorable quantum geometry for realizing FCI phases. We show that in the tractable `chiral' limit, the flat bands wavefunctions are an analytic function of the crystal momentum. This remarkable property fixes the quantum metric up to an overall momentum dependent scale factor, the local Berry curvature, whose variation is itself small. Thus the three conditions associated with FCI stability (i) narrow bands (ii) quantum metric satisfying the isotropic ideal droplet condition and (iii) relatively uniform Berry curvature are all satisfied. Our work emphasizes continuum real space approaches to FCIs in contrast to earlier works which mostly focused on tight binding models. This enables us to construct a Laughlin wavefunction on a real-space torus that is a zero energy ground state of pseudopotential-like interactions. The latter can be realized for Coulomb interaction in the limit of very short screening length. We also map the problem onto an equivalent one of Dirac particles in an inhomogeneous but periodic magnetic field and no crystal potential. Finally we discuss evolution of the band geometry on tuning away from the chiral limit and show numerically that some of the desirable properties continue to hold at a quantitative level for realistic parameter values.


I. INTRODUCTION
The recent discovery of correlated insulators and superconductors in bilayer graphene twisted near a magic angle has brought the physics of moiré superlattices to the fore. A recent striking observation is that of a quantized Hall response in magic angle bilayer graphene. This has been attributed to the presence of an aligned hBN substrate, which breaks the C 2 rotation symmetry relating opposite valleys generating bands with valley chern number [1][2][3][4][5] . Spontaneous spin and valley polarization then leads to a chern insulator with integer quantized Hall conductance. Similar physics is also expected in other Moire platforms where C 2 symmetry is broken such as ABC trilayer graphene on hBN [6][7][8] and twisted double bilayer graphene [9][10][11][12] .
In this work we will be mainly concerned with the ground state at fractional electron filling (assuming spin and valley polarization) of the nearly flat bands of magic angle graphene with aligned hBN. Of course, in fractionally filled Landau levels (e.g. electrons in the lowest Lan-dau level at 1/3 filling) fractional quantum Hall states arise. In fact, in the related problem of strained graphene, where Landau levels of opposite effective magnetic field arise from the opposite valleys 13 , a fractional filling was numerically shown to lead to spontaneous valley polarization and fractional quantum Hall states 14 . This is particularly relevant for twisted bilayer graphene due to the close analogy between the nearly flat bands of magic angle graphene and Landau levels which was highlighted in several recent works [15][16][17][18][19][20] , and is all the more apparent in the presence of a hBN substrate that leads to a staggered sublattice potential and isolates Chern bands in each valley. However, the flat band wavefunctions differ crucially from the Landau levels due to the presence of actual (Moiré) translational symmetry instead of continuous magnetic translations. As a result, the relevant problem to TBG is that of fractional Chern insulators (FCIs) which investigates the necessary ingredients for realizing the fractional quantum hall phases in systems with discrete translation symmetry.
Previous studies of FCI physics have performed numerical and theoretical analysis primarily of tight-binding models, see the reviews 21,22 . Tight binding models can be related to lowest Landau level (LLL) physics by mapping the LLL Hilbert space and pseudopotential interactions to that of a lattice model 23,24 . This "ideal" lattice model will necessarily host an FCI state, and one may test adiabatic continuity between it and a realistic model. 25 However, the construction of the mapping depends on an arbitrary k-dependent phase choice and the ideal lattice models one obtains are often unrealistic.
While many numerical studies (see numerics discussed and cited in the comprehensive reviews 21,22 ) have successfully demonstrated that FCI states exist and agree with Chern-Simons effective field theory, microscopic theoretical analysis on FCI stability has been hampered by the fact that tight binding models lack continuous magnetic translation symmetry and rotation symmetry due to the discreteness of the lattice. For this reason, most analytic work has focused on identifying ideal characteristics of FCI models in momentum-space.
While flat bands and non-zero Chern number are obvious conditions for an ideal FCI state, band geometry also plays a role. It was noticed early on that flat Berry curvature was important for realizing FCI states; for example Ref. 26 showed its importance for the long wavelength limit of the Girvin-MacDonald-Platzmann (GMP) algebra 27 , a Hamiltonian formulation of the quantum Hall problem presented entirely in terms of projected density operators that obey said algebra. However, higher Landau levels do not display the same fractional quantum hall effect (FQHE) as the LLL despite having flat bands and flat Berry curvature. This suggests that other band geometric conditions involving relationships between the Fubini-Study metric and the Berry curvature play a crucial role in stabilizing FCIs by mimicking the band geometry of the LLL. Under these stricter band-geometric conditions 21,28 , the GMP algebra is satisfied to all orders. One may also construct a first quantized description of Chern bands in momentum space by assuming these geometric conditions 29 . Exact diagonalization studies have shown that band geometry correlates strongly with the stability of FCI phases 30 . Motivated by this, some models with ideal band geometry have been categorized and fine-tuned to have nearly-flat Berry curvature 31 .
Here, we will show that chiral magic angle TBG (chiral-MATBG) has nearly ideal band geometry for FCI phases. Due to the unique properties of the TBG bands which we will explain in detail in this work, models like it were missed in previous discussions of ideal FCI bands. More specifically, we will show that chiral-MATBG has relatively uniform Berry curvature and satisfies the so-called isotropic ideal droplet condition (also referred to as the trace condition), discussed and defined below.
The ideal criteria above are phrased in momentum space. Indeed, most toy FCI models are built on real space tight binding lattices whose discreteness makes geometric arguments in real space difficult. However, we will use the continuum-nature of the model to show that a first quantized description in real space is possible. We will construct a Laughlin state on a realspace torus, built from chiral-MATBG wavefunctions, that is a zero-energy ground state for pseudopotential-like interactions. To our knowledge, such analysis has not appeared in the FCI literature due to the focus on tight-binding models, but see Ref. 32 for the FQHE in a periodic potential. We will also show that the chiral-MATBG wavefunctions may be interpreted as wavefunctions of a generalized LLL, in which the magnetic field is modulated in space, in the sense that density-density interactions have identical matrix elements for both systems. This provides a natural pathway to test adiabatic continuity from the LLL to potential realistic FCI states in TBG and a potential new tool for theoretical FCI analysis.
In the next section we review the chiral model for TBG 33 to set notation and results that will be convenient for the rest of the paper. We also discuss the introduction of a sublattice potential that separates the two flat bands so that the lower (or upper) band can host a FCI state. If the sublattice potential is chosen to be the same on both layers the chiral-MATBG wavefunctions are unchanged but the two bands are gapped from each other. In section III, we re-view band-geometric conditions for FCIs and show that chiral-MATBG has near-ideal band geometry; the non-ideality being inhomogeneous Berry curvature. In section IV, we construct a Laughlin ground state on a real space torus 34 which is an exact zero-energy ground state for pseudopotential-like interactions 35 . In section V we show that the chiral-MATBG wavefunctions are very similar to wavefunctions for a Dirac particle in an inhomogeneous magnetic field. Indeed, the only difference is a position dependent phase and layer polarization for the chiral-MATBG wavefunctions which does not matter for matrix elements of density operators. Finally, in section VI we show quantitatively using numerics how the momentum space energetics and geometry change as one moves from the chiral model to more realistic models for TBG. The chiral model, characterized by a parameter κ = w 0 /w 1 = 0, has the most promising band energetics and geometry. As κ increases to realistic values, ≈ 0.7−0.8 36,37 , these characteristics worsen, though the chiral model bands are still adiabatically connected to realistic TBG bands. Because FCI states are stable to small perturbations and the interaction screening length is highly tunable 38,39 , we argue that FCI states could be realizable in realistic models of TBG for some screening lengths. Nonetheless, it is worth noting that if κ can be lowered by modifying the system (it can be raised by adding dopants 40 ) then this would likely have a positive impact on FCI stability. Pressure may have this effect because it could increase lattice relaxation.

II. THE CHIRAL MODEL FOR TBG WITH A SUBLATTICE POTENTIAL
The continuum model for a single spin and valley of TBG is given by the Hamiltonian 33,41 which acts on a spinor Ψ(r) = (ψ 1 , χ 1 , ψ 2 , χ 2 ).
The lattice vectors and triangular lattice coordinates and reciprocal lattice analogues are The chiral model (κ = 0) enjoys many simplifications. First, the θ dependence in the Pauli matrices can be removed by a gauge transformation such that σ ±θ/2 → σ, leading to an exact particle-hole symmetry 15,20 . Second and most importantly, the single-particle Hamiltonian anticommutes with the chirality operator σ z which denotes the graphene sublattice. Additionally, the Hamiltonian can now be written in a dimensionless form H chiral = E 0 H chiral where E 0 = vk θ and The Hamiltonian now acts on the spinor Ψ(r) = (ψ 1 , ψ 2 , χ 1 , χ 2 ) T . The anti-holomorphic derivative is given by∂ = 1 2 (∂ x + i∂ y ), corresponding to the complex coordinate z = (x + iy). The parameter α = w/E 0 ∼ 1/ sin θ measures the strength of tunneling.
At certain values of α, corresponding to magic angles for θ, the Hamiltonian (4) exhibits two exactly flat bands at zero energy. This occurs for infinitely many α, with a quasiperiodicity α n ≈ α n−1 + 3/2 33 . We now obtain the wavefunctions for the exactly flat bands at these values of α and show that they exhibit very interesting properties.
We are interested in zero modes of the Hamiltonian (4). Because the "chirality" matrix anticommutes with the Hamiltonian (4), we can label the zero energy solutions with their σ z eigenvalue. We therefore seek solutions to for the positive chirality, σ z = +1, eigenspace. The solutions for the negative chirality eigenspace are then χ(r) = ψ * (−r). We pause to recall that the Hamiltonian (4) is not in traditional Bloch form; the two layers have a crystal momentum that differs by q 1 . The "Bloch" states, of positive chirality, are then of the form ψ k (r) = e ik·r u 1k (r) u 2k (r)e iq 1 ·r .
We take k = 0 to correspond to the K point in the moiré Brillouin zone (BZ). There is guaranteed to be a zero-mode at the K point, k = 0. This is because there is one for α = 0 and it is guaranteed to exist at all α because it is protected and pinned by symmetries 33 . Another way to understand this is by noting that the two different wavefunctions at K transform under conjugate representations of C 3 which are glued together by C 2 T symmetry 42 . The chiral symmetry ensures that this pair of levels remains at zero (since it interchanges positive and negative energy states). We therefore have D(r)ψ K (r) = 0, where ψ K (r) is the K point wavefunction for positive chirality.
To find the other zero modes, we use that D(r)(F (z)ψ(r)) = 0 if D(r)ψ(r) = 0 for any holomorphic function F (z) 33 . Since the K-point wavefunction is a zero mode, we may generate other zero modes by multiplying it by some F (z). The wavefunctions F (z)ψ K (r) are generically not normalizable since Liouville's theorem implies that a non-constant F (z) cannot be bounded. However, at the magic angles, we have ψ K (r 0 ) = 0 (both components), for r 0 = (a 1 − a 2 )/3. 33 We can therefore let F (z) be meromorphic, provided it only has poles when z = z 0 = (a 1 − a 2 )/3 and translates by lattice vectors z → z + a 1 , z → z + a 2 . Here we use a non-bold letter to denote the complex version of a vector, e.g.: The chiral-MATBG wavefunctions are then of the form where ω = e iφ = a 2 /a 1 . Here f (z) is a holomorphic function that is not and does not need to be bounded. The function 10) is the Jacobi theta function of the first kind, is odd in z, and has zeros when z = m + nω where m, n are integers 43 . Under translations it satisfies 43 For Bloch-periodic states we choose where k = k x + ik y . The Bloch-periodicity then follows from the behavior of the theta functions under translations (11) 43 . The conventions here and that of Ref. 33 agree up to k-dependent normalizations. We choose this gauge because the moiré-periodic functions (13) are holomorphic in k. The existence of such a gauge originates from the fact that D(r) only has antiholomorphic derivatives. The zero mode equation D(r)ψ k (r) = 0 can then be recast as where∂ k =∂ + k. The function u k (r) can then be chosen to only depend on k and notk, as it is the only zero mode for a matrix that only depends on k.
Introducing a sublattice potential that is equal on both layers then creates a promising situation for a fractional Chern insulator state without changing wavefunctions (9). Indeed, such a potential has the effect where δ = ∆/E 0 = α∆/w is the dimensionless potential strength. Because the flat band states can be chosen to be eigenstates of σ z , introducing this potential splits the two bands without changing the wavefunctions. A potential on a single layer would change the wavefunctions, but it would still not introduce mixing between the flat bands because it commutes with σ z . It would therefore only change the wavefunctions via mixings with higher bands. Such effects are suppressed by factors of δ/E 0 . For the rest of this paper, we therefore work with split flat bands and focus on the lower band with wavefunctions (9), (12). The positive and negative chirality bands have Chern number −1 and +1 respectively. As they are perfectly flat and separated in the presence of a sublattice potential (15), we attain ideal band energetics for hosting a FCI state by fractionally filling the positive chirality band. We assume here that the system is spin and valley polarized for simplicity, which in realistic systems is achieved through interaction-induced flavor polarization [1][2][3][4][5]44 , though spin unpolarized states can be included in the same way as in the LLL. In the next section we discuss the near-ideal band geometry of the flat bands.

III. BAND GEOMETRY OF THE FLAT BANDS
In this section we discuss the band geometry of the wavefunctions (9). In particular, we will show that chiral-MATBG has near-ideal band geometry due to holomorphicity in momentum space (13). We emphasize that band geometry is a strong predictor of FCI stability: it appears in many analytical arguments for ideal descriptions and its deviation from ideality strongly correlates with the many body gap 30 .
Before we move forward, we review the various band-geometric conditions. All geometric information is encapsulated in the complex tensor (16) which does not depend on the k-dependent norm or phase of the wavefunctions u k . We use moiréperiodic wavefunctions in these overlaps as is standard. The imaginary and real components of this tensor respectively capture the phase and magnitude of Bloch overlaps under varying k and correspond to the Berry curvature Ω(k) and the Fubini-Study metric g ab (k) respectively: To motivate these geometrical considerations, we review the GMP algebra in the LLL 27 . The many-body quantum hall Hamiltonian projected to the LLL can be written as where ρ q is obtained from the density operator ρ q = i e iq·ri by LLL projection: The coordi-nateR i is the guiding center coordinate for the i'th particle and the pre-factor of e −q 2 2 is known as the LLL form factor. Because the guiding center coordinates do not commute, R x ,R y = −i 2 B , the projected densities do not either and can be shown to obey the GMP algebra 27 It was noticed early on that to mimic Landau levels, Berry curvature should be flat. In the context of the GMP algebra, this allows one to reproduce the algebra at small q, in particular one obtains , are the projected densities in the Chern band such that The algebra (20) matches (19) in the small q limit, with Ω as the analogue of 2 B . The algebra no longer closes, however.
One must add conditions on the metric g(k) to reproduce the algebra in full. The full algebra (19) is recovered, with Ω taking the role of 2 B , if and only if the metric is k-independent and satisfies for each k in the BZ. In particular, it is recovered for the operators ρ q = e q 2 Ω/2ρ q , such that LLL physics is reproduced 21,28 . Other Landau levels crucially have a different form factor.
The condition (21) holds for all two band tight binding models, although the metric is usually not constant and its eigenvectors typically vary between different values of k. The condition (21) is referred to as the determinant condition or the ideal droplet condition. If the metric varies only in magnitude then one may diagonalize it, which leads to a stronger version of (21), This condition is considerably harder to satisfy and is referred to as the trace condition, as it saturates the inequality tr g ≥ |Ω(k)| 28 . It is also known as the isotropic ideal droplet condition.
The real-space periodic Bloch states (12) have the remarkable property that they are holomorphic in k x + ik y . It was noticed in Ref. 29 that this type of holomorphicity implies (22) in tight-binding models. The proof is similar in our continuum model, which we show now for completeness. Holomorphicity in k x + ik y implies ∂ kx u k = (∂ +∂)u k = ∂u k , where the holomorphic derivatives are now in k-space, 2∂ = from which (22) follows.
Refs. 29 and 31 previously considered ideal models satisfying (22) by writing down wavefunctions that are holomorphic in k x + ik y . If one imposes periodic boundary conditions in kspace, these are elliptic functions which have |C| ≥ 2 because they have two or more poles in the BZ 45 . The authors of Refs. 29 and 31 thus restricted their attention to elliptic function models with |C| ≥ 2. However, periodic boundary conditions in k-space are not required for a well defined Berry curvature or Chern number. The chiral-MATBG Bloch states (12) are not periodic, and they were therefore missed in the above studies of ideal models satisfying (22).
The only non-ideal aspect of chiral-MATBG is inhomogeneous Berry curvature. Because Berry curvature spreads out in the chiral limit (and outside of the chiral limit in the presence of a sublattice potential), this is not a large effect. We can quantify Berry curvature fluctuations through where A BZ is the area of the BZ and C = −1 is the Chern number for the positive chirality band in this valley. For the first through fourth magic angles we have F = 0.373, 0.197, 0.173, 0.213 respectively. This is illustrated in Fig. 1 which shows the distribution of Berry curvature in momentum space for the first four magic angles showing relatively uniform Berry curvature, particularly for the second and third magic angles.
To compare with other models it will be useful to quantify violation of (22). We use the BZ average T = T (k) BZ of the (non-negative 28 ) quantity Because chiral-MATBG exactly satisfies (22) and has small F , it is one of the most promising models to host FCI states, even compared to fine-tuned models. The study of ideal FCI models has shown that the confluence of flat bands, constant Berry curvature, and the ideal metric (22) are very rare. These geometric qualities are important: Ref. 30 finds that the many body gap for FCI states is strongly negatively correlated with F and T . Because chiral-MATBG has small F and exactly satisfies (22), it is one of the most promising systems for FCI states. To argue this, we compare with models in Ref. 30. The results are shown in Table I; chiral-MATBG is just as promising as fine-tuned lattice models with only nearest neighbor hoppings, and next nearest neighbor hoppings for the Haldane model, included 30 . If arbitrary length hoppings with > 3 bands are chosen then both F and T can be made arbitrarily small 31 . It is remarkable that chiral-MATBG compares well to fine-tuned models without any parameter optimization.
While the above analysis suggests that chiral-MATBG is promising for FCI ground states, to make analytical progress that shows chiral-MATBG can host an FCI ground state one must pursue pseudopotential-like arguments. For the FQHE, Haldane 35 used rotation and magnetic translation symmetry to reduce all interaction matrix elements to a discrete collection of numbers V m , where m is the (conserved) relative angular momentum of a pair of particles. Break- ing these symmetries splits the degeneracy between these matrix elements 22 . One expects that the FQHE can be recovered, say at ν = 1/3, if the gap between V 1 and V 3 is much larger than the splitting of each. In the next section we pursue an argument similar to this: we construct a Laughlin state on a real-space torus geometry that is a zero-energy ground state of pseudopotential-like interactions. These pseudopotential-like interactions will not be described by a set of discrete matrix elements and will be split into bands 22 . For a screened coulomb interaction with screening length much less than a, we will argue that these bands are well separated.

IV. LAUGHLIN STATE ON A MOIRÉ TORUS
In this section we investigate the chiral model on a moiré torus and derive a Laughlin-like wavefunction. We closely follow Haldane and Rezayi's derivation for the FQHE 34 . In fact, we obtain an exact mapping onto their derivation which can be anticipated as follows. Previous studies of FCIs have identified the dictionary 2 B → Ω for translating from the FQHE to FCIs, where 2 B is the squared magnetic length and Ω is the average Berry curvature (for example: Refs. 21, 26, and 29). One may under-stand this in the context of Harper-Hofstadter models, where both a magnetic length and average Berry curvature is present, such that 2π 2 B is equal to the area of an effective unit cell where translation operators commute.
Applying this to TBG, we replace 2π 2 B with the area of the moiré cell: where a = a 1,2 is the moiré lattice constant. We now apply this dictionary to obtain TBG on a torus by suitably modifying each step of Ref. 34 . The parallelogram that defines our torus will be spanned by the vectors L 1 = N 1 a 1 and L 2 = N 2 a 2 , where N 1 and N 2 are positive integers so that there is an integer number of unit cells in each direction. This is the equivalent of magnetic flux quantization through a torus; indeed under the replacement (26) the flux quan- We specify generic twisted boundary conditions ψ(r + N 1 a 1 ) = e iφ1+iπN1 ψ(r), ψ(r + N 2 a 2 ) = e iφ2 ψ(r), up to a relative phase factor between the layers that we ignore (and can be removed with a gauge transformation). The factor e iπN1 is unimportant, and can be absorbed into φ 1 , but it will conveniently absorb the sign (−1) N1 (11) from z → z + N 1 a 1 in the upcoming discussion. For the wavefunctions (9), the twisted boundary conditions imply where we used the behavior of ϑ 1 ((z − z 0 )/a 1 |ω) under translations (11) 43 . These conditions exactly match those in Ref. 34 provided we identify N s = N 1 N 2 and τ = N 2 ω/N 1 . This is because our 1/ϑ 1 factor has similar properties under translations to the factor e −y 2 /2 2 B that appears in the QHE, a fact that is exploited further and made more precise in the next section. Indeed, everything from here onward is the same as in Ref. 34 if we use these identifications. In particular, we write a Jastrow ansatz for the many-body wavefunction of N e electrons where Z = i z i . The boundary conditions (28) imply 34 for some z-independent constants η 1,2 . The integral of d dz ln(g(z)) around the parallelogram defining the moiré torus is 2πi times the number of zeros of g(z). The boundary conditions above then imply that the number of zeros is m = N 1 N 2 /N e , which must then be an integer. If we demand that all the zeros are at the same point, we obtain the Laughlin factor g(z) = (ϑ 1 (z/L 1 |τ )) m .
We must have m odd for fermion antisymmetry to hold. The function (31) vanishes as z m when two particles approach each other. This ensures that this wavefunction will be a zero energy eigenstate of any "pseudopotential" where m < m and is a delta function that respects the periodicity of the torus. In particular, because all bandprojected states are smeared on the scale of the moiré length, a Coulomb interaction screened at distances much shorter than the moiré length will be well approximated as V (r) ≈ v 0 δ(r) + v 1 ∇ 2 δ(r), and (31) suffers no energy cost from this potential for m = 3 (though we cannot identify the coefficients with matrix elements of relative angular momentum eigenstates).
We note that angular momentum conservation or magnetic translation symmetry are not necessary for these arguments; they simply make the discussion on a plane or sphere more streamlined as symmetry arguments alone dictate that one only needs to care about Haldane pseudopotentials 35 . Matrix elements here in contrast will have dependence on the center of mass coordinate (due to lack of magnetic translation symmetry), but this should not matter as for sufficiently short range potentials this dependence will be dwarfed by the gap between v 1 and v 3 .
As in Ref. 34, one may proceed to derive boundary conditions for the center of mass piece of the wavefunction From these expressions, one finds that there are m linearly independent solutions which characterizes the m-fold degeneracy of Laughlin-like FCI states on the torus. Furthermore, these solutions braid into eachother upon adiabatically inserting a flux quantum through slowly increasing φ 1 or φ 2 . One may go further and discuss quasihole excitations as well, as in Ref. 34.

V. FLAT TBG BAND AS A GENERALIZED LLL
In this section, we interpret the previous results by showing that the chiral-MATBG wavefunctions (9) also arise in a generalized LLL where the magnetic field is inhomogeneous.
To do this we note that only the magnitude of ψ K (r)/ϑ 1 ((z − z 0 )/a 1 |ω) enters in bandprojected matrix elements of density-density operators: the position dependent layer polarization and phase has no effect. We will find generalized LLL wavefunctions that match the chiral-MATBG wavefunctions (9), except with ψ K (r)/ϑ 1 ((z − z 0 )/a 1 |ω) replaced by its magnitude.
Generally an inhomogeneous magnetic field leads to Landau level dispersion and wavefunctions that are not holomorphic. There are two situations that we found where this is not the case and states in the LLL can be defined as zero modes of the antiholomorphic canonical momentum Π = (Π x + iΠ y )/2 where Π = p − eA and ∇ × A = B(r).
One situation is that of a Dirac particle with the Hamiltonian The zeroth Landau level is then sublattice polarized and a zero mode of Π. The other situation is that of a nonrelativistic particle in a Kähler geometry. Here, one has both a curved space and an inhomogeneous magnetic field, such that the magnetic field is proportional to the Kähler form [46][47][48][49] . However, this case comes with the additional stipulation that B cannot vary in sign (else √ det g would), and so the Dirac model is more general. We will see that the inhomogeneous B(r) of interest to us does vary in sign within the moiré cell.
We now need to discuss the zero modes of Π = −i∂−eĀ. We takeĀ = iB 0∂ K/4 for some K(r). The magnetic field is then B(r) = B 0 ∂∂K. The states in the generalized LLL are where B0 = /eB 0 . We will match (38) with (9) for a suitable K that we now define. We set 2π 2 B0 = a 2 √ 3/2 so that there is one flux quantum per moiré cell, as in (26). We also decompose K = K 0 + K λ so that FIG. 2. Effective magnetic field B(r) such that TBG states are LLL states. The magnetic field has moiré-periodicity stemming from its relation to TBG wavefunctions (42) .
satisfies ∂∂K 0 = 1 and K λ generates mean-zero spacial variation B(r) − B 0 = B 0 ∂∂K λ . The numerical shift in the Landau gauge for K 0 was chosen to match the quasiperiodicity of the magnitude of the wavefunctions (9) under transla- For later convenience we define The function |ψ 0 | 2 is positive and moiré-periodic. We will see that |ψ 0 | 2 represents an "ideal" version of |ψ K | 2 . The K λ term generates the intra-unit-cell variations of B and breaks continuous translation symmetry down to translations by moiré vectors. It is given by such that the inhomogeneities are driven by ψ K = ψ 0 . Substituting (39), (40), and (41) into (38) one obtains the TBG wavefunctions (9) up to the previously discussed difference in phase and layer polarization. The inhomogeneous magnetic field is It is plotted in Fig. 2. We note that while B(r) is strongly inhomogeneous and changes sign within a moiré cell, much of the variation is averaged over in matrix elements as all band projected states are delocalized over areas ∼ a 2 . These results explain the success of the arguments in the previous section, where we constructed a Laughlin state on a torus using Haldane's arguments in Ref. 34. Indeed, one would expect this for any generalized LLL where the magnetic field is periodic, as the only difference in the wavefunctions from the ordinary LLL is a periodic factor e −K λ /2 2 B 0 . This extra periodic factor does not change anything about the boundary conditions on the torus and so all arguments go through unaffected.
We may relate the above approach to the momentum space geometry described in section III by carrying out this process in reverse. I.e., one may consider any periodic magnetic field B(r) and the corresponding LLL wavefunctions (38) with K 0 given by (39) and ∂∂K λ = (B(r) − B 0 )/B 0 so that |ψ K | 2 = |ψ 0 | 2 exp −K λ /2 2 B0 . One may then use Bloch periodic states similar to (12), with the overall position dependent phase and layer polarization removed. Because the periodic states are holomorphic in k x + ik y the condition (22) is exactly satisfied. The Berry curvature will in general be inhomogeneous, as it is for chiral-MATBG.
It is interesting to ask whether every FCI model that satisfies (22) can be described this way. Indeed, there has been a lot of work on curved-space Landau levels including both firstquantized real-space Laughlin approaches and effective field theories, for example Refs. 46-49. Despite their success in the FQHE, these techniques so far have seen little use in the FCI literature due to the discreteness of the lattice. We argue that our work may change this state of affairs for models satisfying (22); other models may be approached from these special models through adiabatic continuity.

VI. DEVIATIONS FROM CHIRAL
MODEL: κ > 0 In this section we explore deviations from the chiral model: κ > 0. The parameter κ in (2) describes the relative strength between tunneling at AA sites and AB/BA sites, such that if κ = 0 AA tunneling is absent. Lattice relaxation effects reduce AA tunneling and increase AB/BA tunneling so that κ is smaller than the value of 1 which was assumed in the original Bistritzer-Macdonald model 41 . Refs. 36 and 37 estimates κ ≈ 0.7 − 0.8. While the approximation κ = 0 may seem drastic, the sublatticepolarized flatband wavefunctions for the realistic model are adiabatically connected to those of the chiral model 20 Combined with the fact that the topological order of FCI states is stable to deformations up to some critical value, this suggests that the chiral model has the potential to describe the physics of FCIs described here.
We proceed to investigate the effects of κ > 0 on band energetics and band geometry. The results are plotted in Fig.3. We quantify band energetics through inverse flatness ratios R = Bandwidth/Bandgap to avoid divergences at κ = 0. There are bands above and below the band hosting the FCI state (the lower of the two flat bands) and so we have two different ratios: R upper is the ratio of the bandwidth to the band gap to the upper band (upper flat band) and R lower is the ratio of the bandwidth to the gap to the lower band (non-flat band). Both ratios increase with κ. Larger values of δ have the adverse effect of lowering the band gap to the lower band so that R lower increases. The band gap to the upper band increases with δ, though the bandwidth does as well such that for δ 0.05 the ratio R upper is nearly δ-independent.
We do not consider points in parameter space where the bandgap to the lower band is less than 0.03E 0 , as the near-crossing poses numerical issues and these points in parameter space are likely not suitable for FCI phases anyway.
We now discuss the influence of κ > 0 on the wavefunctions. We plot the distance from the chiral model wavefunctions, quantified by The distance increases with κ but barely depends on δ. We also track the Berry curvature inhomogeneity F . For fixed δ, the Berry curvature distribution becomes less uniform as κ is increased and begins to diverge as the band gap to the lower band closes. This is illustrated in Fig. 4 showing that the Berry curvature gets more concentrated at the Γ point with increasing κ. As δ increases F does too, as expected since increasing δ decreases the lower band gap. As κ increases T does too. Similar to the distance (43) it is nearly independent of δ. One conclusion we arrive at is that a small but not too small value δ ≈ 0.05 seems to be optimal for FCI phases. This δ is large enough to avoid large R upper , and keeps R lower and F smaller than if δ were larger. Sublattice potentials from hBN are somewhat higher: ∆ ≈ 17 − 30 meV, depending on alignment angle, translates to δ = 0.10 − 0.18, where we used α = 0.568 for the first magic angle 33 and w = 100 meV. It is also worth noting that we kept α constant and equal to the chiral magic angle value, whereas it is likely advantageous to tune it as a function of κ and δ. For example, lowering α by 0.02 at κ = 0.73 and δ = 0.1 lowers Berry curvature fluctations by more than 20%. We will expand on the effects of changing α in the future.
In the presence of Coulomb interaction, interaction-driven band renormalization effects coming from normal ordering terms as well as the effect of remote bands become important 5,20,44,50,51 . Due to the small "bare" Upper Band < l a t e x i t s h a 1 _ b a s e 6 4 = " c H q R I C a e n T b z f E 3 T 7 w T x i H q 5 L Y k = " > A A A C N n i c b V B N S 8 Q w E E 3 9 d v 1 a 9 e g l u A g e p L Q q q L d F L 1 6 E F V w V t s u S p t M 1 m K Y l m Y p L 6 a / y 4 u / w t h c P i n j 1 J 5 i t e / B r S I b H m z e Z y Q s z K Q x 6 3 t C Z m J y a n p m d m 6 8 t L C 4 t r 9 R X 1 y 5 N m m s O b Z 7 K V F + H z I A U C t o o U M J 1 p o E l o Y S r 8 P Z k V L + 6 A 2 1 E q i 5 w k E E 3 Y X 0 l Y s E Z W q p X P w s i i G 1 v 9 V K h I Y r K Q v f D s v D c o x 3 q u X 6 V y l q A c I / f R Y F O i n a W g Q 5 2 q D 3 H T E V l 2 a s 3 P N e r g v 4 F / h g 0 y D h a v f p T E K U 8 T 0 A h l 8 y Y j u 9 l 2 C 2 Y R s E l 2 L G 5 g Y z x W 9 a H j o W K J W C 6 R b V H S b c s E 9 E 4 1 f Y q p B X 7 v a N g i T G D J L T K h O G N + V 0 b k f / V O j n G h 9 1 C q C x H U P x r U J x L i i k d e U g j o Y G j H F j A u B Z 2 V 8 p v m G Y c r d M 1 a 4 L / + 8 t / w e W u 6 + + 5 u + f 7 j e b x 2 I 4 5 s k E 2 y T b x y Q F p k l P S I m 3 C y Q M Z k h f y 6 j w 6 z 8 6 b 8 / 4 l n X D G P e v k R z g f n x P y q 4 M = < / l a t e x i t >  h P T h s f q 0 U 4 S F k g y V D z i l F g n t X t a Z u H F t F + u e F V v D r x K / J x U I E e j X / 7 q D W K a S q Y s F c S Y r u 8 l N s i I t p w K N i 3 1 U s M S Q s d k y L q O K i K Z C b L 5 u V N 8 5 p Q B j m L t S l k 8 V 3 9 P Z E Q a M 5 G h 6 5 T E j s y y N x P / 8 7 q p j W 6 C j K s k t U z R x a I o F d j G e P Y 7 H n D N q B U T R w j V 3   (22) T are plotted in (b-f) respectively. All become less ideal as κ increases, though overall less so for sublattice potentials small but not too small δ ≈ 0.05. Points for which the lower band gap is less than 0.03E0 are excluded. dispersion, these effects alter the band dispersion significantly while keeping the flatband wavefunctions largely unchanged. Thus, we expect our conclusions in regard to the geometry of the flat band wavefunctions to be robust to such effects. On the other hand, we expect band dispersion to be strongly renormalized 20,44 . This is also supported by STM data 52-55 which find a separation of around 10-15 meV between the van Hove peaks at empty filling which is significantly larger than the expected 3-5 meV based on the non-interacting band structure. We note, however, that despite the substantial bandwidth renormalization, it still remains a factor of 2-3 smaller than the main interaction scale 20,44 . This is also supported by STM data [52][53][54][55] where the interaction-induced separation between the van Hove peaks at charge neutrality is significantly larger ∼35-50 meV. This suggests that the realization of an FCI is feasible in the realistic setting even when intraction induced bandwidth renormalization is taken into account.
We also note that the interaction potential is a highly tunable quantity in a continuum model. This has been beautifully illustrated in two recent experimental works 38,39 which managed to tuned the interaction potential by changing the distance to the metallic gate which controls the screening of the interaction. The strong dependence of the obtained phase diagram 38,39 on the gate distance suggests the sensitivity of the different interaction parameters, including band renormalization effects, to screening. This makes it reasonable to expect the existence of some parameter regime where the FCI phases discussed here are stable.

VII. CONCLUSION
Our results provide a good starting point for analyzing the possibility for FCI states in TBG through numerical tests of adiabatic continuity. Indeed, one may start from LLL wavefunctions where the FQHE is well established, turn on the inhomogeneous part of the magnetic field until one reaches the chiral model wavefunctions, and then gradually turn on κ until one reaches realistic wavefunctions. If the many body gap does not close during this process, one has established adiabatic continuity. One may also explore the continuous space of paths, parameterized by δ and the strength and range of the interaction potential, and see if any lead to realizable FCI states.
Independent of relevance to current TBG ex-periments, our model is a promising idealized model for numerical studies of FCI stability. Previous studies have learned much from toy lattice models 30 and Harper-Hofstadter models 56 but were restricted by lack of tunable parameters that produced distinct effects on band energetics and geometry. Twisted bilayer graphene is highly tunable; it may be modified through (small) changes in angle, different strength sublattice splittings on each layer, a perpendicular electric field, pressure, screening length, and more. It can add many more data points to these studies. Note Added: After completion of this work, reference 57 appeared, which contains a numerical study of FCIs in MATBG . This is in agreement with our analytical results regarding the stability of Laughlin FCIs.
Acknowledgements: This work was partly supported by the Simons Collaboration on Ultra-Quantum Matter, which is a grant from the Simons Foundation (651440, A. V.) and A. V., E. K., were supported by a Simons Investigator grant. E. K. was supported by the German National Academy of Sciences Leopoldina through grant LPDS 2018-02 Leopoldina fellowship. P.J.L. was partly supported by nsfdmr 1411343 and the Department of Defense (DoD) through the National Defense Science & Engineering Graduate Fellowship (NDSEG) Program. G.T. acknowledges support from the MURI grant W911NF-14-1-0003 from ARO and by DOE grant de-sc0007870.