Landau Level Degeneracy in Twisted Bilayer Graphene: Role of Symmetry Breaking

TL;DR

This paper investigates how symmetry breaking, particularly $C_3$ and $C_2$ symmetries, influences Landau level degeneracy in twisted bilayer graphene near the magic angle, supported by numerical analysis.

Contribution

It proposes $C_3$ symmetry breaking as the key mechanism for degeneracy lifting in Landau levels, supported by numerical results, and explores the impact of $C_2$ symmetry breaking.

Findings

Degeneracy reduces from eight-fold to four-fold near the magic angle.

$C_3$ symmetry breaking explains robust degeneracy lifting.

$C_2$ symmetry breaking also affects Landau level structure.

Abstract

The degeneracy of Landau levels flanking charge neutrality in twisted bilayer graphene is known to change from eight-fold to four-fold when the twist angle is reduced to values near the magic angle of $\approx 1.05^\circ$. This degeneracy lifting has been reproduced in experiments by multiple groups, and is known to occur even in devices which do not harbor the correlated insulators and superconductors. We propose $C_3$ symmetry breaking as an explanation of such robust degeneracy lifting, and support our proposal by numerical results on the Landau level spectrum in near-magic-angle twisted bilayer graphene. Motivated by recent experiments, we further consider the effect of $C_2$ symmetry breaking on the Landau levels.

Figures

Equations

$$\begin{split}H=\left(\begin{array}[]{cc}\hbar v(-i\bm{\nabla}+\mathbf{q}_{3})\cdot\bm{\sigma}_{-\theta/2}&T(\mathbf{r})\\ T^{\dagger}(\mathbf{r})&\hbar v(-i\bm{\nabla}-\mathbf{q}_{2})\cdot\bm{\sigma}_{\theta/2}\end{array}\right)\end{split}$$

$$\begin{split}T^{j}=t_{M}\left(\begin{array}[]{cc}\alpha&e^{-i\frac{2\pi}{3}(j-1)}\\ e^{i\frac{2\pi}{3}(j-1)}&\alpha\end{array}\right),\end{split}$$

$$H_M=M_t{\psi }_t^{†}{\sigma }_z{\psi }_t+M_b{\psi }_b^{†}{\sigma }_z{\psi }_b$$

$$H=H_0+H_M$$

$$H_0=\sum _{\mathbf{k}}{h}_{\alpha \beta }(\mathbf{k})f_{\alpha }^{†}(\mathbf{k})f_{\beta }(\mathbf{k})$$

$$H_M=\sum _{\mathbf{k}}\sum _{{\mathbf{Q}}_{\mathbf{j}}}{T}_{\alpha \beta }({\mathbf{Q}}_{\mathbf{j}})f_{\alpha }^{†}(\mathbf{k})f_{\beta }(\mathbf{k}+{\mathbf{Q}}_{\mathbf{j}})$$

$$[{\pi }_x,{\pi }_y]=i\frac{1}{l_B^2}$$

$$[R_x,R_y]=-i{l}_B^2$$

$$a$$

$$a^{†}$$

$$a|n⟩=\sqrt{n}|n⟩$$

$$f_{\alpha ;k_x,k_y}=\sum _n{\tilde{\phi }}_{n,k_x}(k_y){\psi }_{\alpha ;n,k_x}$$

$${\phi }_{n,k_x}(y)$$

$${\tilde{\phi }}_{n,k_x}(k_y)$$

$$H_0=\sum _{k_x}\sum _{nm;\alpha \beta }{\psi }_{\alpha ;n,k_x}^{†}{h}_{\alpha \beta ;nm}{\psi }_{\beta ;m,k_x}$$

$$H_M=\sum _{k_x}\sum _{\mathbf{Q}}\sum _{nm}{T}_{\alpha \beta }(\mathbf{Q})F_{nm;k_x}(\mathbf{Q}){\psi }_{\alpha ;n;k_x}^{†}{\psi }_{\beta ;m;k_x+Q_x}$$

$$F_{nm;k_x}(\mathbf{Q})$$

$$=F_{nm}^0(\mathbf{Q})e^{i(k_x{Q}_y+\frac{1}{2}{Q}_x{Q}_y)l_B^2}$$

$$F_{nm}^0(\mathbf{Q})$$

$$=\sqrt{\frac{m!}{n!}}(z_x+i{z}_y{)}^{n-m}{e}^{-\frac{|z{|}^2}{2}}{L}_m^{n-m}(|z{|}^2)$$

$$H_0$$

$$\displaystyle+e^{-i\frac{\theta}{2}}\sqrt{n+1}f^{\dagger}_{A_{2};n+1;j}(\mathbf{k})f_{B_{2};n;j}(\mathbf{k})+h.c.\big{)}$$

$$+\sum _{\mathbf{k}}\sum _n{M}_t(f_{A_1;n,j}^{†}(\mathbf{k})f_{A_1;n,j}(\mathbf{k})-f_{B_1;n,j}^{†}(\mathbf{k})f_{B_1;n,j}(\mathbf{k}))$$

$$+\sum _{\mathbf{k}}\sum _n{M}_b(f_{A_2;n,j}^{†}(\mathbf{k})f_{A_2;n,j}(\mathbf{k})-f_{B_2;n,j}^{†}(\mathbf{k})f_{B_2;n,j}(\mathbf{k}))$$

$$H_M$$

$$T^{j}=t_{M}\left(\begin{array}[]{cc}\alpha&e^{-i\varphi(j-1)}\\ e^{i\varphi(j-1)}&\alpha\end{array}\right)$$

$$H_{C_2}=M_t{\psi }_t^{†}{\sigma }_z{\psi }_t+M_b{\psi }_b^{†}{\sigma }_z{\psi }_b$$

$$H_{eff}=\sum_{\mathbf{k}}\psi^{\dagger}(\mathbf{k})\left(\begin{array}[]{cc}M_{K,K^{\prime}}&\upsilon_{eff}(k_{x}\mp ik_{y})\\ \upsilon_{eff}(k_{x}\pm iky)&-M_{K,K^{\prime}}\end{array}\right)\psi(\mathbf{k})$$

Full text

I Introduction

The discoveryCao et al. (2018a, b) of correlated insulators and superconductivity in magic-angle twisted bilayer graphene (TBG) has sparked tremendous experimental and theoretical activity. The original results of Refs. Cao et al., 2018a, b have been confirmed and extended in Refs. Yankowitz et al., 2019; Lu et al., 2019. Ferromagnetism, accompanied by an anomalous (possibly quantized) Hall effect, has also been observed at 3/4 filling of the conduction band for some devicesSharpe et al. (2019); Lu et al. (2019). In addition, gate tunable correlated insulatorsChen et al. (2019) as well as signs of superconductivityChen et al. have been demonstrated in the moire bands of ABC trilayer graphene aligned with a hexagonal boron nitride (h-BN) substrate. Very recently, twisted double bilayers of graphene have also been studied and are shown to host spin-polarized correlated insulatorsShen et al. (2019); Liu et al. ; Cao et al. and superconductivityShen et al. (2019); Liu et al. .

In spite of the experimental progress, theoretically there is very little understanding of the many-body physics of these systems. In this paper, we address one aspect of the phenomenology of the normal metallic state of near-magic-angle twisted bilayer graphene that has resisted explanation so far. Specifically, we will focus on understanding the Landau level degeneracy near charge neutrality. As we discuss below, the same pattern of Landau fan is consistently observed across samples with varying twist angles and in a reasonably wide range of densities and temperatures. This robust experimental observation, however, has resisted a clear theoretical explanation so far.

More specifically, it is well known that the band structure of twisted bilayer graphene features two Dirac points at charge neutrality within each microscopic valley. This band touching is protected by a $C_{2}\mathcal{T}$ symmetry111$C_{2}$ refers to 180 degree spatial rotation and $\mathcal{T}$ to time reversal. We will refer to these as mini-Dirac points. With two mini-Dirac points per valley and per spin, we have a total of eight Dirac points, which is double of that in monolayer graphene. Upon doping away from neutrality, Landau levels will form out of the mini-Dirac cones. The corresponding Landau fan sequence will then be expected to be double of that of monolayer graphene, namely, $\pm 4,\pm 12,\pm 20,\dots$. Indeed, precisely this degeneracy pattern is seen in experiments on devices with relatively large twist angle of $1.8^{\circ}$ Cao et al. (2016), which is far away from the magic angle $\approx 1.1^{\circ}$.

In the vicinity of the magic angle, however, there is a surprise. The Landau fan emerging from neutrality has the Hall conductance sequence $\pm 4,\pm 8,\pm 12,\dots$. In other words, the eight-fold Landau-level degeneracy is reduced to four-fold. Experimentally, this degeneracy lifting is always seen in samples showcasing the correlated insulators and superconductivity. However, the converse is not true, i.e., the same sequence is observed even in some devices which do not show the other correlated phenomena 222 For instance, see supplementary data in Ref.  Yankowitz et al., 2019 on the $1.27^{\circ}$ device (D2) at the ambient pressure of [math] GPa. . Furthermore, the Landau fans are found to terminate once the half-filling correlated insulators set in, and the degeneracy is further reduced to two-fold on the other side of the insulator. As superconductivity is seen on both sides of the insulators, understanding the Landau fan, which conveys information on the nature of the charge carriers and the possible patterns of symmetry breaking, may shed light on the nature of the correlated phenomena closer to the magic angle.

Here, we show that the experimentally observed Landau fan can be reproduced within a free-fermion model of twisted bilayer graphene, for various angles close to the magic angle. Our main result is that the sequence is stabilized by a weak breaking of $C_{3}$ symmetry. As is well-known, $C_{3}$ symmetry breaking splits the Van Hove singularity (VHs) in the density of states Choi et al. (2019); Bi et al. . For concreteness, let us focus on the conduction band. In the energy range between the two split Van Hove singularities, the equal-energy contours in the momentum space (within each valley) take a qualitatively different shape compared to those allowed in the symmetric case: there is a single electron-like orbit which encloses both of the Dirac points. In contrast, without $C_{3}$ symmetry breaking, the equal energy contours either consist of two disjoint pockets each enclosing a Dirac point, or a single hole-like pocket enclosing the band maximum at $\Gamma$. In the presence of a magnetic field, these new equal-energy contours give rise to four-fold degenerate Landau levels stemming out from charge neutrality. We find that even for a relatively small degree of the $C_{3}$ breaking, the degeneracy in the vicinity of the magic angle is four-fold at the magnetic field strengths at which the experiments are done ($\gtrapprox 1$ T) .

We trace the origin of this effect to the large susceptibility of near-magic-angle TBG to such a $C_{3}$ breaking perturbation, which leads to a significant modification to the electronic spectrum even when the “bare” $C_{3}$ breaking is weak. Such symmetry breaking may arise due to strain effects (known to be generally present in TBG), or due to spontaneous symmetry breaking driven by interactions. While we leave open the physical origin of the $C_{3}$ breaking, we remark that the smallness of symmetry breaking required suggests that the presence of strain in the sample is a plausible explanation.

In some devices of twisted bilayer graphene, for example in Ref. Sharpe et al., 2019, it is known that there is very likely also a breaking of $C_{2}$ symmetry. For the device in Ref. Sharpe et al., 2019, this is due to the near alignment with a h-BN substrate. Such $C_{2}$ breaking gaps out the Dirac points and results in an insulator at charge neutrality. In the device studied in Ref. Lu et al., 2019, it is observed that the system is insulating at neutrality. A possible explanation is that $C_{2}$ is broken in this system as well, although it is unclear if the symmetry breaking is again due to alignment with h-BN, or is interaction driven In light of these experimental results, we also study the expected Landau levle degeneracy in devices with $C_{2}$ symmetry breaking, both with and without an additional $C_{3}$ breaking.

We remark that multiple theoretical attempts have already been made to address the neutrality Landau fan in small-angle TBG Bistritzer and MacDonald (2011a); de Gail et al. (2011); Moon and Koshino (2012); Lian et al. (2018); Hejazi et al. . While all of the emergent symmetries of TBG are preserved in these earlier works, some of them did identify specific choice of parameters for which the experimental sequence of $\pm 4,\pm 8,\pm 12,\dots$ can be observed. For instance, Ref. Hejazi et al., found numerically that the experimental sequence can be reproduced in a narrow range of parameters near the magic angle. However, the sequence is also found to depend very sensitively on the choice of model parameters Hejazi et al. , and might not explain the experimental robustness (across groups and samples) of the Landau level sequence.

Alternatively, Ref. de Gail et al., 2011 proposed a VHs-induced Landau level splitting mechanism that is essentially the same as the one we discuss below. As emphasized above, we find that this picture is valid only if $C_{3}$ symmetry breaking is invoked. This is consistent with the more microscopic calculation presented in Ref. Moon and Koshino, 2012, which did not observe the experimental sequence for the parameter range suggested in Ref. de Gail et al., 2011 when all the emergent symmetries are kept. Finally, we also note that Ref. Bi et al., suggested that layer-dependent strain could lead to the observed Landau level sequence near neutrality, although they did not perform an explicit calculation on this point. The mechanism they considered, however, is different from the one discussed in the present work. In particular, in our model the $C_{3}$ symmetry does not differentiate the two layers. Very recently Ref. Wu and Das Sarma, studied the effects of such strain on the pairing structure of the superconductor in TBG but that work did not address the Landau fans which are our interest here.

II Hofstadter Butterfly Calculations

Our starting point is the single-valley continuum model for twisted bilayer graphene Lopes dos Santos et al. (2007); Bistritzer and MacDonald (2011b). Let $q\equiv\frac{4\pi}{3a}\theta$ with $a=2.46$ Å being the lattice constant of monolayer graphene, and define the wave vectors $\mathbf{q}_{1}=-q\hat{\mathbf{y}}$, $\mathbf{q}_{2}=q(\sin\frac{2\pi}{3}\hat{\mathbf{x}}-\cos\frac{2\pi}{3}\hat{\mathbf{y}})$, and $\mathbf{q}_{3}=q(-\sin\frac{2\pi}{3}\hat{\mathbf{x}}-\cos\frac{2\pi}{3}\hat{\mathbf{y}})$. The Hamiltonian is given by

[TABLE]

where $\bm{\sigma}_{\phi}\equiv e^{i\phi\sigma_{3}/2}\bm{\sigma}e^{-i\phi\sigma_{3}/2}$, $T(\mathbf{r})=\sum_{j=1}^{3}T^{j}e^{-i\mathbf{b}_{j}\cdot\mathbf{r}}$ with

[TABLE]

and $\mathbf{b}_{j}=\mathbf{q}_{j}-\mathbf{q}_{1}$ are reciprocal lattice vectors for the moiré potential. We use $v=10^{6}$ m/s, $t_{M}=110$ meV Bistritzer and MacDonald (2011b), and $\alpha=0.8$ to incorporate the effect of lattice relaxation Koshino et al. (2018).

A perpendicular magnetic field can be incorporated by substituting $-i\bm{\nabla}\mapsto-i\bm{\nabla}-e\mathbf{A}/\hbar$, with the vector potential $\mathbf{A}=-By\,\hat{\mathbf{x}}$ in the Landau gauge. As in Ref. Bistritzer and MacDonald, 2011a, the spectrum for non-zero $B$ can be solved by first going into a Landau level basis for two decoupled layers, and then re-expressing the inter-layer potential $T(\mathbf{r})$ in this basis (more details can be found in Appendix A). In practice, in the numerical calculations we keep only finite Landau levels (typically several hundreds) close to the charge neutrality. Note that in the limit $B\rightarrow 0$ an infinite number of Landau levels should to be kept, and so we only consider fields $B>1.2$ T.

Eq. (1) has all the effective symmetries of twisted bilayer graphene. In the following, we will also consider the effect of symmetry breaking (at the single-particle level). In particular, we will set $(1-\beta)T^{1}=T^{2}=T^{3}$, such that $\beta\neq 0$ controls the degree of $C_{3}$ symmetry breaking. We will also introduce $C_{2}$ symmetry breaking through $H_{M}=M_{t}\sigma_{3}\oplus M_{b}\sigma_{3}$, which corresponds to a staggered chemical potential that would be present if the coupling to the hBN substrate is significant.

We focus on small twist angles $\theta>1.1^{\circ}$, slightly away from the true magic angle of $\theta\approx 1.05^{\circ}$ 333 Defined here as the angle at which the active band has a minimum band width, for the parameters we used. , for the following reasons: (1) In the experiment, the four-fold degeneracy is observed in a large range of angles, including in devices where the correlated insulators and superconductors are not observed. For instance, the sequence $\pm 4,\pm 8,\pm 12,\dots$ has been seen in the device at $1.27^{\circ}$ in Ref. Yankowitz et al., 2019 under ambient pressure. This suggests that the key physics leading to the four-fold degeneracy should be in effect over a relative broad range of angles. (2) Due to the very small bandwidth exactly at the magic angle, the Landau level is not well-formed within our free-electron treatment, which complicates the analysis.

III Landau Fan for the Symmetric Case

First, we show the results for the fully symmetric model, which retains the $C_{3}$ and $C_{2}$ symmetries. We will show that the Landau fan close to charge neutrality should have eight-fold degeneracy, in agreement with the experiments at twist angle $\theta=1.8^{\circ}$ but in disagreement with the results near the magic angle Cao et al. (2018a); Yankowitz et al. (2019). The spectrum is shown in Fig. 1. From neutrality, only Landau levels $n=0,\pm 1,\pm 2$ are well-developed, and the higher levels are superseded by the Landau levels coming from the top and bottom of the active bands. In a semi-classical picture, such blending of Landau levels occurs near the zero-field Van Hove singularity separating the equal-energy contours enclosing the two mini Dirac points from those enclosing the band extrema located at $\Gamma$.

Importantly, each Landau level emanating from charge neutrality is eight-fold degenerate. This degeneracy can be explained by noting the doubling due to spin, valley and layer degrees of freedom. We remark that the layer index is not generally a good quantum number; however, as these low-lying Landau levels can be understood as coming from the mini Dirac points, there is still a two-fold degeneracy arising from the mini-valley degeneracy. In contrast, in the experiments the Landau fan is only four-fold degenerate at such small twist angles.

IV Landau Fan with $C_{3}$ Breaking

In this section, we consider the effect of $C_{3}$ symmetry breaking by introducing an anisotropy $\beta$ such that $(1-\beta)T^{1}=T^{2}=T^{3}$. We will not specify the physical origin of this anisotropy in this work. Our main purpose is to demonstrate that $C_{3}$ breaking can qualitatively change the Landau fan sequence, and so long as a mean-field description remains valid, we expect our conclusions to hold even if the $C_{3}$ breaking is interaction-driven.

The Landau fan spectrum at $\beta=0.07$ is shown in Fig. 2. Compared to the $C_{3}$ invariant case, there are several clear differences. First, more Landau levels from neutrality can be seen, suggesting that the Van-Hove singularity is pushed further away from charge neutrality. Such shift of Van-Hove singularity has also been observed in STM experiment at neutralityChoi et al. (2019). Second, the previously eight-fold degenerate Landau level is split into two groups when magnetic field is larger than a critical field $B_{c}(n)$ for each Landau level $n$. $B_{c}(n=0)=6$ T and $B_{c}(n=\pm 1)=1.2$ T. For $|n|\geq 2$, $B_{c}$ is smaller than the lowest field we can reach in the calculation. Later we will argue that $B_{c}(n)\sim\frac{1}{|n|}$, which is confirmed for a smaller $C_{3}$ breaking parameter $\beta=0.01$ in the Appendix. From the splitting of the Landau levels, one expects the Landau fan sequence $\pm 4,\pm 8,\pm 12,...$ as observed in the experiments. This is confirmed in the Wannier plot in Fig. 3.

Here we try to give a simple, qualitative explanation of the splitting of the eight-fold degenerate Landau level based on a semi-classical picture derived from the band structure at zero magnetic field. With $C_{2}$ symmetry preserved, the four-fold degeneracy from spin and valley is robust. In our numerics, we have explicitly checked that the splitting happens within each fixed spin-valley sector, and so for the discussion here we can focus on a fixed sector. In the fully symmetric case (i.e., for $\beta=0$), there are two Dirac cones at the $K$ and $K^{\prime}$ points of the mini Brillouin zone (MBZ). At doping slightly above charge neutrality, the Fermi surface consists of independent electron pockets, each enclosing one of the Dirac points (Fig. 4a). As discussed in the previous section, the Landau levels arising from such orbits exhibit a two-fold mini-valley degeneracy, which, when combined with valley and spin degeneracy, leads to eight-fold degenerate Landau levels. When the doping is increased beyond the Van Hove singularity, the Fermi surface becomes a single hole pocket enclosing the $\Gamma$ point, which leads to four-fold degenerate Landau levels coming from the superlattice gap at the top of the active bands.

Importantly, the breaking of the $C_{3}$ symmetry unpins the Dirac cones from the $K$ and $K^{\prime}$ points and distorts the dispersion. Similar to before, at doping slightly above charge neutrality, the Fermi surface consists of two independent electron pockets, and one expects mini-valley degeneracy in the corresponding Landau levels. However, at higher doping the mentioned distortion opens up the possibility of a new form of Fermi surface, which is electron-like and encloses both of the Dirac points (Fig. 4b). The corresponding Landau level would appear to emerge from charge neutrality, but would no longer showcase mini-valley degeneracy. In other words, the eight-fold degenerate Landau levels split. Equivalently, the emergence of such new form of Fermi surfaces can be seen from the splitting of the Van Hove singularity as $\beta$ increases (Fig. 4c). Such splitting from $C_{3}$ symmetry breaking is also discussed in Ref. Choi et al., 2019. Note that the preceding discussion assumes $\beta\geq 0$; while the Van Hove singularity still splits for negative values of $\beta$, the new equal energy contours are open orbits and do not lead to four-fold degenerate Landau levels in the semi-classical picture. Consistently, we find that the Landau levels do not split for $\beta<0$ (Appendix.B).

As a complementary point of view, one can consider the Landau levels associated with the mini Dirac cones, which have energies $E_{n}=v_{\rm eff}\sqrt{2e\hbar nB}$ for $n\geq 0$. In this picture, the lifting of the mini-valley degeneracy originates from the hybridization between the two set of mini-Dirac Landau levels. Qualitatively, such hybridization is expected to become significant when $E_{n}\gtrapprox E_{c}$ de Gail et al. (2011), where $E_{c}$ denotes the energy of the first split Van Hove singularity arising from $C_{3}$ breaking. As such, we expect the Landau level splitting to be observable for $B\geq B_{c}(n)$ where $B_{c}(n)\propto 1/n$. This trend is in qualitative agreement with the numerics shown in Fig. 2 and Appendix C. However, we caution that such a simple picture does not explain, for instance, the splitting of the zeroth Landau level shown in Fig. 2.

V Landau Fan with $C_{2}$ Breaking

Recent experiments have seen signatures of $C_{2}$ breakingSharpe et al. (2019); Lu et al. (2019). A typical $C_{2}$ breaking term is the staggered potential:

[TABLE]

where $M_{t},M_{b}$ is the staggered potential strength for top and bottom layer. $\sigma$ denotes Pauli matrices in the $A,B$ sublattice subspace.

If the sample is well-aligned with only the top h-BN substrate, one expects $M_{t}\neq 0,M_{b}\approx 0$Zhang et al. (2019); Bultinck et al. (2019). In this case, the zeroth Landau level is always split into four groups, each of which has only two-fold spin degeneracy (see Append. D). Therefore the first two Landau fan sequences are always $\pm 2,\pm 4$, which has indeed been observed in both Refs. Sharpe et al., 2019 and Lu et al., 2019. In this case, the activation gap for the $\nu=2$ quantum Hall sequence is proportional to the magnitude of $|M_{t}|-|M_{b}|$ and therefore is almost a constant when changing magnetic field. This feature may be a useful test for this scenario in the experiment. Note that the valley splitting for zeroth Landau level is from the valley-sublattice locking and should not be thought as a simple valley Zeeman coupling in the semi-classical picture. To get a strong $\nu=2$ sequence in the single particle level, we find it necessary to invoke a strong mirror reflection breaking through a finite $|M_{t}|-|M_{b}|$. Therefore the sample in Ref. Lu et al., 2019 should also have a strong mirror reflection symmetry breaking to explain the observed $\nu=\pm 2$ sequence. The origin of the mirror reflection symmetry breaking is not clear because no obvious signature of hBN alignment is reported in Ref. Lu et al., 2019.

Next, we turn to the higher Landau levels. For $|n|\geq 1$ the valley splitting is quite small, as shown in Fig. 5. Without $C_{3}$ breaking, the $|n|\geq 1$ Landau levels from mini $K$ and $K^{\prime}$ overlap with each other and therefore we still expect eight-fold degeneracy if $C_{3}$ is a good symmetry. $C_{3}$ breaking is necessary to get the observed four fold degeneracy in higher Landau levels. One example is illustrated in Fig. 5.

In the above, we just assume the order parameters without asking about their origins. The simplest source to get $M_{t}=M$ and $M_{b}=0$ is just alignment with the top h-BN substrate. For spontaneous symmetry breaking, the most natural ansatz of spontaneous symmetry breaking in Hartree Fock calculation is actually $M_{t}=M_{b}$Xie and MacDonald (2018) and cannot reproduce the $\nu=\pm 2$ sequence. It remains a question whether the $C_{2}$, mirror reflection and $C_{3}$ breakings needed to match the Landau fan in Ref. Lu et al., 2019 are from external sources or from correlation effects.

VI Discussion

Since the first experiment on near magic angle twisted bilayer graphene, the Landau fan near neutrality have been a mystery. The observed 4-fold Landau fan is in striking contrast with the naive theoretical expectation of an 8-fold degeneracy which is also what is seen at larger twist angles. The 4-fold degeneracy is seen in devices that approach but need not be too close to magic angle. In this paper we have shown how to understand these observations through explicit calculations of the Landau fan within the continuum model of the band structure. While there exist previous such calculations with mixed success in explaining the experiments, we find that a key ingredient is to include $C_{3}$ symmetry breaking into the continuum Hamiltonian. Upon approaching the magic angle the decreasing bandwidth (and the concomitant increase in density of states) leads to an amplified response of the electronic structure to a weak $C_{3}$ breaking anisotropy. This makes it easier for the equal energy contours near the mini-Dirac points to touch and reconstruct. Beyond a low energy scale the equal energy contours then enclose both mini-Dirac points thereby reducing the Landau fan degeneracy from 8 to 4. We also consider the effect of $C_{2}$ breaking.

With this understanding let us now revisit the experimental results on various devices. The samples in Ref. Cao et al., 2018a, b; Yankowitz et al., 2019 behave as a semimetal at neutrality, suggesting a good $C_{2}\mathcal{T}$ symmetry. To explain the four fold degeneracy of Landau fan, we require the existence of $C_{3}$ breaking, either from strain or from interaction driven symmetry breaking. Such $C_{3}$ breaking has been seen in STM experimentsChoi et al. (2019); Kerelsky et al. (2018) on other near magic angle twisted bilayer graphene samples.

The sample in Ref. Sharpe et al., 2019 is aligned with the h-BN substrate on one side. This leads to an explicit breaking of both $C_{2}$ and the effective mirror (a $180^{0}$ rotation in 3d about an axis that is parallel to and bisects the 2 layers) symmetries. The $C_{2}$ breaking leads to to a gap at neutrality, consistent with the high resistivity observed at neutrality compared to the devices in Ref. Cao et al., 2018a, b; Yankowitz et al., 2019. Extensive Landau fan data is not available in Ref. Sharpe et al., 2019 though the reported Landau levels $\pm 2,\pm 4$ are reproduced by our analysis that takes into account both $C_{2}$ and mirror breaking.

In the recent study of Ref. Lu et al., 2019, a number of additional features are reported. In contrast to Ref. Cao et al., 2018a, b; Yankowitz et al., 2019 but somewhat similar to Ref. Sharpe et al., 2019 charge neutrality is insulating with a sizeable gap $\approx 0.8$ meV. Further, a Landau fan sequence $\pm 2,\pm 4,\pm 8,......$ is reported near neutrality. The insulating behavior at neutrality could be due to $C_{2}\mathcal{T}$ breaking either from explicit alignment with at least one of the h-BN substrates, or spontaneously due to interactions 444Alternately, it could be due to an interaction-driven inter-valley coherence orderPo et al. (2018).. In Ref. Lu et al., 2019 the observed Landau fan sequence near neutrality is $\pm 2,\pm 4,\pm 8,\dots$. In our calculations if we further retain the effective mirror symmetry, then we are unable to reproduce the observation of the $\pm 2$ level. Thus we conclude that the mirror symmetry must be strongly broken. An additional $C_{3}$ breaking is needed to further reduce the degeneracy of higher Landau levels to be four fold. Thus we think it likely that $C_{2},C_{3}$ and effective mirror are broken in the device in Ref. Lu et al., 2019.

The $C_{3}$ and $C_{2}$ breaking that we have crucially invoked in understanding the Landau fan clearly sets the stage on which the other correlated phenomena happen. An immediate open question is whether it is an explicit or spontaneous symmetry breaking driven by interactions. We hope to study this in the future.

Acknowledgements.

We thank Z. Bi, Y. Cao, P. Jarillo-Herrero, D. Mao, M. Yankowitz, and M. P. Zaletel for helpful discussions. TS is supported by a US Department of Energy grant DE-SC0008739, and in part by a Simons Investigator award from the Simons Foundation. HCP is supported by a Pappalardo Fellowship at MIT and a Croucher Foundation Fellowship.

Appendix A Calculation of Hofstadter butterfly

In this section we describe the procedure of our calculation of Hofstadter butterfly spectrum. The method we used is essentially the same as in Ref. Bistritzer and MacDonald, 2011a.

We first describe a general scheme for any graphene moiré superlattice systems. Then calculation for the twisted bilayer graphene is straightforward.

For moiré system with large lattice constant, the two valleys always form separated spectrum and here we only treat valley $+$ as an example. At zero magnetic field, a moiré system is described by a continuum model:

[TABLE]

$H_{0}$ is just the effective model for the valley:

[TABLE]

where $\alpha,\beta$ is the combination of layer and sublattice index.

The $H_{M}$ is the moiré hopping or superlattice potential term. It involves the scattering of momentum from $\mathbf{k}$ to $\mathbf{k}+\mathbf{G_{M}}$, where $\mathbf{G_{M}}$ is the reciprocal vector of the small mini Brillouin zon (MBZ) folded by $H_{M}$. A general $H_{M}$ term looks like

[TABLE]

When the moireé superlattice is formed by twisting two layers, the momenta of the Dirac cones (or more general band crossing points) from the two layers are different. Here we follow the notation of Ref. Bistritzer and MacDonald, 2011b and define the momentum $\mathbf{k}$ in $\psi_{\alpha}(\mathbf{k})$ as the momentum relative to the Dirac point in the corresponding layer fixed by $\alpha$. In this convention $\mathbf{Q_{j}}$ does not belong to the reciprocal vector of the MBZ.

The key idea to calculate the spectrum under magnetic field is the following. We first ignore $H_{M}$ and add magnetic field to $H_{0}$. We can easily get the new eigenstate in the Landau level basis. Then we express $H_{M}$ in the Landau level basis, which generically mix different Landau levels. Numerically we can still solve the resulting Hamiltonian and get the spectrum.

First let us solve $H_{0}$ plus magnetic field. The Hamiltonian is labeled by dynamical momenta $\pi_{x}=p_{x}-eA_{x}$ and $\pi_{y}=p_{y}-eA_{y}$ with the commutation relation:

[TABLE]

where $l_{B}^{2}=\frac{\hbar}{eB}$.

We can define guiding center coordinates $R_{x}=x+l_{B}^{2}\pi_{y}$ and $R_{y}=y-l_{B}^{2}\pi_{x}$ with commutation relation:

[TABLE]

We also have $[\pi_{a},R_{b}]=0$. $\pi_{a}$ generates the Landau level index $n$ and $R_{a}$ generates the orbital within each Landau level.

Define the ladder operator:

[TABLE]

One can easily check that $[a,a^{\dagger}]=1$. Therefore $a$ can be viewed as an annihilation operator. Landau level $n$ eigenstate satisfy:

[TABLE]

We use the Landau gauge: $\mathbf{A}=(-By,0)$, then $R_{x}=x+p_{y}l_{B}^{2}$ and $R_{y}=-p_{x}l_{B}^{2}$. We can use $R_{y}=-p_{x}l_{B}^{2}$ to label the basis of each Landau level. Each basis is then labeled by Landau level index $n$ and $k_{x}$.

We label the funciton $\varphi_{n,k_{x}}(y)$ as the wavefunction in the real space for the basis $\ket{n,k_{x}}$. The corresponding annihilation operator is labeled as $\psi_{\alpha;n,k_{x}}$, where $\alpha=A_{1},B_{1},A_{2},B_{2}$ is the sublattice index.

The continuum model can be written in terms of the operator $f_{\alpha;k_{x},ky}$, which is related to $\psi$ through

[TABLE]

where $\tilde{\varphi}_{n,k_{x}}(k_{y})$ is the Fourier transformation of $\varphi_{n,k_{x}}(y)$ with respect to $y$.

We have the following translation property:

[TABLE]

$H_{0}$ can be expressed in the following form

[TABLE]

Here $k_{x}$ is defined in the whole range of $R=[-\infty,\infty]$. In the graphene problem, $h_{nm}$ is usually very simple. For example, in monolayer graphene, it is just an off diagonal term from $B,n$ to $A,n+1$.

$H_{M}$ generically scatters momentum. In terms of $\psi_{\alpha;n,k_{x}}$, $H_{M}$ in Eq. 6 can be rewritten as

[TABLE]

where

[TABLE]

where we have used Eq. 12. We have

[TABLE]

where $(z_{x},z_{y})=\frac{1}{\sqrt{2}}(Q_{x},Q_{y})$. $L^{n}_{m}(z)$ is the associated laguerre polynomials.

If we view $k_{x}=k^{0}_{x}+x\Delta$ with $x\in Z$ as site in a 1D chain, Eq. 14 is basically a hopping term in this 1D chain. The lattice constant is $\Delta=Q_{x}$ for the smallest non-zero $Q_{x}$. Here $\alpha,n$ is just the orbital label of each site labeled by $k_{x}$. At each site $k_{x}$, the hopping term $t(k_{x})=F_{nm;k_{x}}(Q_{x})$ is a function of the site $k_{x}$. Generically we have a quasi-periodic 1D model because $t(x)\sim e^{ik_{x}Q_{y}l_{B}^{2}}\sim e^{iQ_{x}Q_{y}l_{B}^{2}x}$ is oscillating according to Eq. 15.

With the commensurate condition $\Delta Q_{y}l_{B}^{2}=2\pi\frac{q}{p}$, we have a unit cell with $p$ sites. Then we can solve the model with good quantum number $k_{1}=k_{x}^{0}\in[0,\Delta]$ and $k_{2}\in[-\frac{\pi}{p},\frac{\pi}{p}]$ as the crystal momentum of this 1D chain. For TBG, the corresponding flux per moiré unit cell is $\Phi=\frac{2\pi}{6}\frac{p}{q}\frac{h}{e^{2}}$.

Now we apply the above general scheme to the TBG. We keep $M$ Landau levels. Typical $M$ is several hundreds. For the flux corresponding to $\Phi=\frac{2\pi}{6}\frac{p}{q}\frac{h}{e^{2}}$, we can construct a $4pM$ matrix for $H_{0}+H_{M}$ corresponding to each fixed $(k_{1},k_{2})$. We use $f_{\alpha,n,j}(k_{1},k_{2})$ to label the state with momentum $k_{2}$ in the 1D chain fixed by $k_{1}$. $j=1,2,...,p$ is the sublattice index for each unit cell.

The Hamiltonian corresponding to $H_{0}$ is:

[TABLE]

Note that here $\mathbf{k}=(k_{1},k_{2})$ is just a label of the eigenstate of the 1D chain. It is not the true 2D momentum.

For TBG, $H_{M}$ is

[TABLE]

where $T^{j}_{\alpha\beta}$ is the inter-layer moiré hopping matrix. We have $T^{j}_{tb}=T^{j\dagger}_{bt}=T^{j}$. Here $t,b$ means top and bottom layer. The momentum transfer is $\mathbf{Q_{1}}=(0,\frac{4\pi}{3a_{M}})$, $\mathbf{Q_{2}}=(-\frac{2\pi}{\sqrt{3}a_{M}},-\frac{2\pi}{3a_{M}})$ and $\mathbf{Q_{3}}=(\frac{2\pi}{\sqrt{3}a_{M}},-\frac{2\pi}{3a_{M}})$.

[TABLE]

where $\varphi=\frac{2\pi}{3}$. $\alpha\leq 1$ is a parameter to incorporate lattice relaxation. We use $\alpha=0.8$ in our calculation.

Appendix B Landau levels when $\beta<0$

In the main text we show that a $C_{3}$ breaking with $\beta>0$ can reduce the degeneracy to four fold. Here we discuss the case with $\beta<0$. As is shown in Fig. 6, each Landau level from the charge neutrality is still eight fold degenerate for $\beta<0$. This implies that very specific form of $C_{3}$ is required to explain the experimental result of four fold degenerate Landau fan.

Appendix C Landau fan at larger twist angle

For $\theta=1.24^{\circ}$, the spectrum and the Wannier plots are shown in Fig. 7 and Fig. 8. Clearly even at this larger angle $C_{3}$ breaking can lead to four fold degenerate Landau fan at finite magnetic field, even though the band width is as large as $40$ meV here.

Appendix D Landau Fan with $C_{2}$ Breaking

Recent experiments have seen signatures of $C_{2}$ breakingSharpe et al. (2019); Lu et al. (2019). A typical $C_{2}$ breaking term is just the staggered potential:

[TABLE]

where $M_{t},M_{b}$ is the staggered potential strength for top and bottom layer. $\sigma$ is a Pauli matrix in the $A,B$ sublattice subspace.

For $C_{2}$ breaking coming from alignment with the top (bottom) h-BN layer, we have $M_{t}=M,M_{b}=0$ ($M_{t}=0,M_{b}=M$).

For $C_{2}$ breaking from spontaneous symmetry breaking, a natural choice is $M_{t}=M_{b}=M$Xie and MacDonald (2018). Alignment with both top and bottom h-BN layers also gives this term.

In the following we discuss the Landau fan close to neutrality from the above two different $C_{2}$ breaking ansatz.

D.1 Zeroth Landau levels

We start from the $n=0$ Landau level with an emphasis on the observed $\nu=\pm 2$ sequenceSharpe et al. (2019); Lu et al. (2019). When $C_{2}$ is preserved, there are two Dirac cones at $K$ and $K^{\prime}$ of the MBZ for each spin-valley sector.

Let us suppress the spin index. For the Dirac cone at $K$ (or $K^{\prime}$) of valley $a=\pm$, we can write a $2\times 2$ effective model:

[TABLE]

where $\psi(\mathbf{k})$ is a two component spinor. Note that the pseudospin basis of $\psi(\mathbf{k})$ is a mixture of the sublattice and layer.

When $M_{t}=M_{b}$, mirror reflection guarantees $M_{K}=M_{K^{\prime}}$. When $|M_{t}|>>|M_{b}|$, we also expect $|M_{K}|>>|M_{K^{\prime}}|$. Let us make $M_{t}=M_{b}=0$ first. Now we add magnetic field and get four zeroth Landau levels formed by $(\pm,K/K^{\prime})$ for each spin sector. In total the zeroth Landau level has eight fold degeneracy. Next we treat $M_{t}$ and $M_{b}$ as small perturbation and add them on this eight fold degenerate $n=0$ Landau level subspace. For the zeroth Landau level, it is well known that the pseudospin is polarized and locked to the valley. The $+,K$ and $+,K^{\prime}$ eighenstate only contains $\psi_{1}$ component while the $-K$ and $-K^{\prime}$ eigenstate only contains the component $\psi_{2}$. Because of this special property, $M_{K}$ and $M_{K^{\prime}}$ (resulting from $M_{t}$ and $M_{b}$) just shift the energy of the zeroth Landau level from $(\pm,K)$ ($(\pm,K^{\prime})$) by a constant energy $\pm M_{K}$ ($\pm M_{K^{\prime}}$). The final result is that we have four groups of $n=0$ Landau levels: two groups from valley $-$ at energy $-M_{K}$ and $-M_{K^{\prime}}$, two groups from valley $+$ at energy $M_{K^{\prime}}$ and $M_{K}$. This is clearly shown in Fig. 9. For $M_{t}=5$ meV, $M_{b}=0$, each of the four groups only has two fold spin degeneracy. There should be a $\nu=\pm 2$ Quantum Hall sequence with the activation gap almost constant with the magnetic field. For $|M_{t}|\approx|M_{b}|$, we have $|M_{K}|\approx|M_{K^{\prime}}|$. Therefore there are only two groups of $n=0$ Landau levels. In this case the $\nu=\pm 2$ quantum Hall sequence should be quite weak or even absent without invoking interaction effects.

D.2 Higher Landau Levels

For higher Landau levels, Ref. Lu et al., 2019 observe four fold degeneracy. Here we try to give an explanation with $C_{2}$ breaking ansatz. With $C_{2}$ breaking, the only exact degeneracy at finite magnetic field is the two fold spin degeneracy (we ignore the small spin Zeeman coupling). However, for $n\geq 1$, the valley splitting is quite weak, as shown in Fig. 9. For $M_{t}=5$ meV and $M_{b}=0$, the mirror reflection is strongly broken. However, the $n\geq 1$ Landau levels from $K$ and $K^{\prime}$ for each valley still overlap with each other. Therefore the $n\geq 1$ Landau levels have eight fold degeneracy as long as $C_{3}$ is not broken. With $C_{3}$ breaking, the Landau levels from $K$ and $K^{\prime}$ will couple with each other and reduce the degeneracy to four fold, as discussed for $C_{2}$ symmetric case. We show the Wannier plots for several cases in Fig. 10. One can see that the ansatz with $M_{t}=M_{b}$ can not give the $\nu=\pm 2$ quantum Hall sequence regardless of whether $C_{3}$ is broken. The only ansatz which can reproduce the sequence $\pm 2,\pm 4,\pm 8,\pm 12,...$ in Ref. Lu et al., 2019 is the one shown in Fig. 5 in the main text.