Electron-hole asymmetry and band gaps of commensurate double moire patterns in twisted bilayer graphene on hexagonal boron nitride

Spontaneous orbital magnetism observed in twisted bilayer graphene (tBG) on nearly aligned hexagonal boron nitride (BN) substrate builds on top of the electronic structure resulting from combined G/G and G/BN double moire interfaces. Here we show that tBG/BN commensurate double moire patterns can be classified into two types, each favoring the narrowing of either the conduction or valence bands on average, and obtain the evolution of the bands as a function of the interlayer sliding vectors and electric fields. Finite valley Chern numbers $\pm 1$ are found in a wide range of parameter space when the moire bands are isolated through gaps, while the local density of states associated to the flat bands are weakly affected by the BN substrate invariably concentrating around the AA-stacked regions of tBG. We illustrate the impact of the BN substrate for a particularly pronounced electron-hole asymmetric band structure by calculating the optical conductivities of twisted bilayer graphene near the magic angle as a function of carrier density. The band structures corresponding to other $N$-multiple commensurate moire period ratios indicate it is possible to achieve narrow width $W \lesssim 10$ meV isolated folded band bundles for tBG angles $\theta \lesssim 1^{\circ}$.


I. INTRODUCTION
The electronic properties of magic angle twisted bilayer graphene (tBG) and other graphene moire systems have become an intense focus of research in recent years after observations of strong correlations and superconductivity in transport experiments 1-3,5-7 . One of the striking observations revealing the inherent topological character of the bands in graphene moire systems 8,9 is the spontaneous anomalous Hall effects at 3/4 filling in tBG that appear even in the absence of an external magnetic field [10][11][12] , which can be viewed as a solid state realization of the Haldane-like model in a honeycomb lattice 13 . Theories attempting to explain the spontaneous Hall effects relied on band structures that included a small mass term in the Dirac Hamiltonian contacting the BN substrate 14,15 that give rise to precursor bands of the spontaneous quantum Hall phases. The fact that this effect in tBG has been observed so far only in twisted bilayer graphene for samples that are nearly aligned with hexagonal boron nitride (tBG/BN) [16][17][18] motivates further research on the effects of the BN substrate in modifying the electronic structure of twisted bilayer graphene.
In this manuscript we show that commensurate double moire patterns of tBG on BN can have significant electron-hole asymmetric bands in addition to gaps between the moire bands which provide the conditions for the formation of valley Chern bands prone to strong correlations. Here we focus our discussions on commensurate double moire patterns of tBG/BN whose continuum Hamiltonian can be conveniently formulated in terms in a unique moire Brillouin zone. The electronic structure of double moire patterns were addressed in the past for arbitrary relative twist angles in systems like BN encapsulated graphene [19][20][21] or twisted trilayer graphene 22 giving rise the super-moire patterns. Earlier studies of BN/G/BN double moire systems 19 have shown that the largest double moire interference effects are seen when the patterns are commensurate and in these cases the electronic structure undergo significant modifications when interayer sliding is introduced.
Electronic structure calculations for tBG/BN relying on commensurate double moire patterns carried out recently have shown the opening of gaps between moire bands that approximately maintains electron-hole symmetry 23,24 .
In contrast, in our work we show that the moire bands can have significant electron-hole asymmetry depending on the relative sliding between the moire patterns and twist angles. For type-1 double moire systems we can obtain on average a narrow conduction and a wider valence band, while we predict that type-2 double moire systems can form ordered phases for both conduction and valence bands slightly favoring the valence bands. Our results provides a possible justification that the spontaneously broken symmetry phases have been observed in experiments only for the conduction bands for the type-1 of double moire arrangements 7,10 , while correlated phases are likely for the valence bands for type-2. In both cases we find that isolation of the flat bands required to endow finite valley Chern numbers is facilitated by the opening of primary and secondary gaps. We find that valley Chern numbers are finite with values ±1 over a large parameter space in keeping with experimentally observed quantized Hall conductivities, while the local density of states associated to the flat bands invariably concentrate around the AAstacked regions of tBG regardless of the BN arrangement. We illustrate for a particularly pronounced electron-hole asymmetric band structure the impact of the BN moire pattern by showing the differences in the optical conductivities between tBG and tBG/BN as a function of carrier density. Band structures corresponding to other rational number N -multiple commensurate moire periods indicates that it is possible to find narrow width W 10 meV isolated folded band bundles for tBG angles θ 1 • .
The manuscript is structured as follows. In Sec. II we introduce the model Hamiltonian for the doubly commensurate tBG on BN. In Sec. III we present the exact geometric conditions under which the double moire patterns become commensurate with the same period and orientations. We then move on to analyze the electronic structures in Sec. IV where we present the band structures for different relative sliding vectors paying particular attention to the bandwidth, gaps, and valley Chern numbers, and study the effects of interlayer potential differences due to a perpendicular electric field. In Sec. V we illustrate the effects of the BN substrate in the optical conductivy near the magic angle twisted bilayer graphene for a clearly electron-hole assymetric electronic structure with narrow conduction band. In Sec. VI we discuss the electronic structures for commensurate double moire patterns with different periods. Finally we close the paper in Sec. VII with the conclusions.
For tBG/BN, the moire potential Hamiltonian H M BN (r) has the form where Here, ϕ m are the polar angles of six moire reciprocal lattice vectorsG m of the tGBN interface with a prime symbol to distinguish from the unprimed moire reciprocal lattice vectorsG m for tBG. Here we use the moire potential parameter set given by C AA = −14.88 meV, φ AA = 50.19 • , C BB = 12.09 meV, φ BB = −46.64 • , C AB = C BA = 11.34 meV, φ AB = 19.60 • following Refs. 26 and 30. For tBG, the reciprocal lattice vectors of BZ of graphene are G 1 = 2π(0, 2/ √ 3a G ), G m =R π(m−1)/3 G 1 and the moire reciprocal lattice vectors are given bỹ G i = −θẑ × G i (i = 1, · · · , 6). On the other hand, for Moire patterns and the corresponding energy bands in type-1 for (b) AA-tBG/BN, (c) AB-tBG/BN, (d) BA-tBG/BN, and in type-2 for (g) AA-tBG/BN, (h) AB-tBG/BN, (i) BA-tBG/BN. The first electron and hole bands are emphasized by red and blue colors, respectively. The distinct difference between type-1 and type-2 is in the direction of twist angle θ between the bottom G layer and the hBN layer, which is positive for type-1 and negative for type-2. Underneath each moire pattern figure we illustrate three commensurate stackings with colored circles, blue and green for the hBN layer, red and tender red for the bottom G layer, and tender violet and tender blue for the top G layer. The labeling of stacking configuration starts from the lowest layer (hBN) and the sublattice on the left within a layer. Energy bands were plotted along the high-symmetry lines in (e) for type-1 and (j) for type-2. For type-1, mBZ of tBG and tGBN have an angle difference φ − φ = ∆φ = −120 • , while ∆φ = −60 • for type-2 (see Appendix A for more details). tGBN we haveG i = G i − θẑ × G i . We assume that for the commensurate angle sets (see the next section for further details), the mBZ of tBG and tGBN have the same size and share their two corners in momentum-space as shown in Fig. 1 (e), (j) (see Appendix A for further details).
We label each commensurate stacking in real-space with three characters, A, B, and C to label from bottom to top layers the relative position of the left sublattice of the two atoms unit cell of G or BN. For instance, when the local AA-stacking region in tBG is AA stacked on top of hBN (AA-tBG/BN), a type-1 system gives rise to three AAA, ABC, and ACB local stacking at the symmetry points, while AAA, ACA, and ABA local stacking appear in type-2 systems. We plot each of the commensurate stacking configurations of tBG/BN with six small colored circles under the moire patterns for type-1 in Fig. 1 (b)-(d) and for type-2 in (g)-(i).

III. COMMENSURATE TWIST ANGLE SETS
Among the many possible commensurate double moire patterns [19][20][21] here we want to focus on the subset that can be obtained by requiring two constraints: (1) moire length: the moire lengths of tBG and tGBN should be the same, and (2) moire angle: the difference of moire angle should be 0 (mod 60 • ). The moire length of tBG is and the moire length of tGBN is where α = a G /a BN = 1+ and = (a G −a BN )/a BN 26 are parameters that represent the lattice mismatch between two layers. We use for the lattice constant of graphene a G = 2.461Å and boron nitride a BN = 2.504Å, which results in = −0.017172. Here, we use the notation for angles measured with respect to the bottom graphene (G) layer for convenience since both moire interfaces are sharing it; we define θ = θ tBG as the twist angle of tBG which is measured with respect to the bottom G layer and θ = −θ tGBN is the twist angle of BN measured with respect to the contacting bottom G layer at the tGBN interface and hence the minus sign. The condition for coincident moire lengths L M tBG = L M tGBN requires the following relation between θ and θ θ = cos −1 (1 + α)(1 − α) + 2α cos θ 2 .
Accordingly, the rotation angle φ = φ tBG of the moire pattern in tBG with respect to the bottom G layer is 26 This expression reduces to φ = −sgn(θ) π 2 + θ 2 (12) indicating that the moire pattern orientation rotates with respect to the bottom G by one half of the total twist angle between the graphene layers. On the other hand, the tGBN moire pattern rotation angle φ = φ tGBN measured with respect to the bottom G layer is given by The condition for the two equal period moire patterns to be commensurate is where n is an integer number 19 . There is no analytical solution for n = 0 but we have for n = 1 the type-1 solution at ∆φ = 60 • or equivalently, −120 • because of π periodicity of the tangent function, and have for n = 2 the type-2 solution at ∆φ = −60 • or equivalently 120 • . The rotation of the moire patterns and associated moire Brillouin zones are shown graphically in Fig (8,9) leads to 4 sin 2 (θ/2) = (α − cos(θ )) 2 + sin 2 (θ ) (15) that when we use the θ = 2θ condition 4,23 is equivalent to Eq. (14) with n = 1, that yields same angle signs for type-1 solutions (θ, θ ) = (1.1329 • , 0.5664 • ). The equations for unequal sign angles for type-2 solutions can be formulated by combining both Eqs. (14,15) for n = 2 and their numerical solutions are (θ, θ ) = (1.1463 • , −0.5932 • ). Geometric discussions of the double moire commensuration conditions in real and momentum space are discussed in Appendix A. We show the real-space moire patterns in three different stackings between the top and bottom G layers for both type-1 in Fig. 1 (b-d) and type-2 in Fig. 1 (g-i) where we can clearly see the differences in both real-space patterns and momentum space band structures.

IV. BAND STRUCTURES
In this section we present the band structures of tBG/BN corresponding to the two types of commensurate stacking arrangements introduced in earlier sections as a function of moire pattern sliding and perpendicular external fields. We will be paying attention to the bandwidth, gaps, and the electron-hole asymmetry. Knowledge of the band structure allows to estimate regimes of strong correlations that satisfy U eff /W 1 where the bandwidths W are sufficiently narrow versus the effective Coulomb interaction where we used the dielectric constant of graphene as r = 4, the moire length l M = a G /(2 sin(θ/2)) a G /θ and the screening length λ D = 2 0 /e 2 D(δ p , δ s ) that depends on the two-dimensional density of states D(δ p , δ s ) = 4(|δ p |u(−δ p ) + |δ s |u(−δ s ))/(W 2 A M ) defined in terms of the Heaviside step function u(x) that includes screening due to overlap between neighboring bands. 28 The valley Chern number of the n th energy band C n = mBZ d 2 k Ω n (k)/(2π) is obtained integrating the Berry curvature 31 given by Ω n (k) = −2 n =n Im n| ∂H ∂kx |n n | ∂H ∂ky |n /(E n − E n ) 2 . Additional information on the electronic structure can be found in appendix C where we present the normalized local density of states (LDOS) whose maxima are expected at the AA stacking sites of tBG.

A. Type-1 and type-2 electronic structures
We can distinguish two arrangements of double moire patterns that we classify as type-1 and type-2 depending on the relative twist of the layers in the tBG and tGBN substrate that in turn impacts the relative rotation of the associated moire patterns ∆φ = φ − φ. For double moire patterns of type-1 we require a relative rotation of ∆φ = −120 • or equivalently ∆φ = 60 • between the mBZ of tBG and tGBN as shown in Fig. 1(e). On the other hand, for type-2 we require ∆φ = −60 • (see Appendix A) or equivalently ∆φ = 120 • as shown in Fig. 1(j). The impact in the band structure of tBG/BN is represented along the high-symmetry line of mBZ of tBG (K −K −Γ −Γ −K), we show the band structures for in Fig. 1(b-d) for type-1 and in Fig. 1 (g-i) for type-2 for different displacement of the top G layer for an electronic structure model with unequal interlayer tunneling value that incorporates an effective vertical relaxation 28,29 . In our analysis the twist angles of tBG compatible with the double commensuration defined by Eqs. (13)- (14) essentially depends on the value of the lattice constant mismatch α = a G /a BN . The double commensuration angles θ 1.13 • and θ 1.14 • result from our choice of graphene and BN lattice constants. We can expect the bandwidths will remain narrow provided that these angles remain close to the experimental 32 θ exp MA 1.08 • and theoretical 25,28,29 θ th MA 1.05 • values where flat bands emerge. We estimated the changes in the electronic structure for twist angles θ of tBG away from the magic angle by using an approximation scheme that preserves commensuration with the hBN moire pattern potential which shows a similar qualitative behavior in the vicinity of the magic twist angle for variations δθ on the order of ∼ 10% or smaller, see appendix C Fig. A3.
The electronic structures depend strongly on the relative sliding between the moire patterns of tBG and tGBN. We use a vector τ = (τ x , τ y ) to indicate the sliding of the top G layer prior to rotation and with this notation the local AA-, AB-, and BA-stackings at the origin correspond to (τ x , τ y ) = (0, 0), a G (0, 1/ √ 3), and a G (0, 2/ √ 3). The local stacking maps between the two graphene sheets and BN subsrate are traced by assuming that the BN substrate is locked at the AA stacking configuration with respect to the bottom graphene layer. For instance, the AA-tBG/BN case where tBG is AA stacked prior to twisting corresponds to Fig. 1(b) and (g). For type-1 we can find a small gap of ∼ 2.5 meV atK while for type-2 the lowest energy bands are quasi-flat with a large gap of ∼ 35 meV at K and ∼ 20 meV atΓ. For AB-tBG/BN corresponding to Fig. 1(c) and (h), type-1 bands consist of mutually overlapping bands with large secondary gaps δ s for both conduction and valence bands that separate the lowest energy bands with the next band that is father in energy from charge neutrality, whereas for type-2 we find isolated low energy conduction and valence bands with a δ p gap of size ∼ 15 meV atΓ. For BA-tBG/BN corresponding to Fig. 1 (d), we see a large electronhole asymmetric electronic structure separated with a small primary gap δ p 3.5 meV consisting of a nearly flat low energy conduction band with W ∼ 5 meV and wider valence band W ∼ 30 meV that reminds of the experimental observations of spontaneous Hall effects for electrons but not for holes 7,10 . Conversely for type-2 we have wider conduction and narrower valence bands that are separated by a finite gap of δ p 5.6 meV for the entire mBZ.
For both types of solutions, the secondary gap δ s for the conduction and valence bands remain open for the three AA, AB, BA symmetric stacking geometries, and we will show shortly that this is the case for the entire range of sliding geometries. This is in contrast to the behavior of the primary gap δ p that clearly closes when the top G layer is placed at (τ x , τ y ) = a G (0, 1/ √ 3) in a type-1 system.
The Berry curvatures for discrete three commensurate slidings of the top G layer for conduction and valence band and the normalized local density of statesD(r, E) = D(r, E)/max(|D(r, E)|) at the van Hove singularity (vHs) can be found in appendix C Fig. A4 that illustrates the impact of the BN substrate in altering those quantities.
Having identified the distinct behavior of the bands for type-1 and type-2 moire pattern arrangements at three different symmetric stacking configurations we discuss in the following the dependence of the band structure to other continuous stacking sliding geometries and interlayer potential differences. Here we extend our analysis of the lowest energy bands to other continuous stacking geometries as a function of τ = (τ x , τ y ) and explore the changes in the electronic structure due to an interlayer potential difference η introduced by a perperdicular electric field. For all possible local stacking geometries we have obtained the maxima, minima, average and the standard deviations of the bandwidths, which we succinctly summarize in Table I. From this data we can observe a more pronounced electron-hole asymmetry for type-1 systems with narrower conduction bands on average and smaller bandwidth fluctuations that posits the preference of the conduction bands to form ordered phases with respect to the valence bands. On the other hand, we expect a higher likelihood of observing ordered phases for both conduction and valence bands for type-2 systems because the electron-hole asymmetry is generally less pronounced, while narrower valence bands with smaller fluctuations are found for the average results.
The phase diagram maps for the bandwidth W , the U eff /W ratios, the δ p , δ s gaps, and valley Chern number maps for continuous stacking geometries are shown in Fig. 2. The stacking configurations with the narrowest bandwidth and largest gaps have greater chances of enhanced Coulomb interactions that can lead to ordered phases. For roughly half of the stacking configuration space we can observe favorable U eff /W 1 condition and they largely coincide with the regions that have finite valued valley Chern numbers. The regions with suppressed U eff /W are mainly due to the closure of the primary gap δ p that overlaps the two low energy flat bands and enhances the Coulomb screening.
The phase diagram maps for the bandwidth W , the gaps and resulting U eff /W ratios in Fig. 3 shows that a perpendicular electric field directed in the appropriate sense can enhance the narrowing of the bands and has an overall effect of enhancing the secondary gap δ s , while the opening and closure of the primary gap δ p follows a more complex behavior as a function of stacking and applied electric fields.  3. (Color online) Bandwidth W, secondary δs and primary δp gaps, ratio of the effective Coulomb interaction to W for type-1 (upper row) and type-2 (lower row) as a function of continuous sliding of the top G layer in the y-direction (τy) prior rotation obtained for various values of the interlayer potential difference η. We see that both stacking configuration and electric fields can substantially alter the bandwidth and band isolation of the low energy moire bands.

C. Band gap analysis atK point
As discussed earlier, distinct electronic low energy bands are expected depending on the moire pattern stacking types that we classify as type-1 and type-2. For example, in type-1 of AA-tBG/BN there is a small gap ∼ 2.5 meV atK in contrast to the large gap ∼ 35 meV in type-2. A band gap near charge neutrality in tBG/BN double moire systems can be expected to appear due to moire pattern interference effects. We illustrate this behavior by obtaining the analytical expressions for the size of the gaps atK for different moire pattern types and relative stackings arrangements as a function of the model parameters that define the moire patterns. For this purpose we consider a truncated 8 × 8 Hamiltonian as a simple model to describe tBG system that captures the scattering of the electrons from the three Dirac cones of the bottom layer to one Dirac cone on the top layer 25 . AtK the eigenstates can be represented by Ψ = (ψ 1 , ψ 2 , ψ 3 , ψ 4 ) T consisting of four two-component sublattice pseudospin spinors 33 . Here the ψ 2 , ψ 3 spinors can be rewritten in terms of ψ 4 by multiplying with a unitary matrix ψ 2 = U ψ 4 and ψ 3 = U † ψ 4 , and thus allowing to further reduce our Hamiltonian to a 3×3 matrix to obtain a simpler characteristic equation. Then, the expression for the band gap atK is, where Q = υ F θk D . This expression agrees up to 2 significant digits with the gap size obtained by diagonalizing numerically the eight-band model. The above equation holds for both type-1 and type-2 solutions and for the three symmetric AA, AB, BA stackings but the analytical coefficients A 3 and B 3 are different for each case. Below we present the expression for the A 3 and B 3 coefficients for type-1 BA-tBG/BN for the electronic structure showing the largest electron-hole asymmetry that incorporates the effect of the interlayer potential difference η that can be induced by a perpendicular electric field and the moire pattern potentials generated by the BN substrate that affect the magnitude of the gap. For sufficiently small magnitudes of η the A 3 and B 3 coefficients in Eq. 18 can be expressed in terms of the moire pattern coefficients stemming from the hBN substrate as follows In Appendix B we show a more detailed derivation of the analytic expressions for other types of solutions and stackings.

V. OPTICAL CONDUCTIVITY
Electronic structure reconstruction in the energy ranges of a few tens of meV in the vicinity of the Dirac cone of graphene can be conveniently explored by means of mid-infrared and terahertz spectroscopy. Here we show the numerical results for the optical absorption in type-1 BA-tBG/BN as well as tBG without a BN substrate in an effort to estimate the differences in the optical signals that can be expected in these two systems. The real part of the longitudinal optical conductivity normalized by the universal optical conductivity of grapheneσ 0 = πe 2 /2h is given by [34][35][36][37][38] Re[σ xx (ω)]/σ 0 = 16 ω where J α = −∂H/∂k α is the current operator, f ( ) is the Fermi-Dirac distribution function. k,i is the ith eigenstate energy at k = (k x , k y ). In Fig. 4(a) we show the energy bands for type-1 BA-tBG/BN (black solid line) and tBG (red dotted line) without a BN substrate where we achieve double commensuration of equal period moire patterns by using the twist angles of (θ, θ ) = (1.13291 • , 0.5664 • ) as we introduced earlier.
The horizontal lines are the three different chemical potentials µ = −0.021 (magenta circle), 0.012 (blue triangle), and 0.035 eV (green square) considered, and we plot the real part of the linear optical absorption in the longitudinal direction for tBG in Fig. 4(b) and type-1 BA-tBG/BN in Fig. 4(c). The presence of the BN substrate leads to markedly different asymmetric optical absorptions for µ = −0.021 eV (magenta circle) and 0.035 eV(green square), and one can observe small gaps less than 10 meV for µ = −0.021 eV (magenta circle), 0.012 eV (blue triangle) in type-1 BA-tBG/BN as shown in Fig. 4(c), unlike without BN.

VI. COMMENSURATE DOUBLE MOIRE ANGLE SETS
So far we have considered doubly commensurate systems with equal moire pattern periods and here we extend our band structure analysis restricting to commensurate double moire systems whose periods are multiples of each other and have modulo 60 • alignment. The moire period L M tBG a G /θ for tBG and L M tGBN a G /(θ 2 + 2 ) 1/2 for tGBN as a function of twist angles θ, θ are given in Eq. (8) and (9). The set of twist angles that satisfy the commensuration conditions are listed in Table II for type-1 and type-2 systems when L M tGBN = N ×L M tBG for integer and rational values of N considering that longer moire periods can be achieved by using smaller twist angles.
Our analysis focuses on the L M tBG > L M tGBN cases when the moire periods of tBG are greater than those of tBG for N = 1/2, 1/3, and 1/4 ratios. These commensurate patterns can be achieved in the limit of small tBG twist angles θ 1 • and gives rise to band folding as the moire Brillouin zone reduces in size. The resulting band structures are presented in Fig. 5 for the three different starting stackings, AA, AB, and BA where we can observe bundles of narrow bands that we highlight with blue-shades that are separated from each other forming bands of moire bands. Likewise we have also highlighted with red shades the isolated nearly flat moire bands that are separated from each other through avoided gaps. In particular the band bundles we observe for N = 1/4, 1/3 commensurate ratios are as narrow as ∼ 10 meV and are comparable to the bandwidth of the lowest energy flat bands that are formed near the magic twist angle θ ∼ 1 • of N = 1.
These folded nearly flat bands are weakly dispersive in the reduced moire Brillouin zones in the limit of small tBG twist angles. We expect that the moire pattern features that appear at low energies for long moire periods due to the BN substrate can introduce gaps that separate the folded band bundles that enhance their isolation and reduce their screening. The notable changes in the band structures that we can observe for FIG. 5. Band structures for the twist angle sets of N -multiple moire periods of N = 1/4, 1/3, 1/2 for type-1 and type-2 systems listed in Table II. The AA, AB, and BA stacking labels indicate the starting stacking of the top G layer prior to rotation. The isolated folded band bundles with narrow width 10 meV are indicated by the blue-shaded regions, while the isolated single flat bands of width 10 meV are indicated by red-shaded regions.
different relative sliding configurations for the doubly commensurate geometries indicate that small deviations of the twist angles can in principle lead to significant changes to the electronic structures, while in most cases we can also expect the formation of gaps between the folded moire band bundles due to the BN substrate.

VII. CONCLUSION
We have investigated the electronic structure of commensurate double moire patterns of twisted bilayer graphene (tBG) on hexagonal boron nitride (BN) in an effort to understand the effects introduced by a BN substrate when it is brought to near alignment with the graphene layers, and identify the conditions that are favorable for the appearance of gapped moire bands with finite valley Chern numbers that ultimately lead to the anomalous Hall effects observed in experiments.
The band structure of tBG has been calculated based on a continuum model that effectively accounts for vertical interlayer relaxations by using unequal interlayer tunneling parameters. The effects of a second moire pattern produced by the BN substrate have been modeled assuming a moire pattern that is perfectly commensurate with tBG which conveniently allows us to study the physics of the system based on a single moire Brillouin zone continuum Hamiltonian model. For these doubly commensurate geometries we have highlighted the importance of aligning properly the moire reciprocal lattice upon rotation to properly construct the continuum moire Hamiltonian, and have illustrated the role of the BN moire patterns and interlayer potential differences in opening the primary band gaps by obtaining the analytical form of the gaps nearK for a truncated model.
When the periods of each moire patterns are equal, we have identified two types of commensurate double moire patterns depending on the relative twist angles of the tBG layers with respect to the BN substrate that we classified as type-1 and type-2 and whose band structures have revealed markedly different features in their electron-hole asymmetry, the formation of band gaps, and dependence to local interlayer sliding arrangement. In type-1 double moire systems we find narrower overall conduction bands than valence bands with the most prominent asymmetry appearing in the vicinity of the BA-tBG/BN stacking configuration where the conduction band has the smallest bandwidth ∼ 5 meV suggesting that conduction bands are more prone for broken symmetry phases than the valence bands, whereas the electron-hole symmetry is partially restored in type-2 solutions that slightly favor on average the narrowing of the valence bands suggesting that ordered phases will likely appear for both conduction and valence bands. The impact of the BN substrate on the electronic structure should be observable not only through transport measurement but also through the optical absorption probes as we illustrated for the specific BA-tBG/BN example versus tBG without a substrate.
We have extended our electronic structure study to other commensurate double moire systems with unequal moire pattern periods L M tBG > L M tGBN focusing on the smaller bilayer graphene twist angles and verified that they often lead to separated folded moire bands bundles with reduced widths ∼ 10 meV. Hence, our results suggest that tBG twist angles in the range of 0.3 • ∼ 0.5 • clearly below ∼ 1 • can also be candidate systems for exhibiting correlated phases.
The analysis presented in this work for commensurate double moire patterns for a variety of different sliding geometries shows how the BN substrate impacts the electronic structre of tBG, specifically leading to band gap openings that isolate the moire bands. Additionally, these band structures provide a useful starting point for understanding the electronic structure of incommensurate double moire tBG/BN systems that can be pictured as a coherent combination of the commensurate states for different stacking. A more detailed study of the incommensurate double moire patterns as they depart from perfect commensuration as well as the impact of lattice relaxations in the electronic structure is desirable to understand in greater detail the precursor states leading to the experimentally observed anomalous Hall orbital magnetization in tBG/BN.

VIII. ACKNOWLEDGMENTS
We gratefully acknowledge discussions with J. H. Sun for the derivation of the analytic band gaps. This work was supported by Samsung Science and Technology Foundation under project no. SSTF-BA1802-06 for J. S., the Korean National Research Foundation (NRF) grant NRF-2020R1A5A1016518 for Y. P., the Basic Science Research Program of the NRF-2018R1A6A1A06024977 for B. L. C., and by grant NRF-2020R1A2C3009142 for J. J. We acknowledge computational support from KISTI through grant KSC-2020-CRE-0072.

APPENDIX A: GEOMETRICAL DERIVATION OF COMMENSURATE TWIST ANGLES AND MOIRE PATTERN ROTATIONS
There are two classes of solutions for same period commensurate tBG/BN as shown in Fig. 1. Here we present a geometrical proof for both types of solution in momentum-space and a more detailed analysis about how type-1 and type-2 moire potentials in real-space lead to relative moire pattern rotation angles of ∆φ = −120 • and ∆φ = −60 • respectively.
We assume that the bottom G and top G layer are rotated by −θ/2 and θ/2, respectively. In Fig. A1 we plot the BZ of the bottom G layer with a thick blue solid line and the BZ of the top G layer with a thick green solid line, generating the mBZ of tBG which plotted by a black thin solid line. The thick brown solid line represents BZ of hBN. For tBG/BN to be commensurate, mBZ of tBG and of tGBN should have the same size and share their two corners. The geometrical understanding for type-1 is straightforward since the simple relation θ = 2θ holds as shown in Fig. A1(a). For type-2, on the other hand, if we assume the length from the center to the corner of BZ of hBN layer is lα, then the length from the center to the corner of BZ of G layer is l, and the two lengths indicated by the dashed purple line in Fig. A1 must be the same to share the two corners of the mBZ of tBG and tGBN. Then we get an equation which is equivalent to L M tBG = L M tGBN condition in Eq. (10). The shaded areas in red and blue in the upper row in Fig. A1(a) and (b) correspond to the mBZ depicted in Fig. 1(e) and (j), respectively.
We now move on to discuss the rotation of the moire patterns. The moire potential in Eq. (6) can be rewritten for the bottom G layer with a 2×2 matrix as follows. in tGBN for type-1 (upper two rows) and type-2 (lower two rows). The red arrow indicates the lattice vector of tGBN. The rest three columns represent the mass term and the real part of the in-plane gauge term in tBG/BN for type-1 (upper two rows) and type-2 (lower two rows) for three different ∆φ values 0 • , −60 • , and −120 • . The black arrows correspond to the lattice vector of tBG and the red arrows represent that of tBG/BN. Note that the moire patterns generated from tBG/BN coincide with those from tGBN (on the leftmost column) when ∆φ (the angle difference between red and black arrows) is −120 • for type-1 and ∆φ = −60 • for type-2.
whereG m are the moire reciprocal vectors of tGBN and APPENDIX B: BAND GAP ANALYSIS ATK Here derive the analytic expressions for the energy gap at a specific high-symmetry pointK that often is a useful estimate of the actual band gap at charge neutrality. We first consider truncating the range ofGvectors for an eight-bands model of tBG 25 possessing three Dirac Hamiltonians for the bottom layer and one Dirac Hamiltonian for the top layer. We use the fact that among the four two-component pseudospin eigenstates in the basis Ψ = (ψ 1 , ψ 2 , ψ 3 , ψ 4 ) T , 33 the two spinors ψ 2 and ψ 3 can be written atK in terms of ψ 4 by a unitary transformation as follows, for example, for AA-tBG/BN for conduction (valence) band in type-1 (type-2), and for valence (conduction) band in type-1 (type-2), we reduced our 8×8 Hamiltonian into 3×3 matrix to find a secular equation. In type-1 for the conduction band, and for the valence band, where the eigenstates of the conduction and valence band for eight-band model are given as We get characteristic equations from the above 3 × 3 matrices. For example, and where and (A12) for t = 1 (type-1) and t = 2 (type-2). For type-1, Θ and Ξ in the order of AA-, AB-, BA-tBG/BN are given as follows.
Also, A 3 and B 3 are given as follows, for type-1 AA-tBG/BN, for type-1 AB-tBG/BN, and and for type-1 BA-tBG/BN, and Since A 3 and B 3 are negligibly small, we can get three roots for each, E C = 0, E C = A 1 ± A 2 1 + A 2 , and E V = 0, E V = B 1 ± B 2 1 + B 2 , presuming A 3 = B 3 = 0. By retrieving A 3 , B 3 and plugging the lowest energies E C = 0, E V = 0, whose eigenvectors remain the same, into the original characteristic equations, we get a much more simplified expressions for the energy gap as written in Eq. (18). In the same manner, for other stackings of tBG the characteristic equations have the same form but the coefficients are different since the unitary matrix U connecting the spinor bases are different, resulting in different characteristic equations.
For AA-tBG/BN of type-2, the 3 × 3 matrix for a characteristic equation is given as where the eigenstates of the conduction and valence band for eight-bands model are given The A 3 and B 3 parameters for type-2 systems are given by for AA-tBG/BN, for AB-tBG/BN and The analytical values of the gaps in the presence of an interlayer potential difference η will be shown only for the case of BA-tBG/BN in type-1. The associated 3 × 3 matrix for the characteristic equation is for conduction band, and Here we supplement information on electronic structure presented in the main text for the bandwidth W , the secondary band gaps δ s , primary band gaps δ p , the ratio of the screened Coulomb interaction strength to the bandwidth U/W , and valley Chern numbers by considering the effects of twist angles under the assumptions of commensurate double moire patterns and also provide information of the local density of states at the van Hove singularities and the Berry curvatures of the nearly flat bands.
We begin by presenting in Fig. A2 through 1D plots various electronic structure information for variable τ y and fixed τ x = 0 that had been discussed in the main text. The topmost row (a) panel shows the fluctuations in the bandwidth for the low energy conduction and valence bands as a function of τ y where we can find a relatively narrower conduction band versus valence band for type-1 systems and more electron-hole symmetric bands for type-2 systems. Subsequent rows represent (b) the secondary band gaps and (c) the primary band gaps.
Information from (a)-(c) panels are used to determine (d) the U eff /W ratio where we can clearly observe its relative enhancement for conduction bands in type-1 but maintaining comparable magnitudes for electrons and holes for type-2 systems. The valley Chern numbers in (e) showing well quantized values or oscillating with noninteger values as a function of τ y are in keeping with the fact that the bands are isolated when we have positive primary and secondary gaps.
In Fig. A3 we provide additional information about how the electronic structure will be modified when we allow the twist angle θ of tBG to change from the specific values determined by the double moire commensuration condition. For this purpose we make the approximation that double commensuration is maintained for every value of θ which implies that the moire pattern period and angle also changes continuously together with the pattern of tBG. The dashed vertical line indicates the commensuration angle θ for tBG consistent with the lattice constants of graphene a G = 2.461Å and hexagonal boron nitride a BN = 2.504Å. We may expect that the double commensuration condition can still hold approximately in the vicinity of this vertical dashed line thanks to lattice relaxation effects that rotate and globally relax the lattices. The information in Fig. A3 allows to distinguish the expected electronic structures features in the lowest energy conduction and valence bands for different initial interlayer sliding stacking geometries AA-tBG/BN [τ = (0, 0)], AB-tBG/BN [τ = a G (0, 1/ √ 3)], and BA-tBG/BN [τ = a G (0, 2/ √ 3)] prior to rotation. While the relative sliding of the layers are shown to considerably modify the low energy electronic structure, as we illustrate in Fig. A4 through normalized local density of states (LDOS) we also find that the electron localization largely concentrates at the AA stacking sites of tBG regardless of the arrangement of the BN substrate. Quantitative differences in the normalized LDOS are still expected to be observable near AB or BA local stacking regions when relative sliding between layers are introduced. Likewise the Berry curvature maps undergo changes in their hotspot distributions depending on which stacking configuration is chosen. These differences could impact the Hall transport and optical measurements especially when valley contrasting circularly polarized light is used.