Dirac cones and magic angles in the Bistritzer–MacDonald TBG Hamiltonian

Simon Becker ETH Zurich, Institute for Mathematical Research, Rämistrasse 101, 8092 Zurich, Switzerland simon.becker@math.ethz.ch , Solomon Quinn School of Mathematics, University of Minesota, Minneapolis, MN, 55455, USA srquinn@umn.edu , Zhongkai Tao Department of Mathematics, University of California, Berkeley, CA 94720, USA ztao@math.berkeley.edu , Alexander Watson School of Mathematics, University of Minesota, Minneapolis, MN, 55455, USA abwatson@umn.edu and Mengxuan Yang Department of Mathematics, University of California, Berkeley, CA 94720, USA mxyang@math.berkeley.edu

Abstract.

We demonstrate the generic existence of Dirac cones in the full Bistritzer–MacDonald Hamiltonian for twisted bilayer graphene. Its complementary set, when Dirac cones are absent, is the set of magic angles. We show the stability of magic angles obtained in the chiral limit by demonstrating that the perfectly flat bands transform into quadratic band crossings when perturbing away from the chiral limit. Moreover, using the invariance of Euler number, we show that at magic angles there are more band crossings beyond these quadratic band crossings. This is the first result showing the existence of magic angles for the full Bistritzer–MacDonald Hamiltonian and solves Open Problem No. 2 proposed in the recent survey [Zw23].

1. Introduction

The discovery of superconductivity and strong electronic interactions in twisted bilayer graphene (TBG) has established TBG as an important platform for investigating strongly correlated physics within a system characterized by a topologically non-trivial band structure.

The Bistritzer-MacDonald (BM) Hamiltonian serves as the standard model for describing the effective one-particle band structure of TBG [BiMa11]. Considerable efforts in both mathematical and physics literature have been devoted to analyzing specific properties of a notable limit of this model, known as the chiral limit [TKV19]. Within this chiral limit, a discrete set of parameters has been identified [TKV19, Be*22], at which the bands closest to zero energy become entirely flat. This phenomenon is often interpreted as a key factor contributing to the system’s strongly correlated electronic properties.

It has been shown in [TKV19, Be*22, Zw23] that in the chiral limit, the band structure exhibits either a Dirac cone (i.e., a conic singularity of the band structure) at zero energy or a flat band. In short, the model exhibits Dirac cones if and only if the angle is not magic.

Since completely flat bands at zero energy are not observed away from the chiral limit, the absence of Dirac cones is studied to identify magic angles in the band structure of the self-adjoint BM Hamiltonian [BiMa11]. The BM Hamiltonian is defined as $H(\alpha,\lambda):H^{1}(\mathbb{C};\mathbb{C}^{4})\subset L^{2}(\mathbb{C};\mathbb{C}^{4})\to L^{2}(\mathbb{C};\mathbb{C}^{4})$ given by

H(\alpha,\lambda):=\begin{pmatrix}\lambda C&D(\alpha)^{*}\\ D(\alpha)&\lambda C\end{pmatrix}

(1.1)

with

D(\alpha):=\begin{pmatrix}2D_{\bar{z}}&\alpha U(z)\\ \alpha U(-z)&2D_{\bar{z}}.\end{pmatrix}\text{ and }C:=\begin{pmatrix}0&V(z)\\ \overline{V(z)}&0\end{pmatrix}

and tunnelling potentials $U,V$ satisfying symmetries for $\gamma\in\Lambda={\mathbb{Z}}+\omega{\mathbb{Z}}$

\begin{split}U(z+\gamma)=e^{i\langle\gamma,K\rangle}U(z),\quad U(\omega z)=\omega U(z),\text{ and }\overline{U(\overline{z})}=-U(-z)\\ V(z)=V(\overline{z})=\overline{V(-z)},\quad V(\omega z)=V(z),\text{ and }V(z+\gamma)=e^{i\langle\gamma,K\rangle}V(z),\end{split}

where $\omega=e^{2\pi i/3}$ and $K=\frac{4}{3}\pi$ with inner-product $\langle a,b\rangle:=\operatorname{Re}(a\bar{b})$ for $a,b\in{\mathbb{C}}.$ Honeycomb lattices consist of two non-equivalent vertices per fundamental domain that we denote by $A$ and $B$ . Then parameters $\alpha$ and $\lambda$ tune the strengths of the tunnelling potentials of the $A/B$ and $A/A,B/B$ regions, respectively and are inversely proportional to the twisting angle $\theta$ for $\theta\ll 1.$ The chiral limit is obtained by setting $\lambda=0$ in the Hamiltonian (1.1), so $\lambda$ also represents how far the model is from this limit. We refer the reader to Section 2.1 for a detailed discussion of the special form of this Hamiltonian and additional definitions on the tunnelling potentials. The Hamiltonian (1.1) commutes with the translation operator

\mathscr{L}_{\gamma}\mathbf{w}(z):=\operatorname{diag}(\omega^{\gamma_{1}+\gamma_{2}},1,\omega^{\gamma_{1}+\gamma_{2}},1)\mathbf{w}(z+\gamma),\ \ \ \gamma\in\Lambda,\quad\gamma=\gamma_{1}+\omega\gamma_{2},\ (\gamma_{1},\gamma_{2})\in{\mathbb{Z}}^{2}.

(1.2)

Thus, by Bloch–Floquet theory, the Hamiltonian is equivalent to a family of Bloch transformed operators

	$\displaystyle H_{k}(\alpha,\lambda)=\begin{pmatrix}\lambda C&D(\alpha)^{*}+\bar{k}\\ D(\alpha)+k&\lambda C\end{pmatrix}:L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})\to L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4}),$		(1.3)
	$\displaystyle L^{2}_{k}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4}):=\{\mathbf{v}\in L^{2}_{\operatorname{loc}}(\mathbb{C},{\mathbb{C}}^{4}):\mathscr{L}_{\gamma}\mathbf{v}=e^{i\langle k,\gamma\rangle}\mathbf{v},\ \ \gamma\in\Lambda\},\quad k\in{\mathbb{C}}/\Lambda^{*}$		(1.4)

such that

\operatorname{Spec}_{L^{2}(\mathbb{C};\mathbb{C}^{4})}(H(\alpha,\lambda))=\bigcup_{k\in{\mathbb{C}}/\Lambda^{*}}\operatorname{Spec}_{L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}(H_{k}(\alpha,\lambda)),

where $\Lambda^{*}$ is the dual lattice associated with $\Lambda.$ The spectrum of $H_{k}(\alpha,\lambda)$ on $L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})$ is discrete, and we label it as follows

\begin{gathered}\{E_{\pm j}(\alpha,\lambda;k)\}_{j\in\mathbb{N}_{+}},\ \ E_{\pm 1}(\alpha,\lambda;-K)=E_{\pm 1}(\alpha,\lambda;0)=0,\\ \ \cdots\leq E_{-2}(\alpha,\lambda;k)\leq E_{-1}(\alpha,\lambda;k)\leq 0\leq E_{1}(\alpha,\lambda;k)\leq E_{2}(\alpha,\lambda;k)\leq\cdots,\end{gathered}

(1.5)

forming the Bloch bands. The points $k=-K,0$ are high-symmetry points and are typically denoted by $K$ and $K^{\prime}$ in the physics literature¹¹1We shift them by $K$ for notational convenience in the computation.. At points $k=0,-K$ , there exist two protected states (cf. [BeZw23-2, Proposition 2])

\varphi_{k}(\alpha,\lambda)=(u,v)^{T},\ \psi_{k}(\alpha,\lambda)=(\tilde{u},\tilde{v})^{T}\in\ker_{L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}(H_{k}(\alpha,\lambda)),\ k=0,-K

(1.6)

for all $\lambda,\alpha\in{\mathbb{C}}$ explaining $E_{\pm 1}(\alpha,\lambda;-K)=E_{\pm 1}(\alpha,\lambda;0)=0$ in (1.5). As we shall argue, the protected states play an essential role in analyzing the presence of Dirac cones in this model.

As mentioned above, we shall study the band structure of the BM Hamiltonian (1.1) near the points $0,-K$ in the Brillouin zone ${\mathbb{C}}/\Lambda^{*}$ . In particular, we investigate the existence of conic singularities (see Figure 1) in $E_{\pm 1}(\alpha,\lambda;k)$ for $(\alpha,\lambda)\in{\mathbb{R}}^{2}$ . First, we introduce the following

Definition 1.1 (Simple Dirac cone).

Assume $\dim\ker_{L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}(H_{0}(\alpha,\lambda))=2$ . We say that the BM Hamiltonian (1.1) exhibits a simple Dirac cone at $k=0$ at $(\alpha,\lambda)$ if and only if

E_{\pm 1}(\alpha,\lambda;k)=\pm|v_{F}(\alpha,\lambda)||k|+\mathcal{O}(|k|^{2})\text{ with }v_{F}(\alpha,\lambda)\neq 0.

We call such parameters $(\alpha,\lambda)$ non-magic. When $v_{F}(\alpha,\lambda)=0$ , i.e. the dispersion relation is not linear, we call the parameters $(\alpha,\lambda)$ magic.

Remark 1.

The Dirac cone at $k=-K$ is defined analogously. To simplify the presentation, we prove all results near $k=0$ . In fact, the two points are connected by the symmetries $\mathscr{S}$ and $\mathscr{M}$ that will be defined in (2.14).

Refer to caption — Figure 1. The two bands $E_{\pm 1}(\alpha,\lambda;k)$ closest to zero energy of BM Hamiltonian exhibiting two Dirac points for parameters $\alpha=0.7$ and $\lambda=0.3.$

We start by proving a theorem on the generic existence (cf. Figure 2) of Dirac cones for the bands $E_{\pm 1}(\alpha,\lambda,k)$ near $k=0$ .

Theorem 1 (Generic existence of Dirac cones).

There is a locally finite family of points and analytic curves $\{S_{i}\}$ such that the BM Hamiltonian (1.1) exhibits a simple Dirac cone at $k=0$ for all $(\alpha,\lambda)\in{\mathbb{R}}^{2}\setminus\bigcup_{i}S_{i}$ .

Next we show that there exist real parameters $(\alpha,\lambda)\in{\mathbb{R}}^{2}$ at which the Dirac cones disappear, as suggested by Figure 2 numerically.

Recall that in the chiral limit of the BM Hamiltonian $H(\alpha,0)$ where $\lambda=0$ in (1.1), it has been proven (cf. [TKV19, Be*22, BHZ22-2]) that there exists a discrete set $\mathcal{A}\subset{\mathbb{C}}$ such that

\alpha\in\mathcal{A}\iff E_{\pm 1}(\alpha,0;k)=0\ \ \forall k\in{\mathbb{C}},

which means that there is a flat band at zero energy for $\alpha\in\mathcal{A}$ . If $E_{2}(\alpha,0,k)$ does not vanish for all $k$ as well, then $\alpha$ is called simple. The existence of the first real magic angle $\alpha\in\mathcal{A}\cap{\mathbb{R}}$ is shown in [WaLu21, BHZ22-1] for the BM potential (2.3) used in the physics literature [BiMa11]. It is shown in [Zw23, Appendix] and [TKV19] that $(\alpha,0)$ is non-magic in the sense of Definition 1.1 if and only if $\alpha\notin\mathcal{A}$ . This means that for the chiral Hamiltonian $H(\alpha,0)$ , the Dirac cone at $k=0$ disappears precisely for magic $\alpha\in\mathcal{A}$ at which the Hamiltonian exhibits a flat band. The next theorem generalizes this result to the BM Hamiltonian (1.1). The theorem shows that magic angles of the chiral model, persist in the full BM Hamiltonian. This is illustrated in Figure 3.

Theorem 2 (Persistence of chiral magic angles).

Assume $\partial_{\alpha}v_{F}(\alpha_{0},0)\neq 0$ (cf. Proposition 3.2) for some $\alpha_{0}\in\mathcal{A}\cap\mathbb{R}$ simple. There is $\varepsilon>0$ and a real valued real analytic function $\lambda\mapsto\alpha(\lambda)$ for $|\lambda|<\varepsilon$ with $\alpha(0)=\alpha_{0}$ and $\alpha^{\prime}(0)=0$ such that $(\alpha,\lambda)$ is magic along the curve $(\alpha(\lambda),\lambda)$ for $\lambda\in(-\varepsilon,\varepsilon)$ .

On the other hand, if $\alpha_{0}\notin\mathcal{A}$ , then there exists a neighborhood $U_{\alpha_{0}}\subset{\mathbb{R}}^{2}$ of $(\alpha_{0},0)$ such that $(\alpha,\lambda)$ is non-magic on $U_{\alpha_{0}}$ .

Remark 2.

For the standard BM potential (2.3), the condition $\partial_{\alpha}v_{F}(\alpha_{0},0)\neq 0$ in Theorem 2 is numerically checked for the first six magic angles (see Table 1).

In Section 5, we study band touching near magic parameters. We show that when transitioning from conic singularities to quadratic band touching at Dirac points $-K,0$ , the change of winding number implies the existence of other band crossings at points other than $-K,0$ points as shown in Figure 7. In fact, we prove the following

Theorem 3.

For magic parameters $(\alpha,\lambda)$ near simple chiral magic angles $\alpha_{0}\in\mathcal{A}$ and $\lambda=0$ , suppose the first two bands $E_{\pm 1}(\alpha,\lambda,k)$ are isolated from other bands and exhibit quadratic band crossing in the sense that near $k=0$

E_{\pm 1}(\alpha,\lambda,k)=c_{\pm}(\alpha,\lambda)|k|^{2}+\mathcal{O}(|k|^{3}),\quad c_{\pm}(\alpha,\lambda)\neq 0.

Then the two bands touch at some additional points other than $k=0,-K$ .

Remark 3.

The assumption of the theorem can be weaken: it holds as long as there exists a path $\gamma$ in ${\mathbb{R}}^{2}$ connecting $(\alpha_{0},0)$ and $(\alpha,\lambda)$ such that the first two bands $E_{\pm 1}(\alpha,\lambda,k)$ are gapped from the rest of bands along $\gamma$ . See Section 5 for more details.

To understand the transition from conical intersection in the band structure to quadratic ones, we need to develop the study of topological phases for our model. The mechanism behind this involves a winding number known as the Euler number which classifies the triviality of real vector bundles of rank two. The reality condition appears in this model due to the $C_{2z}T$ symmetry.

Here we only discuss the BM model (1.1) but our argument can be generalized to studying band touchings of two-band systems gapped from the rest of the bands.

See also [ChWe24] for discussions of the splitting of quadratic band touchings into Dirac cones for Schrödinger operators and [BeZw23-1] for the formation and splitting of quadratic band touchings by in-plane fields in twisted bilayer graphene.

The analysis of conical intersections in band structures has received a considerable attention and has been studied for honeycomb lattice structures [FW12], describing materials such as graphene, as well as optical and acoustic analogues [Am*20]. The first results on the dispersion relation of tight-binding models for honeycomb structures have been obtained more than half a century ago (see [Wa47, SlWe58]). A mathematical model of honeycomb quantum graphs with even potential on edges for graphene is studied by Kuchment–Post [KuPo07]. The Schrödinger operator $-\Delta+V$ with smooth real valued potentials $V$ exhibiting honeycomb symmetry has been studied by Grushin [Gr09] for small potentials using perturbation theory method. For a generic set of smooth potentials, the existence of Dirac cones has been established by Fefferman–Weinstein [FW12, FLW16] and later by Berkolaiko–Comech [BeCo18]. This is generalized to potentials with singularities at honeycomb lattice points in the work of Lee [Le16]. For generalizations to different class of elliptic operators defined on the honeycomb lattice, see [CaWe21, LWZ19, LLZ23, Ca24]. See also the survey by Kuchment [Ku23]. The omnipresence of conical singularities in the context of topological insulators has been analyzed by Drouot [Dr21]. For twisted bilayer graphene, in the chiral limit of the BM Hamiltonian, the existence of Dirac cone has been proved by some of the authors of this paper in [Zw23, Appendix] and in [TKV19].

Structure of the paper

•

In Section 2, we review properties of the BM Hamiltonian which includes a basic derivation and symmetries of the BM Hamiltonian and the existence of protected states.
•

In Section 3, we prove the generic existence of Dirac cones.
•

In Section 4, we extend the notion of magic angles, show the stability of magic angles obtained from the chiral limit to the general BM Hamiltonian, and prove the existence of magic angles.
•

In Section 5, we discuss Euler numbers and winding numbers of bands and use them to study the behavior and crossing of bands near magic parameters.

Acknowledgements

The authors are very grateful to Maciej Zworski for many helpful discussions and proposing this project. The authors would also like to thank Gregory Berkolaiko, Patrick Ledwith, Qiuyu Ren and Oskar Vafek for helpful discussions. Their insights were key to the development of the paper. We would like to thank Jens Wittsten for allowing us to use the left figure in Figure 4 that he created. ZT and MY were partially supported by the NSF grant DMS-1952939 and by the Simons Targeted Grant Award No. 896630. AW’s research was supported by the NSF grant DMS-2406981. SB acknowledges support by the SNF Grant PZ00P2 216019.

2. Symmetries of the Hamiltonian and protected states

In this section we briefly review the physical origin, symmetries, and Bloch-Floquet theory of the Hamiltonian (1.1).

2.1. BM Hamiltonian

Twisted bilayer graphene consists of two stacked graphene layers with a relative interlayer twist. The atomic structure of each graphene layer is the union of two offset Bravais lattices of carbon atoms, referred to as the $A$ and $B$ sublattices. The single-particle electronic properties of twisted bilayer graphene can be modeled by the following Bistritzer-MacDonald Hamiltonian $\mathcal{H}:L^{2}(\mathbb{R}^{2};\mathbb{C}^{4})\rightarrow L^{2}(\mathbb{R}^{2};\mathbb{C}^{4})$ ,

\mathcal{H}(\beta,\zeta)=\begin{pmatrix}-i\boldsymbol{\sigma}\cdot\nabla&T(\beta,\zeta\mathbf{r})\\ T(\beta,\zeta\mathbf{r})^{\dagger}&-i\boldsymbol{\sigma}\cdot\nabla\end{pmatrix},

(2.1)

acting on functions $\psi(\boldsymbol{r})=(\psi_{1}^{A},\psi_{1}^{B},\psi_{2}^{A},\psi_{2}^{B})^{\top}$ , $\boldsymbol{r}=(x,y)$ , where $|\psi_{i}^{\sigma}|^{2}$ represents electronic probability densities on layer $i$ and sublattice $\sigma$ .

In the diagonal (intralayer) terms, we have $\boldsymbol{\sigma}=(\sigma_{x},\sigma_{y})$ , where $\sigma_{\bullet}$ are the Pauli matrices, and $\boldsymbol{\sigma}\cdot\nabla=\sigma_{x}\partial_{x}+\sigma_{y}\partial_{y}$ . These terms capture the effective Dirac dispersion for electrons in monolayer graphene near to the Fermi level. The off-diagonal terms describe interlayer tunnelling, which is approximated by the interlayer tunnelling matrix potential

T(\beta,\zeta\mathbf{r})=\begin{pmatrix}\beta_{AA}V(\zeta\mathbf{r})&\beta_{AB}\overline{U(-\zeta\mathbf{r})}\\ \beta_{BA}U(\zeta\mathbf{r})&\beta_{AA}V(\zeta\mathbf{r})\end{pmatrix}.

Here, $\beta=(\beta_{AA},\beta_{AB})$ are the tunnelling potential amplitudes between like and unlike sublattices, $\zeta$ is the relative twist angle (here we assume $\zeta$ is small so that we can use the approximation $\sin(\zeta)\approx\zeta$ ), and $U$ and $V$ are given by

\begin{gathered}U(\mathbf{r})=\sum_{i=0}^{2}\omega^{i}e^{-iq_{i}\cdot\mathbf{r}},\ \ V(\mathbf{r})=\sum_{i=0}^{2}e^{-iq_{i}\cdot\mathbf{r}},\ \ q_{i}:=R^{i}\begin{pmatrix}\ \ 0\\ -1\end{pmatrix},\ \ R:=\tfrac{1}{2}\begin{pmatrix}-1&-\sqrt{3}\\ \sqrt{3}&-1\end{pmatrix},\end{gathered}

where $\omega=e^{2\pi i/3}$ . The model (2.1) was introduced in [BiMa11]; for its mathematical justification and detailed discussion of the various approximations involved in deriving (2.1), see [Wa*22, Ko*24, Q*24, Ca23]. If we denote the lattice of periodicity of the Hamiltonian (2.1) by $\Gamma$ , the Hamiltonian (2.1) commutes with translations in the moiré lattice $\frac{1}{3}\Gamma$ up to cubic root of unity phases; see e.g. [Wa*22] for details.

The Hamiltonian (1.1) is obtained from (2.1) as follows. Let $K=\operatorname{diag}(1,\sigma_{x},1)$ . Then, upon making the change of variables $\zeta\mathbf{r}\mapsto\mathbf{r}$ , $\lambda=\beta_{AA}/\zeta$ , and $\alpha=\beta_{AB}/\zeta$ , we find

\zeta^{-1}K\mathcal{H}(\beta,\zeta)K=\begin{pmatrix}\lambda C&D(\alpha)^{*}\\ D(\alpha)&\lambda C\end{pmatrix},

where

\begin{split}D(\alpha)&=\begin{pmatrix}D_{x}+iD_{y}&\alpha U(\mathbf{r})\\ \alpha U(-\mathbf{r})&D_{x}+iD_{y}\end{pmatrix}\text{ and }C=\begin{pmatrix}0&V(\mathbf{r})\\ V(-\mathbf{r})&0\end{pmatrix},\end{split}

which now acts on functions $\psi(\boldsymbol{r})=(\psi_{1}^{A},\psi_{2}^{A},\psi_{1}^{B},\psi_{2}^{B})^{\top}$ . Changing coordinates from $\boldsymbol{r}=(x,y)$ to $z$ so that $x+iy=\frac{4}{3}\pi iz$ we obtain (1.1), where $D(\alpha)$ and $C$ are given by (we abuse notation to write functions of $\boldsymbol{r}$ and $z$ by the same letters)

D(\alpha)=\begin{pmatrix}2D_{\bar{z}}&\alpha U(z)\\ \alpha U(-z)&2D_{\bar{z}}\end{pmatrix}\text{ and }C=\begin{pmatrix}0&V(z)\\ V(-z)&0\end{pmatrix},

(2.2)

where

U(z)=\sum_{k=0}^{2}\omega^{k}e^{\frac{1}{2}(z\bar{\omega}^{k}-\bar{z}\omega^{k})}\text{ and }V(z)=\sum_{k=0}^{2}e^{\frac{1}{2}(z\bar{\omega}^{k}-\bar{z}\omega^{k})}.

(2.3)

The moiré lattice in these coordinates is $\Lambda=a_{1}{\mathbb{Z}}+a_{2}{\mathbb{Z}},$ where $a_{\ell}=\omega^{\ell-1}$ . This lattice and its associated reciprocal lattice $\Lambda^{*}$ are related to those of the physical space moiré by

\tfrac{1}{3}\Gamma=\tfrac{4}{3}\pi i\Lambda,\quad 3\Gamma^{*}=\frac{3}{4\pi i}\Lambda^{*}.

As discussed above, throughout this work we make the more general assumption that $U$ and $V$ are arbitrary smooth functions satisfying the symmetries

\begin{gathered}\overline{U(z)}=U(\bar{z}),\ \ U(\omega z)=\omega U(z),\ \ U(z+\gamma)=\bar{\omega}^{\gamma_{1}+\gamma_{2}}U(z)\\ V(z)=V(\bar{z})=\overline{V(-z)},\ \ V(\omega z)=V(z),\ \ V(z+\gamma)=\bar{\omega}^{\gamma_{1}+\gamma_{2}}V(z),\end{gathered}

(2.4)

with $\gamma=\gamma_{1}a_{1}+\gamma_{2}a_{2}\in\Lambda.$ Since our proofs rely only on symmetries of the model (1.1), this generalization does not complicate our proofs, but we use (2.3) for our numerical computations for simplicity.

Remark 4.

The assumption (2.3) on the potentials $U,V$ amounts to an approximation where momentum-space hopping is truncated to nearest-neighbors in the moiré reciprocal lattice. This is generally an excellent approximation since the non-zero spacing between the graphene sheets causes the magnitude of momentum-space hops to decay exponentially with distance [Wa*22, BiMa11]. Note that mechanical relaxation [Ca18, Ca20] complicates this picture by enhancing the strength of longer-range momentum hops; see e.g. [Ma23].

Remark 5.

It is natural to ask whether our results could be generalized beyond the other simplifying approximations implicit in the model (2.1). One of these approximations is neglecting the rotation of the monolayer Dirac cones, see, e.g. [BiMa11, equation (8)]. Another is the neglect of derivative terms in the interlayer tunneling, see e.g. [Ca19]. Since these terms preserve the translation and rotation symmetries of the model, and the $\mathcal{P}\mathcal{T}$ symmetry defined below, all of our results should apply to models including these terms. However, we do not consider these generalizations in the present work because the statement of our results would become much more complicated. This is because these terms break the particle-hole symmetry $\mathscr{S}$ defined below, so that there is no need for Dirac cones to occur at energy $0$ . Our results would then have to allow for Dirac cones to occur in a suitable range of energies.

2.2. Symmetries revisited

We start by recalling the basic translational and rotational symmetries of (1.1). The modified translation operator $\mathscr{L}_{\gamma},$ defined in (1.2), commutes with the Hamiltonian $H(\alpha,\lambda)$ . We also define the rotation operator

\Omega u(z):=u(\omega z),\ \ u\in\mathscr{S}^{\prime}(\mathbb{C};\mathbb{C}^{2}),\ \ \Omega D(\alpha)=\omega D(\alpha)\Omega.

We then extend it to a commuting action with $H(\alpha,\lambda)$ by introducing

\mathscr{C}:=\begin{pmatrix}\Omega&\ 0\\ 0&\bar{\omega}\Omega\end{pmatrix}:L_{\rm{loc}}^{2}(\mathbb{C};\mathbb{C}^{4})\to L^{2}_{\rm{loc}}(\mathbb{C};\mathbb{C}^{4})

(2.5)

such that

\mathscr{C}H(\alpha,\lambda)=H(\alpha,\lambda)\mathscr{C},\ \ \mathscr{L}_{\gamma}\mathscr{C}=\mathscr{C}\mathscr{L}_{\omega\gamma},\ \ \mathscr{C}\mathscr{L}_{\gamma}=\mathscr{L}_{\bar{\omega}\gamma}\mathscr{C}.

From these two symmetries, it is possible [BeZw23-1, (3.10),(3.15)] to obtain the following orthogonal decomposition

\begin{gathered}L^{2}_{k}(\mathbb{C}/\Lambda)=\bigoplus_{p\in\mathbb{Z}_{3}}L^{2}_{k,p}(\mathbb{C}/\Lambda),\ \ \ k\in\mathcal{K}/\Lambda^{*},\end{gathered}

(2.6)

where

\ \mathcal{K}:=\{k\in\mathbb{C}:\omega k\equiv k\!\!\mod\Lambda^{*}\}=\{K,-K,0\}+\Lambda^{*},\quad K=\frac{4}{3}\pi

(2.7)

and

\begin{gathered}L^{2}_{k,p}(\mathbb{C}/\Lambda;\mathbb{C}^{4}):=\{u\in L^{2}_{\rm{loc}}(\mathbb{C};\mathbb{C}^{4}):\mathscr{L}_{\gamma}\mathscr{C}^{\ell}u=e^{i\langle k,\gamma\rangle}\bar{\omega}^{\ell p}u\},\\ k\in(\tfrac{1}{3}\Lambda^{*})/\Lambda^{*}\simeq\mathbb{Z}^{2}_{3},\ \ p\in\mathbb{Z}_{3}.\end{gathered}

(2.8)

2.3. Additional symmetries

We recall additional symmetries of the BM Hamiltonian using the notation of [BeZw23-2] that are central in proving the existence of Dirac cones. We start with the symmetries of $D(\alpha)$ , the non-self-adjoint building block of the BM Hamiltonian (1.1) for $U$ satisfying the conditions in (2.4).

We recall an anti-linear symmetry,

{\mathscr{Q}}v(z):=\overline{v(-z)},\ \ \ \ {\mathscr{Q}}(D(\alpha)+k){\mathscr{Q}}=(D(\alpha)+k)^{*},

(2.9)

and two linear symmetries,

\mathscr{H}\begin{pmatrix}u_{1}(z)\\ u_{2}(z)\end{pmatrix}:=\begin{pmatrix}-iu_{2}(-z)\\ iu_{1}(-z)\end{pmatrix},\ \ \ \ \mathscr{H}(D(\alpha)+k)\mathscr{H}=-(D(\alpha)-k),

(2.10)

\mathscr{N}\begin{pmatrix}u_{1}(z)\\ u_{2}(z)\end{pmatrix}:=\begin{pmatrix}u_{2}(-\bar{z})\\ u_{1}(-\bar{z})\end{pmatrix},\ \ \ \ \mathscr{N}(D(\alpha)+k)\mathscr{N}=-(D(\bar{\alpha})-\bar{k})^{*}.

(2.11)

All these symmetries are involutions, i.e. satisfy $A^{2}=I$ with $A:L^{2}_{0}\to L^{2}_{0}$ . For potentials satisfying (2.4) these symmetries extend to symmetries of the BM Hamiltonian (1.1)

$\mathcal{P}\mathcal{T}$ symmetry:	$\displaystyle\mathcal{P}\mathcal{T}=\begin{pmatrix}0&\mathscr{Q}\\ \mathscr{Q}&0\end{pmatrix},$	(2.12)
particle-hole symmetry:	$\displaystyle\mathscr{S}=\begin{pmatrix}\mathscr{H}&0\\ 0&\mathscr{H}\end{pmatrix},$
mirror symmetry:	$\displaystyle\mathscr{M}=\begin{pmatrix}\ 0&i\mathscr{N}\\ -i\mathscr{N}&0\end{pmatrix}.$

We recall the following proposition from [BeZw23-2] summarizing basic properties of the aforementioned symmetries

Proposition 2.1.

For the BM Hamiltonian given in (1.1) with potentials satisfying (2.4),and $\alpha\in\mathbb{C}$ , $\lambda\in\mathbb{R}$ , we have

\begin{gathered}\mathcal{P}\mathcal{T}H(\alpha,\lambda)=H(\alpha,\lambda)\mathcal{P}\mathcal{T},\ \ \ \ \mathscr{M}H(\alpha,\lambda)=H(\bar{\alpha},\lambda)\mathscr{M},\\ \mathscr{S}H(\alpha,\lambda)=-H(\alpha,\lambda)\mathscr{S}.\end{gathered}

(2.13)

Additionally, for $k=\{0,\pm K\}+\Lambda^{*}$ , $p\in\mathbb{Z}_{3}$ , and $K=\frac{4}{3}\pi$ ,

\mathcal{P}\mathcal{T}:L^{2}_{k,p}\to L^{2}_{k,-p+1},\ \ \mathscr{S}:L^{2}_{k,p}\to L^{2}_{-k-K,p},\ \ \mathscr{M}:L^{2}_{k,p}\to L^{2}_{-k-K,1-p},

(2.14)

where $L^{2}_{\bullet,\bullet}:=L^{2}_{\bullet,\bullet}(\mathbb{C}/\Lambda;\mathbb{C}^{4})$ .

As a final ingredient from [BeZw23-2], we recall the existence of protected states of the BM Hamiltonian at $k=0$ . (See [GZ23] for a more abstract formulation.)

Proposition 2.2.

In the notation of (1.3) and (2.8) and for all $\alpha,\lambda\in\mathbb{R}$ ,

\dim\ker_{L^{2}_{0,0}}H_{0}(\alpha,\lambda)\geq 1,\ \ \dim\ker_{L^{2}_{0,1}}H_{0}(\alpha,\lambda)\geq 1.

(2.15)

3. Existence of Dirac cones

We prove Theorem 1 in this section. We start with a general argument on perturbations of self-adjoint operators with ${\mathbb{Z}}/3{\mathbb{Z}}$ symmetry. The proof of this theorem relies on the Schur complement formula, which in the version that we use, is often referred to as a Grushin problem. For a general discussion of Grushin problems, we refer to [SjZw07] and [TaZw23, Section 2.6].

Proposition 3.1.

Suppose $H_{k}=H_{0}+kA+\bar{k}A^{*}:\mathcal{D}\to\mathcal{H}$ is a family of (unbounded) self-adjoint operators on a Hilbert space $\mathcal{H}$ with domain $\mathcal{D}$ . In addition, we assume that

•

$\mathcal{H}$ has an orthogonal decomposition $\mathcal{H}=\mathcal{H}_{0}\oplus\mathcal{H}_{1}\oplus\mathcal{H}_{2}$ such that $H_{0}:\mathcal{H}_{j}\cap\mathcal{D}\to\mathcal{H}_{j}$ ;
•

$H_{0}$ has discrete spectrum at $0$ and there exist $\varphi\in\mathcal{H}_{0}$ , $\psi\in\mathcal{H}_{1}$ (normalized) such that $\ker H_{0}={\mathbb{C}}\varphi+{\mathbb{C}}\psi$ ;
•

$A$ is bounded and $A:\mathcal{H}_{j}\to\mathcal{H}_{j+1}$ (with the convention $\mathcal{H}_{3}=\mathcal{H}_{0}$ ).

Then $H_{k}$ has two eigenvalues $E_{\pm 1}(k)$ near $0$ satisfying

E_{\pm 1}(k)=\pm|\langle A\varphi,\psi\rangle||k|+\mathcal{O}(|k|^{2}).

(3.1)

Proof.

Let $R_{+}:\mathcal{D}\to{\mathbb{C}}^{2}$ and $R_{-}:{\mathbb{C}}^{2}\to\mathcal{H}$ be defined by

R_{+}u=(\langle u,\varphi\rangle,\langle u,\psi\rangle),\quad R_{-}(a,b)=a\varphi+b\psi.

(3.2)

Consider the following Grushin problem

\mathcal{P}_{k}:=\begin{pmatrix}H_{k}-z&R_{-}\\ R_{+}&0\end{pmatrix}:\mathcal{D}\times{\mathbb{C}}^{2}\to\mathcal{H}\times{\mathbb{C}}^{2}

(3.3)

as a perturbation of the Grushin problem

\mathcal{P}_{0}=\begin{pmatrix}H_{0}-z&R_{-}\\ R_{+}&0\end{pmatrix}:\mathcal{D}\times{\mathbb{C}}^{2}\to\mathcal{H}\times{\mathbb{C}}^{2}.

(3.4)

For sufficiently small $\varepsilon>0$ , the operator $\mathcal{P}_{0}$ in (3.4) is invertible for $z\notin(\operatorname{Spec}(H_{0})\setminus\{0\})+(-\varepsilon,\varepsilon)$ with inverse given by

\mathcal{P}_{0}^{-1}=:\begin{pmatrix}E&E_{+}\\ E_{-}&E_{-+}\end{pmatrix}

(3.5)

and operators

E=(\Pi^{\perp}(H_{0}-z)\Pi^{\perp})^{-1}\Pi^{\perp},\quad E_{+}=R_{-},\quad E_{-}=R_{+},\quad E_{-+}=z.

(3.6)

Here, $\Pi=R_{-}R_{+}$ is the projection to $\ker H_{0}$ and $\Pi^{\perp}=I-\Pi$ . The operator $\mathcal{P}_{k}$ in (3.3) is invertible for $|k|$ sufficiently close to zero with inverse

\mathcal{P}_{k}^{-1}=:\begin{pmatrix}F&F_{+}\\ F_{-}&F_{-+}\end{pmatrix}.

(3.7)

Here, for $|k|$ small enough

	$\displaystyle F_{-+}$	$\displaystyle=E_{-+}+\sum\limits_{j=0}^{\infty}(-1)^{j}E_{-}(kA+\bar{k}A^{})(E(kA+\bar{k}A^{}))^{j-1}E_{+}$		(3.8)
		$\displaystyle=z-R_{+}(kA+\bar{k}A^{})R_{-}+R_{+}(kA+\bar{k}A^{})E(kA+\bar{k}A^{*})R_{-}+\mathcal{O}(\|k\|^{3}),$		(3.8)

is invertible if and only if $H_{k}-z$ is invertible; see [TaZw23, Lemma 2.10]. In other words, for all $|k|$ sufficiently small, $\det F_{-+}=0$ if and only if $z$ is an eigenvalue of $H_{k}$ . We compute

R_{+}(kA+\bar{k}A^{*})R_{-}=\begin{pmatrix}0&\bar{k}\langle A^{*}\psi,\varphi\rangle\\ k\langle A\varphi,\psi\rangle&0\end{pmatrix}

(3.9)

and

R_{+}(kA+\bar{k}A^{*})E(kA+\bar{k}A^{*})R_{-}=\begin{pmatrix}|k|^{2}\langle(AEA^{*}+A^{*}EA)\varphi,\varphi\rangle&k^{2}\langle AEA\psi,\varphi\rangle\\ \bar{k}^{2}\langle A^{*}EA^{*}\varphi,\psi\rangle&|k|^{2}\langle(A^{*}EA+AEA^{*})\psi,\psi\rangle\end{pmatrix}.

(3.10)

We have

F_{-+}-R_{+}(kA+\bar{k}A^{*})R_{-}=\mathcal{O}(|k|^{2}).

Since $F_{-+}$ and $R_{+}(kA+\bar{k}A^{*})R_{-}$ are both self-adjoint, their eigenvalues near zero differ by $\mathcal{O}(|k|^{2})$ . Since the eigenvalues of (3.9) are $\pm|\langle A\varphi,\psi\rangle||k|$ , the eigenvalues of $E_{\pm 1}(k)$ are by [Ka13, Th. 4.10 p. 291] given by

E_{\pm 1}(k)=\pm|\langle A\varphi,\psi\rangle||k|+\mathcal{O}(|k|^{2}).\qed

Remark 6.

We can compute more terms in the asymptotic expansion of the inverse of the Grushin problem (3.3) as illustrated in (3.8). In our application to $H_{k}=H_{0}+k\begin{pmatrix}0&0\\ I_{2\times 2}&0\end{pmatrix}+\bar{k}\begin{pmatrix}0&I_{2\times 2}\\ 0&0\end{pmatrix}$ , we have an antilinear $\mathcal{PT}$ -symmetry:

\mathcal{PT}:\mathcal{H}\to\mathcal{H},\quad\mathcal{PT}^{2}=I,\quad\mathcal{PT}(\mathcal{H}_{0})=\mathcal{H}_{1},\,\,\mathcal{PT}(\mathcal{H}_{1})=\mathcal{H}_{0},\,\,\mathcal{PT}(\mathcal{H}_{2})=\mathcal{H}_{2}

(3.11)

where the Hilbert spaces are the respective invariant subspaces of the Hamiltonian, so that

H_{0}\mathcal{PT}=\mathcal{PT}H_{0},\quad\mathcal{PT}A=A^{*}\mathcal{PT}.

Consequently, we can choose $\mathcal{PT}\varphi=\psi$ . This immediately implies that

\langle EA^{*}\varphi,A^{*}\varphi\rangle+\langle EA\varphi,A\varphi\rangle=\langle EA\psi,A\psi\rangle+\langle EA^{*}\psi,A^{*}\psi\rangle

and the energies are given by

\begin{split}E_{\pm 1}(k)&=\pm|k\langle A\varphi,\psi\rangle-\bar{k}^{2}\langle(\Pi^{\perp}H_{0}\Pi^{\perp})^{-1}A^{*}\varphi,A\psi\rangle|\\ &\quad-|k|^{2}\langle(\Pi^{\perp}H_{0}\Pi^{\perp})^{-1}(A+A^{*})\varphi,(A+A^{*})\varphi\rangle+\mathcal{O}(|k|^{3}).\end{split}

(3.12)

Remark 7.

Under the condition $\langle A\varphi,\psi\rangle\neq 0$ , the bands $E_{\pm 1}(k)$ have a conic singularity near $k=0$ . We also notice when the first order term $\langle A\varphi,\psi\rangle$ vanishes, we get

E_{\pm 1}(k)=c_{2,\pm}(\varphi,\psi)|k|^{2}+\mathcal{O}(|k|^{3})

where, with the notation introduced in the proof of Proposition 3.1,

c_{2,\pm}(\varphi,\psi)=\pm|\langle(\Pi^{\perp}H_{0}\Pi^{\perp})^{-1}A^{*}\varphi,A\psi\rangle|-\langle(\Pi^{\perp}H_{0}\Pi^{\perp})^{-1}(A+A^{*})\varphi,(A+A^{*})\varphi\rangle.

We next apply Proposition 3.1 to the BM Hamiltonian (1.1). We also prove that, because of symmetries, we can assume that $\langle A\varphi,\psi\rangle$ in equation (3.1) is real-valued. This fact will be important in the proof of Theorem 2.

Proposition 3.2.

Suppose there is an open set $\Omega\subset{\mathbb{R}}^{2}$ such that

\dim\ker_{L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}\left(H_{0}(\alpha,\lambda)\right)=2,\quad(\alpha,\lambda)\in\Omega.

(3.13)

Then there exists a real analytic function $v_{F}(\alpha,\lambda):\Omega\to{\mathbb{R}}$ such that

E_{\pm 1}(\alpha,\lambda;k)=\pm|v_{F}(\alpha,\lambda)||k|+\mathcal{O}(|k|^{2}).

(3.14)

In particular, when $v_{F}(\alpha,\lambda)\neq 0$ , the bands exhibit a conic singularity at energy zero with slope $|v_{F}|$ .

Proof.

We consider the spectrum of $H_{k}(\alpha,\lambda)$ on $\mathcal{H}=L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})$ . It satisfies the assumptions of Proposition 3.1 with

\mathcal{H}_{j}=L^{2}_{0,j}\text{ and }A=\begin{pmatrix}0&0\\ I_{2\times 2}&0\end{pmatrix},\text{ and }\ker_{L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}\left(H_{0}(\alpha,\lambda)\right)={\mathbb{C}}\varphi(\alpha,\lambda)+{\mathbb{C}}\psi(\alpha,\lambda).

Here $\varphi(\alpha,\lambda)\in L^{2}_{0,0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})$ and $\psi(\alpha,\lambda)\in L^{2}_{0,1}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})$ are normalized protected states.

We first claim that the spectral projection is analytic for $(\alpha,\lambda)\in\Omega$ . For any $(\alpha_{0},\lambda_{0})\in\Omega$ , by assumption (3.13), there is a neighbourhood of $(\alpha_{0},\lambda_{0})$ in which we can choose $\#\operatorname{Spec}_{L^{2}_{0,j}}(H_{0}(\alpha,\lambda))\cap[-\delta,\delta]=1$ , $j=0,1$ . The spectral projection can then be written as a contour integral

\mathbbm{1}_{0}(H_{0}(\alpha,\lambda)|_{L^{2}_{0,j}})=\frac{1}{2\pi i}\int_{\partial B_{0}(\delta)}(z-H_{0}(\alpha,\lambda))^{-1}dz|_{L^{2}_{0,j}},\quad j=0,1.

(3.15)

The analyticity of $H_{0}(\alpha,\lambda)$ implies that the spectral projection $\Pi_{j}=\mathbbm{1}_{0}(H_{0}(\alpha,\lambda)|_{L^{2}_{0,j}})$ is real analytic.

We now show that $\langle A\varphi,\psi\rangle$ can be chosen to be real valued. Suppose

\varphi(\alpha,\lambda)=\begin{pmatrix}u(\alpha,\lambda)\\ v(\alpha,\lambda)\end{pmatrix},\quad\psi(\alpha,\lambda)=\begin{pmatrix}\tilde{u}(\alpha,\lambda)\\ \tilde{v}(\alpha,\lambda)\end{pmatrix},

then $\langle A\varphi,\psi\rangle=\langle u,\tilde{v}\rangle$ . Recall the particle-hole and mirror symmetries (2.12) from Section 2.2. Note that by Proposition 2.1,

	$\displaystyle\mathscr{S}\mathscr{M}:L^{2}_{0,p}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})\to L^{2}_{0,1-p}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4}),\ \ p=0,\pm 1,$
	$\displaystyle\mathscr{S}\mathscr{M}H_{k}(\alpha,\lambda)=-H_{\bar{k}}(\overline{\alpha},\lambda)\mathscr{S}\mathscr{M}.$

Since we are only considering $(\alpha,\lambda)\in{\mathbb{R}}^{2}$ , $\mathscr{SM}$ maps $\ker_{L^{2}_{0,0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}\left(H_{k}(\alpha,\lambda)\right)$ to the other $\ker_{L^{2}_{0,1}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}\left(H_{\bar{k}(\alpha,\lambda)}\right)$ for $k=\{0,\pm K\}+\Lambda^{*}$ . In particular, for $k=0$ , we can choose $\psi(\alpha,\lambda)=\mathscr{SM}\varphi(\alpha,\lambda)$ . In other words, using the notation in (2.12) and the fact that $(\mathscr{HN})^{2}=-I$ ,

\tilde{v}=-i\mathscr{HN}u,\quad\mathscr{HN}\tilde{v}=iu.

Therefore, since $\mathscr{H}$ and $\mathscr{N}$ are unitary, we obtain

\langle u,\tilde{v}\rangle=\langle\mathscr{HN}u,\mathscr{HN}\tilde{v}\rangle=\langle i\tilde{v},iu\rangle=\langle\tilde{v},u\rangle

and $v_{F}(\alpha,\lambda)=\langle u,\tilde{v}\rangle$ is real valued.

Finally we claim that we can choose $v_{F}(\alpha,\lambda)$ to be real analytic in $\Omega$ . For this we write $v_{F}=\langle A\varphi,\psi\rangle$ using the spectral projection $\Pi_{j}=\mathbbm{1}_{0}(H_{0}(\alpha,\lambda)|_{L^{2}_{0,j}})$ . Note $\Pi_{0}(\alpha,\lambda)\phi=\langle\phi,\varphi(\alpha,\lambda)\rangle\varphi(\alpha,\lambda)$ . Using the $\mathscr{SM}$ symmetry $\psi=\mathscr{SM}\varphi$ , we see that

v_{F}(\alpha,\lambda)=\langle A\varphi(\alpha,\lambda),\mathscr{SM}\varphi(\alpha,\lambda)\rangle=-\operatorname{Tr}_{L^{2}_{0}}\mathscr{SM}A\Pi_{0}(\alpha,\lambda),

as $\mathscr{SM}$ is unitary with $(\mathscr{SM})^{2}=-\mathrm{Id}$ . Since $\Pi_{0}$ is analytic in $(\alpha,\lambda)\in\Omega$ with $\mathscr{SM}$ and $A$ independent of $(\alpha,\lambda)$ , we conclude that $v_{F}(\alpha,\lambda)$ is analytic in $\Omega$ . ∎

Remark 8.

Note that only the choice of an analytic real valued Fermi velocity $v_{F}(\alpha,\lambda)\in{\mathbb{R}}$ relies on the symmetry $\mathscr{SM}$ , whereas the form of equation (3.14) with only analytic complex-valued Fermi velocity $v_{F}(\alpha,\lambda)\in{\mathbb{C}}$ does not rely on the $\mathscr{SM}$ -symmetry. In general, $v_{F}(\alpha,\lambda)$ is uniquely defined up to a phase.

Remark 9.

The condition (3.13) is satisfied when $\lambda=0$ and $\alpha\notin\mathcal{A}$ , or when $\lambda=0$ and $\alpha\in\mathcal{A}$ is a simple magic angle. The condition (3.13) is also stable under perturbation, i.e. the set of $(\alpha,\lambda)$ satisfying (3.13) is open.

The next proposition establishes that the nullspace $\ker_{L^{2}_{0}}(H_{0}(\alpha,\lambda))$ is generically two dimensional. This is done by identifying a real analytic function whose zero set coincides with the set where this condition fails.

Proposition 3.3.

There exists a locally finite family of points and analytic curves $S_{i}\subset{\mathbb{R}}^{2}$ such that for $(\alpha,\lambda)\in\mathbb{R}^{2}\setminus\bigcup_{i}S_{i}$ , $\dim\ker_{L^{2}_{0}}H(\alpha,\lambda)=2$ .

Proof.

We define the Fredholm determinant (cf. [S77, (1.2)] or [DyZw19, Appendix B.5.2])

D(\alpha,\lambda;z)=\det(I+i(H(\alpha,\lambda)+i-z)^{-3}).

(3.16)

Note that for any $(\alpha,\lambda)\in{\mathbb{R}}^{2}$ , the function $D(\alpha,\lambda;z)$ has at least a zero of order two at $z=0$ by the existence of the protected state and is real analytic in $(\alpha,\lambda)$ .

Since $F(\alpha,\lambda):=\partial^{2}_{zz}D(\alpha,\lambda;z)\rvert_{z=0}$ is a real analytic function for $(\alpha,\lambda)\in{\mathbb{R}}^{2}$ and $F(0,0)\neq 0$ , by Łojasiewicz’s structure theorem (see [KrPa02, Theorem 6.3.3]) the zeros of $F(\alpha,\lambda)$ is the union of a locally finite family of points and analytic curves $S_{i}\subset{\mathbb{R}}^{2}$ . Our proposition follows since $F(\alpha,\lambda)\neq 0$ implies $\dim\ker_{L^{2}_{0}}H_{0}(\alpha,\lambda)=2$ . ∎

We are now in a position to almost prove Theorem 1 by the following argument. First, recall that the BM Hamiltonian exhibits a simple Dirac cone at $k=0$ at $(\alpha,\lambda)$ if and only if $\dim\ker_{L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}\left(H_{0}(\alpha,\lambda)\right)=2$ and $\langle A\varphi(\alpha,\lambda),\psi(\alpha,\lambda)\rangle\neq 0$ . Applying Proposition 3.3, we have that the first condition fails only on a locally finite family of points and analytic curves $S_{i}\subset{\mathbb{R}}^{2}$ . By Proposition 3.2, we have that the function $v_{F}(\alpha,\lambda)$ is real analytic everywhere outside of the sets ${S_{i}}$ . Applying Łojasiewicz’s structure theorem ([KrPa02, Theorem 6.3.3]) to a neighborhood of any point in $\mathbb{R}^{2}\setminus\bigcup_{i}S_{i}$ , we have that the zero set of $v_{F}(\alpha,\lambda)$ must be a finite collection of points and analytic curves.

However, this argument does not prove Theorem 1 in full, because the zero set of $v_{F}$ could still accumulate at points in $\bigcup_{i}S_{i}$ . In order to prove Theorem 1 in full, we recall the theory of Kurdyka–Paunescu [KuPa08], whose results can be summarized as follows. First, eigenfunctions of operators depending analytically on parameters $(\alpha,\lambda)\in\mathbb{R}^{2}$ can be chosen analytically everywhere except for a locally finite set of points. In particular, eigenfunctions can generically be chosen analytically in neighborhoods of codimension 1 eigenvalue crossings. Second, at points where analytic eigenfunctions do not exist, the eigenfunctions can again be chosen analytically by lifting these points to appropiate “blowup spaces”. We can then prove Theorem 1 by characterizing the zero set of $v_{F}$ in these blowup spaces and projecting this set down to the original parameter space.

For the reader’s convenience, we recall [KuPa08, Example 6.1], which demonstrates how “blowing up” the parameter space makes it possible to choose eigenfunctions analytically at points where this is impossible otherwise.

Example 1.

Let

A(x_{1},x_{2})=\begin{pmatrix}x_{1}^{2}&x_{1}x_{2}\\ x_{1}x_{2}&x_{2}^{2}\end{pmatrix},\quad(x_{1},x_{2})\in{\mathbb{R}}^{2}.

Then the eigenvalues of $A(x_{1},x_{2})$ are given by $\lambda_{1}=0$ and $\lambda_{2}=x_{1}^{2}+x_{2}^{2}$ . The corresponding normalized eigenvectors are

\displaystyle\phi_{1}=\frac{1}{\sqrt{x_{1}^{2}+x_{2}^{2}}}(x_{2},-x_{1}),\qquad\phi_{2}=\frac{1}{\sqrt{x_{1}^{2}+x_{2}^{2}}}(x_{1},x_{2}),

which are not continuous at the origin, even though the eigenvalues are real analytic. However, if we blow up the origin in ${\mathbb{R}}^{2}$ , i.e., if we take $x_{1}=w_{1},x_{2}=w_{1}w_{2}$ , then the corresponding family

A(w_{1},w_{2})=w_{1}^{2}\begin{pmatrix}1&w_{2}\\ w_{2}&w_{2}^{2}\end{pmatrix},\quad(w_{1},w_{2})\in{\mathbb{R}}^{2}.

admits a simultaneous analytic diagonalization.

This motivates us to introduce the blowup space.

Definition 3.4.

The blowup space ${\rm Bl}_{(0,0)}{\mathbb{R}}^{2}$ of ${\mathbb{R}}^{2}$ at the point $(0,0)$ is

\{((x,y),[\xi,\eta])\in{\mathbb{R}}^{2}\times\mathbb{RP}^{1}:x\eta=y\xi\},

where $\mathbb{RP}^{1}$ is the real projective line. The blow down map $\pi:{\rm Bl}_{(0,0)}{\mathbb{R}}^{2}\to{\mathbb{R}}^{2}$ is given by

((x,y),[\xi,\eta])\mapsto(x,y).

Under this definition, for a neighbourhood $\mathcal{U}=B(0,R)$ of zero, its blowup $\widetilde{\mathcal{U}}=\pi^{-1}(\mathcal{U})$ is given by

\widetilde{\mathcal{U}}=\{((x,y),[\xi,\eta])\in B(0,R)\times\mathbb{RP}^{1}:x\eta=y\xi\},

which again is a two dimensional analytic manifold. One can similarly blow up at any other point. Following Kurdyka–Paunescu [KuPa08], we have the following proposition.

Proposition 3.5.

For any point in ${\mathbb{R}}^{2}$ , there is a neighbourhood $\mathcal{U}$ and a blowup space $\widetilde{\mathcal{U}}$ obtained by a finite composition of blowups such that $\pi:\widetilde{\mathcal{U}}\to\mathcal{U}$ is an isomorphism outside finitely many points $P_{1},\cdots,P_{r}\in\mathcal{U}$ , and there exists a normalized real analytic family of functions

\varphi(u)\in L^{2}_{0,0},\quad\psi(u)\in L^{2}_{0,1},\quad u\in\widetilde{\mathcal{U}}

such that for $u\notin\pi^{-1}(\{P_{1},\cdots,P_{r}\})$ , $\varphi(u),\psi(u)\in\ker H_{0}(\pi(u)).$

Proof.

We may just restrict ourselves to $L^{2}_{0,0}$ . The analysis for $L^{2}_{0,1}$ is similar. First we reduce the problem to a finite dimensional problem. Since $H_{0}(\alpha,\lambda)$ is elliptic and self-adjoint, the eigenvalue at $0$ is isolated. Let $\Pi(\alpha,\lambda)$ be the spectral projector to a neighbourhood of $0$ that is analytic in a neighbourhood $\mathcal{U}$ and we consider $H_{0}(\alpha,\lambda)$ on the space $\Pi(\alpha,\lambda)L^{2}_{0,0}$ . Since this vector space will in general depend on $(\alpha,\lambda)$ , we may choose a basis of $\Pi(\alpha_{0},\lambda_{0})L^{2}_{0,0}$ , say $e_{1},\cdots,e_{N}$ for $(\alpha_{0},\lambda_{0})\in\mathcal{U}$ . As long as $\|\Pi(\alpha,\lambda)-\Pi(\alpha_{0},\lambda_{0})\|<1$ , which locally holds by the holomorphic functional calculus, see (3.15), we can use the Kato-Nagy formula [Ka55]

W(P,Q):=(I-(P-Q)^{2})^{-1/2}(PQ+(I-P)(I-Q)),

(3.17)

such that $W(P,Q)\ Q\ W(P,Q)^{-1}=P$ for $P=\Pi(\alpha,\lambda)$ and $Q=\Pi(\alpha_{0},\lambda_{0})$ to define $w(\alpha,\lambda):=W(\Pi(\alpha,\lambda),\Pi(\alpha_{0},\lambda_{0}))$ such that locally $(w(\alpha,\lambda)e_{i})_{i\in\{1,..,N\}}$ is an analytic orthonormal basis of $\operatorname{ran}(\Pi(\alpha,\lambda))$ in a perhaps smaller neighbourhood that we still denote by $\mathcal{U}$ . This naturally defines an analytic section, see e.g. [TaZw23, Corr. 2.17-2.18]. We may then consider the matrix $M(\alpha,\lambda):=(\langle e_{i},w(\alpha,\lambda)^{*}H_{0}(\alpha,\lambda)w(\alpha,\lambda)e_{j}\rangle)_{i,j}$ representing $H_{0}(\alpha,\lambda)$ acting on this basis. Now $M(\alpha,\lambda)$ is an analytic family of self-adjoint matrices. We can then apply [KuPa08, Theorem 6.2] to conclude the proposition. ∎

Proof of Theorem 1.

Recall that by Definition 1.1 the simple Dirac cone exists at $(\alpha,\lambda)$ if and only if the kernel $\ker_{L^{2}_{0}}H_{0}(\alpha,\lambda)$ is two dimensional and $v_{F}(\alpha,\lambda)\neq 0$ . By Proposition 3.3, the vector space $\ker_{L^{2}_{0}}H_{0}(\alpha,\lambda)$ is two dimensional for $(\alpha,\lambda)\notin\cup_{i}S_{i}$ . Using Proposition 3.5, the function $v_{F}(\alpha,\lambda)$ is locally analytic for $(\alpha,\lambda)\notin\{P_{1},\cdots,P_{r}\}$ , and can be turned into an analytic function after lifting to a blowup space $\widetilde{\mathcal{U}}$ . We can then apply Łojasiewicz’s structure theorem to conclude that the zero set of $v_{F}$ is a locally finite union of points and analytic curves. Note that when we apply Łojasiewicz’s structure theorem in the blowup space we have to check that the blowdown of the zero set remains a finite set of points and analytic curves. But this follows directly from the definition of the blowdown map. ∎

Corollary 3.6.

The set

S:=\{(\alpha,\lambda)\in\mathbb{R}^{2};\text{ the Hamiltonian $H(\alpha,\lambda)$ does not exhibit Dirac points}\}

has Hausdorff dimension

\dim_{\mathrm{Haus}}(S)=\inf\{d\geq 0:H^{d}(S)=0\}\leq 1,

where $H^{d}$ is the $d$ -dimensional Hausdorff measure.

Proof.

This follows from the previous theorem by noticing that for $d>1$

H^{d}\bigg{(}\bigcup_{i}S_{i}\bigg{)}=\sum_{i}H^{d}(S_{i})=0.\qed

4. Magic angles: from chiral limit to the BM Hamiltonian

In this section we prove Theorem 2. By Remark 9 and Propostion 3.2, the Fermi velocity $v_{F}(\alpha,\lambda)=\langle u,\tilde{v}\rangle$ is real valued and real analytic in a neighborhood of a simple magic angle $(\alpha_{0},0)$ with $\alpha_{0}\in\mathcal{A}$ . We refer to Table 1, where the assumption $\partial_{\alpha}v_{F}(\alpha_{0},0)\neq 0$ is numerically verified for the first several magic angles. We have the following expansion for the Fermi velocity.

Lemma 4.1.

Near simple magic parameters in the chiral limit $(\alpha_{0},0)$ with $\alpha_{0}\in\mathcal{A}$ , the Taylor expansion of $v_{F}(\alpha,\lambda)$ is given by

v_{F}(\alpha,\lambda)=c_{10}\alpha^{\prime}+c_{02}\lambda^{2}+c_{11}\alpha^{\prime}\lambda+\mathcal{O}(\alpha^{\prime 2})+\mathcal{O}(|(\alpha^{\prime},\lambda)|^{3}),

(4.1)

where $\alpha^{\prime}:=\alpha-\alpha_{0}$ and the coefficients $c_{10},c_{02}$ , and $c_{11}$ are given by (4.6), (4.7), and (4.8), respectively.

Proof.

We write the Taylor expansion

v_{F}(\alpha,\lambda):=c_{0}+c_{10}\alpha^{\prime}+c_{01}\lambda+c_{02}\lambda^{2}+c_{11}\alpha^{\prime}\lambda+\mathcal{O}(\alpha^{\prime 2})+\mathcal{O}(|(\alpha^{\prime},\lambda)|^{3}).

Note that $c_{0}=0$ as $(\alpha_{0},0)$ is a magic parameter and thus, $v_{F}(\alpha,0)=0$ , due to the presence of a flat band. In a neighborhood of a simple magic parameter $(\alpha_{0},0)$ , by the proof of Proposition 3.2, we can choose protected states $\varphi\in L^{2}_{0,0},\psi\in L^{2}_{0,1}$ analytic in $(\alpha,\lambda)$ such that

\displaystyle\varphi(\alpha_{0},\lambda;z)=\begin{pmatrix}u_{0}+\sum_{i=1}^{\infty}\lambda^{i}u_{0,i}\\ \sum_{i=1}^{\infty}\lambda^{i}v_{0,i}\end{pmatrix},\ \ \psi(\alpha_{0},\lambda;z)=\begin{pmatrix}\sum_{i=1}^{\infty}\lambda^{i}\tilde{u}_{0,i}\\ \tilde{v}_{0}+\sum_{i=1}^{\infty}\lambda^{i}\tilde{v}_{0,i}\end{pmatrix}.

(4.2)

Equations

H_{0}(\alpha_{0},\lambda)\varphi=H_{0}(\alpha_{0},\lambda)\psi=0

(4.3)

yield that, for $k\in{\mathbb{N}}_{+}$ , we have

	$\displaystyle D(\alpha_{0})^{*}v_{0,k}+Cu_{0,k-1}=0,\ D(\alpha_{0})u_{0}=D(\alpha_{0})u_{0,1}=0,\ D(\alpha_{0})u_{0,k}+Cv_{0,k-1}=0;$		(4.4)
	$\displaystyle D(\alpha_{0})^{}\tilde{v}_{0}=D(\alpha_{0})^{}\tilde{v}_{0,1}=0,\ D(\alpha_{0})^{*}\tilde{v}_{0,k}+C\tilde{u}_{0,k-1}=0,\ D(\alpha_{0})\tilde{u}_{0,k}+C\tilde{v}_{0,k-1}=0.$		(4.5)

For $\alpha_{0}\in\mathcal{A}$ , by [Zw23, (A.3)] and equations (4.4) and (4.5), we have

c_{01}=0,\ c_{02}=\langle\tilde{v}_{0},u_{02}\rangle+\langle\tilde{v}_{02},u_{0}\rangle

(4.6)

as $u_{01}$ (resp. $\tilde{v}_{01}$ ) is a scalar multiple of $u_{0}$ (resp. $\tilde{v}_{0}$ ).

To get the coefficient $c_{10}$ , for $\lambda=0$ , we consider $H(\alpha_{0}+\alpha^{\prime})\varphi=H(\alpha_{0}+\alpha^{\prime})\psi=0$ for $\alpha$ small with $\alpha_{0}\in\mathcal{A}$ . By analyticity, write $\psi=(0,\tilde{v})$ with $\tilde{v}=\sum_{i=0}^{\infty}\alpha^{\prime i}\tilde{v}_{i,0}$ and $\varphi=(u,0)$ with $u=\sum_{i=0}^{\infty}\alpha^{\prime i}u_{i,0}$ . Using equation $D(\alpha_{0}+\alpha^{\prime})u=D(\alpha_{0}+\alpha^{\prime})^{*}\tilde{v}=0$ and considering the coefficient of $\alpha^{\prime k}$ yield that

D(\alpha_{0})u_{k,0}+Bu_{k-1,0}=0,\quad B:=\begin{pmatrix}0&U(z)\\ U(-z)&0\end{pmatrix}.

We have

c_{10}=\langle\tilde{v}_{10},u_{0}\rangle+\langle\tilde{v}_{0},u_{10}\rangle.

(4.7)

Analogous computations as above show that

c_{11}=\langle\tilde{v}_{0},u_{11}\rangle+\langle\tilde{v}_{01},u_{10}\rangle+\langle\tilde{v}_{10},u_{01}\rangle+\langle\tilde{v}_{11},u_{0}\rangle.

(4.8)

Summarizing the previous computations, we have established the Taylor expansion

v_{F}(\alpha,\lambda)=c_{10}\alpha^{\prime}+c_{02}\lambda^{2}+c_{11}\alpha^{\prime}\lambda+\mathcal{O}(\alpha^{\prime 2})+\mathcal{O}(|(\alpha^{\prime},\lambda)|^{3}),

(4.9)

with coefficients given by equations (4.6), (4.7), and (4.8). ∎

Proof of Theorem 2.

Recall that by Proposition 3.2, $v_{F}(\alpha,\lambda)=\langle u,\tilde{v}\rangle\in{\mathbb{R}}$ is real analytic in $(\alpha,\lambda)$ with $v_{F}(\alpha_{0},0)=0$ for $\alpha_{0}\in\mathcal{A}\cap{\mathbb{R}}$ . By the implicit function theorem and our assumption that $\partial_{\alpha}v_{F}(\alpha_{0},0)\neq 0$ , there exists a real valued analytic function $(-\varepsilon,\varepsilon)\ni\lambda\mapsto\alpha(\lambda)$ with $\alpha(0)=\alpha_{0}$ and $\alpha^{\prime}(0)=0$ such that $(\alpha,\lambda)$ is magic along the curve $(\alpha(\lambda),\lambda)$ for $\lambda\in(-\varepsilon,\varepsilon)$ . The real analyticity of the curve follows from the real analyticity of $v_{F}(\alpha,\lambda)$ . The vanishing of the first derivative $\alpha^{\prime}(0)=0$ follows from the Taylor expansion (4.1) of $v_{F}(\alpha,\lambda)$ , as $\partial_{\lambda}v_{F}(\alpha_{0},0)=0$ .

Conversely, assume $\alpha_{0}\notin\mathcal{A}$ . We then know by [Zw23, Appendix] and [TKV19] that $(\alpha_{0},0)$ is non-magic in the sense of Definition 1.1. This means $\dim\ker_{L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}(H_{0}(\alpha_{0},0))=2$ and $v_{F}(\alpha_{0},0)\neq 0$ . Since $v_{F}(\alpha,\lambda)$ and the eigenvalues of $H_{0}(\alpha,\lambda)$ are continuous in $(\alpha,\lambda)$ , it follows that $v_{F}(\alpha,\lambda)\neq 0$ and $\dim\ker_{L^{2}_{0}({\mathbb{C}}/\Lambda;{\mathbb{C}}^{4})}(H_{0}(\alpha,\lambda))=2$ for all $(\alpha,\lambda)$ in a sufficiently small neighborhood of $(\alpha_{0},0)$ . ∎

$\alpha$	$\|c_{10}\|$	$\|c_{02}\|$
0.58566355838955	1.5641	0.0493
2.2211821738201	0.4130	0.0973
3.7514055099052	0.1291	1.4239
5.276497782985	0.0355	9.8783
6.79478505720	0.0091	52.5993
8.3129991933	0.0021	252.5188

Table 1. Values of

|c_{10}|

and

|c_{02}|

at first six magic angles in the Taylor expansion (4.1).

5. Magic angles in the BM Hamiltonian

In this section we discuss the topology of the bands. In particular, we show that when there is a quadratic band touching, there will always be other band touching points away from the high symmetry points $-K,0$ . This is due to the topological invariance of the Euler number, which we shall introduce next. In the physics literature, the relevance of the Euler number has also been pointed out in the influential article [APY19] in the context of fragile topology which inspired this section.

5.1. Winding number and Euler number

We recall the notion of winding number and Euler number. Let $\mathcal{E}\to{\mathbb{C}}/\Lambda^{*}$ be an oriented rank two real vector bundle over the torus ${\mathbb{C}}/\Lambda^{*}$ . Suppose there is a connection on $\mathcal{E}$ with curvature

\Omega=\begin{pmatrix}0&\Omega_{12}\\ -\Omega_{12}&0\end{pmatrix}.

The Euler number is

e_{2}(\mathcal{E})=\frac{1}{2\pi}\int_{{\mathbb{C}}/\Lambda^{*}}\operatorname{Pf}(\Omega)=\frac{1}{2\pi}\int_{{\mathbb{C}}/\Lambda^{*}}\Omega_{12},

where $\operatorname{Pf}(\cdot)$ denotes the Pfaffian of the matrix. The Euler number $e_{2}(\mathcal{E})$ gives a complete classification of real rank two oriented bundles $\mathcal{E}$ over the torus up to isomorphism.

We give another characterisation of the Euler number when there is a flat connection $\nabla$ on the rank two vector bundle $\mathcal{E}$ over the torus ${\mathbb{C}}/\Lambda^{*}$ defined outside finitely many points. First we define the winding number of the flat connection. Suppose the connection is orthogonal, i.e., it preserves a given metric $g$ and thus satisfies $\nabla g=0$ . Let $\gamma_{p}(t)$ be a loop around $p\in{\mathbb{C}}/\Lambda^{*}$ and $e_{1}(k)$ , $e_{2}(k)$ be an orthonormal local frame in a neighbourhood of $p$ . Let $\theta(t)$ be the connection $1$ -form under the basis $\{e_{1}(k),e_{2}(k)\}$ , consider

\begin{pmatrix}0&a\\ -a&0\end{pmatrix}=\frac{1}{2\pi}\oint_{\gamma_{p}}\theta.

The winding number ${\rm ind}(\nabla,p)$ is defined to be the value of $a$ . Its exponential

\exp 2\pi\begin{pmatrix}0&a\\ -a&0\end{pmatrix}=\begin{pmatrix}\cos 2\pi a&\sin 2\pi a\\ -\sin 2\pi a&\cos 2\pi a\end{pmatrix}

(5.1)

is the monodromy matrix. From this one can then define the Euler number.

Proposition 5.1.

Suppose there are finitely many points $p_{1},\cdots,p_{\ell}\in{\mathbb{C}}/\Lambda^{*}$ such that $\nabla$ is a flat connection on the vector bundle $\mathcal{E}$ over $({\mathbb{C}}/\Lambda^{*})\setminus\{p_{1},\cdots,p_{\ell}\}$ . Then the Euler number is given by the sum of the winding numbers around the points $p_{j}$ :

e_{2}(\mathcal{E})=\sum\limits_{j=1}^{\ell}{\rm ind}(\nabla,p_{j}).

Proof.

We modify the flat connection smoothly in a small neighbourhood $U_{j}$ of the points $\{p_{j}\}$ . Suppose the connection matrix is given by $\theta\in\Omega^{1}$ , then the Euler number is given by

e_{2}(\mathcal{E})=\frac{1}{2\pi}\int_{{\mathbb{C}}/\Lambda^{*}}\operatorname{Pf}(d\theta)=\frac{1}{2\pi}\sum_{j=1}^{\ell}\limits\int_{U_{j}}\operatorname{Pf}(d\theta),

as curvature vanishes outside $U_{j}$ . Let $\gamma_{j}=\partial U_{j}$ be loops around $p_{j}$ , then Stokes’ theorem gives

e_{2}(\mathcal{E})=\frac{1}{2\pi}\sum\limits_{j=1}^{\ell}\oint_{\gamma_{j}}\operatorname{Pf}(\theta),

where the right hand side is the sum of winding numbers ${\rm ind}(\nabla,p_{j})$ . ∎

5.2. Topology of band with $\mathcal{PT}$ symmetry

As we will see, the $\mathcal{PT}$ symmetry gives rise to a non-trivial band topology in our setting that has effects on the nature of the band crossings.

Suppose we have two bands that are separated from other bands, e.g., when we are near a simple magic angle (cf. [BHZ22-2, Theorem 2]). Following the standard construction we define a rank two complex vector bundle $\mathcal{E}_{0}$ over the torus ${\mathbb{C}}/\Lambda^{*}$ :

\begin{gathered}\mathcal{E}_{0}:=\left\{[k,\phi]_{\tau}\in(\mathbb{C}\times L^{2}_{0}(\mathbb{C}/\Lambda;\mathbb{C}^{4}))/\sim_{\tau}:\phi\in\mathbbm{1}_{E_{\pm 1}(\alpha,\lambda,k)}(H_{k}(\alpha,\lambda))\right\},\\ (k,\phi)\sim_{\tau}(k^{\prime},\phi^{\prime})\ \Longleftrightarrow\ \exists\,p\in\Lambda^{*},\ k^{\prime}=k+p,\ \ \phi^{\prime}=\tau(p)\phi,\end{gathered}

(5.2)

where $\tau(p)\phi=e^{i\langle p,z\rangle}\phi$ . We consider the real subbundle $\mathcal{E}\subset\mathcal{E}_{0}$ defined by

\mathcal{E}=\{\varphi\in\mathcal{E}_{0}:\mathcal{PT}\varphi=\varphi\}.

This is a rank two real vector bundle such that the inner product of $\mathcal{E}_{0}$ restricted to $\mathcal{E}$ is real:

\langle\varphi_{1},\varphi_{2}\rangle=\overline{\langle\mathcal{PT}\varphi_{1},\mathcal{PT}\varphi_{2}\rangle}=\overline{\langle\varphi_{1},\varphi_{2}\rangle}\text{ for all }\varphi_{1},\varphi_{2}\in\mathcal{E}.

(5.3)

At the magic angle $\alpha_{0}\in\mathcal{A}$ of the chiral Hamiltonian $H_{0}(\alpha_{0},0)$ , we have flat bands and $\mathcal{E}$ can be computed explicitly. Suppose $\varphi=(u,v)^{T}$ , then

\mathcal{PT}\varphi=\varphi\Longleftrightarrow\mathscr{Q}u=v.

(5.4)

So $\mathcal{E}$ is defined by the equation $(D(\alpha)+k)u=0$ and equation (5.4), which has Chern number $c_{1}(\mathcal{E})=-1$ by [BHZ22-2, Theorem 4]. As the top Chern class $c_{1}(\mathcal{E})$ equals the Euler class $e_{2}(\mathcal{E})$ of the corresponding real bundle, in a neighbourhood of $(\alpha_{0},0)$ , $\mathcal{E}$ is an oriented rank two real bundle with Euler number $-1$ .

Now suppose the two bands only touches at finitely many points $p_{1},\cdots,p_{\ell}$ . We may define a flat connection $\nabla$ on the bundle $\mathcal{E}$ outside the points $\{p_{j}\}$ by declaring that the normalized eigenfunctions $\phi_{\pm 1}(k)$ are related by parallel transport, i.e. $\nabla\phi_{\pm 1}(k)=0$ . Then by Proposition 5.1,

e_{2}(\mathcal{E})=\sum\limits_{j}{\rm ind}(\nabla,p_{j}).

It remains to compute the winding numbers near the Dirac points. For this we recall our Grushin problem (3.3) with $\mathcal{PT}$ symmetry on $\mathcal{H}=L^{2}_{0}$ . There is also a compatible $\mathcal{PT}$ symmetry on ${\mathbb{C}}^{2}$ given by

\mathcal{PT}\begin{pmatrix}a\\ b\end{pmatrix}=\begin{pmatrix}\bar{b}\\ \bar{a}\end{pmatrix}.

Suppose we are in the setting of Proposition 3.1 and we have a conic singularity at the Dirac point, i.e., $c=\langle A\varphi,\psi\rangle\neq 0$ . Then we claim the winding number is $-1/2$ (see [BeCo18, Appendix] for a different argument).

Using conventions in the proof of Proposition 3.1, recall the Schur’s complement formula (cf. [TaZw23, Lemma 2.10])

(H_{k}-z)^{-1}=F(z)-F_{+}(z)F_{-+}(z)^{-1}F_{-}(z).

We first compute the winding number at a Dirac cone. When there are Dirac cones,

F_{-+}(z)=z-\begin{pmatrix}0&\bar{c}\bar{k}\\ ck&0\end{pmatrix}+\mathcal{O}(|k|^{2})

and

F_{-+}(z)^{-1}=(z-E_{1}(k))^{-1}\Pi_{v_{1}(k)}+(z-E_{-1}(k))^{-1}\Pi_{v_{-1}(k)}+{\rm holomorphic\ terms},

where

v_{1}(k)=\frac{1}{\sqrt{2}}\begin{pmatrix}\left(\frac{k}{|k|}\right)^{-1/2}\\ \left(\frac{k}{|k|}\right)^{1/2}\end{pmatrix}+\mathcal{O}(|k|),\quad v_{-1}(k)=\frac{i}{\sqrt{2}}\begin{pmatrix}-\left(\frac{k}{|k|}\right)^{-1/2}\\ \left(\frac{k}{|k|}\right)^{1/2}\end{pmatrix}+\mathcal{O}(|k|)

are normalized eigenvectors of $F_{-+}(z)$ that are invariant under $\mathcal{PT}$ symmetry, and $\Pi_{v_{\pm 1}(k)}$ are projections to the corresponding eigenvectors. Moreover, by the proof of [TaZw23, Proposition 2.12] we have

F_{+}(z)=E_{+}(z)+\mathcal{O}(|k|),\quad F_{-}(z)=E_{-}(z)+\mathcal{O}(|k|),

where $E_{\pm}(z)=R_{\mp}$ is given by equations (3.2) and (3.6). Suppose the two bands corresponding to $v_{\pm 1}(k)$ are separated from other bands and $\Pi$ is the projection to the two bands. Then

\Pi=-\frac{1}{2\pi i}\oint_{\gamma}F_{+}(z)F_{-+}(z)^{-1}F_{-}(z)dz=\Pi_{v_{1}(k)}+\Pi_{v_{-1}(k)}+\mathcal{O}(|k|),

where $\gamma$ is a curve that encloses $z=0$ away from other bands, so that the true eigenfunctions $\phi_{\pm 1}(k)$ corresponding to the two bands are given by $\phi_{\pm 1}(k)=v_{\pm 1}(k)\cdot(\varphi,\psi)+\mathcal{O}(|k|)$ in a neighbourhood of $k=0$ . Let $\varphi(k)$ and $\psi(k)$ be an orthonormal basis of $\mathbbm{1}_{E_{\pm 1}}(H_{k}(\alpha,\lambda))$ near $k=0$ with $\varphi(0)=\varphi$ and $\psi(0)=\psi$ . Then we can write

\begin{pmatrix}\phi_{1}(k)\\ \phi_{-1}(k)\end{pmatrix}=\frac{1}{\sqrt{2}}\begin{pmatrix}\left(\frac{k}{|k|}\right)^{-1/2}&\left(\frac{k}{|k|}\right)^{1/2}\\ -i\left(\frac{k}{|k|}\right)^{-1/2}&i\left(\frac{k}{|k|}\right)^{1/2}\end{pmatrix}\begin{pmatrix}\varphi(k)\\ \psi(k)\end{pmatrix}+\mathcal{O}(|k|).

Recall $\varphi=\mathcal{PT}(\psi)$ , we can take

\varphi(k)=\mathcal{PT}(\psi(k)),\quad\varphi_{+}(k)=\frac{1}{\sqrt{2}}(\varphi(k)+\psi(k)),\quad\varphi_{-}(k)=\frac{i}{\sqrt{2}}(\psi(k)-\varphi(k))

so that $\mathcal{PT}(\varphi_{\pm}(k))=\varphi_{\pm}(k)$ . Then for $k=|k|e^{i\theta}$ ,

\begin{pmatrix}\phi_{1}(k)\\ \phi_{-1}(k)\end{pmatrix}=\begin{pmatrix}\cos(\theta/2)&\sin(\theta/2)\\ -\sin(\theta/2)&\cos(\theta/2)\end{pmatrix}\begin{pmatrix}\varphi_{+}(k)\\ \varphi_{-}(k)\end{pmatrix}+\mathcal{O}(|k|).

In other words,

\begin{pmatrix}\varphi_{+}(k)\\ \varphi_{-}(k)\end{pmatrix}=\begin{pmatrix}\cos(\omega(\theta))&-\sin(\omega(\theta))\\ \sin(\omega(\theta))&\cos(\omega(\theta))\end{pmatrix}\begin{pmatrix}\phi_{1}(k)\\ \phi_{-1}(k)\end{pmatrix},\quad\omega(\theta)=\theta/2+\mathcal{O}(|k|).

Since $H_{k}\phi_{\pm 1}(k)=E_{\pm 1}(k)\phi_{\pm 1}(k)$ , $\phi_{1}|_{\theta=2\pi}$ has to become $\pm\phi_{1}|_{\theta=0}$ after rotation around once. Therefore, $\omega(2\pi)-\omega(0)\equiv 0\mod\pi$ (this can also be seen from (5.1)). Hence $\omega(2\pi)-\omega(0)=\pi$ and the winding number is given by

a=\frac{1}{2\pi}\oint\langle\nabla\varphi_{+}(k),\varphi_{-}(k)\rangle=-\frac{1}{2\pi}(\omega(2\pi)-\omega(0))=-\frac{1}{2}.

Similarly, when there is a quadratic band touching, i.e.,

F_{-+}(k)=z-\begin{pmatrix}a|k|^{2}&bk^{2}\\ \bar{b}\bar{k}^{2}&a|k|^{2}\end{pmatrix}+\mathcal{O}(|k|^{3}).

When $b=0$ , $a\neq 0$ , the eigenfunctions $v_{\pm 1}(k)$ are given by $(1,0)^{T}+\mathcal{O}(|k|)$ and $(0,1)^{T}+\mathcal{O}(|k|)$ . So the winding number would be zero. When $b\neq 0$ , the eigenfunctions are given by

\begin{pmatrix}k/|k|\\ (k/|k|)^{-1}\end{pmatrix}+\mathcal{O}(|k|),\quad\begin{pmatrix}k/|k|\\ -(k/|k|)^{-1}\end{pmatrix}+\mathcal{O}(|k|).

The phase of $k/|k|$ would change by $2\pi$ , so the winding number is $1$ . We conclude the following.

Proposition 5.2.

Suppose the two bands $E_{\pm 1}(k)$ are isolated from other bands and exhibit quadratic crossing at points $k=0,-K$ in the sense that, following (3.12),

\langle A\varphi,\psi\rangle=0,\quad\langle(\Pi^{\perp}H_{0}\Pi^{\perp})^{-1}A^{*}\varphi,A\psi\rangle\neq 0.

Then the two bands $E_{\pm 1}(k)$ have to touch at some addition points. Moreover if they touch at a discrete set of points, then the sum of the winding numbers of the touching points outside $\{0,-K\}$ would be $-3$ or $-1$ .

5.3. Wannier basis

The nontrivial topology of the band will also give obstructions to the existence of exponential localized Wannier functions. While the Chern number of the bands around zero may be zero, the non-zero Euler number still affects the nature of the associated Wannier function, as also observed in [APY19]. We have the following

Proposition 5.3.

Suppose the two bands $E_{\pm 1}(k)$ are separated from other bands. Then there does not exist exponential localized Wannier functions that is invariant under $\mathcal{PT}$ symmetry. More precisely, there does not exist an orthonormal family $\{\varphi_{\gamma}\}$ inside $\Pi_{E_{\pm 1}(k)}L^{2}({\mathbb{C}})$ of the form

\varphi_{\gamma}(z)=\mathscr{L}_{\gamma}\varphi_{0}(z),\quad\mathcal{PT}\varphi_{\gamma}=\varphi_{-\gamma},\quad\int_{{\mathbb{C}}}|z|^{2}|\varphi_{0}(z)|^{2}dz<\infty.

(5.5)

Proof.

This follows from [TaZw23, Theorem 9], see also [P07] and references therein. We first observe the existence of Wannier basis satisfying (5.5) implies there exists an $H^{1}$ orthonormal trivialization of the bundle $\mathcal{E}$ . This follows from

s(k,z)=\sum\limits_{\gamma\in\Lambda}e^{i\langle z+\gamma,k\rangle}\mathscr{L}_{\gamma}\varphi_{0}(z).

One easily checks that $s$ is a unitary section satisfying $\mathcal{PT}s=s$ . Since

\int_{{\mathbb{C}}/\Lambda}\int_{{\mathbb{C}}/\Lambda^{*}}|\nabla_{k_{j}}s(k,z)|dkdz=C_{\Lambda}\int_{{\mathbb{C}}}|x_{j}\varphi_{0}(z)|^{2}dz,

the section we get is $H^{1}$ . Now think of the oriented rank $2$ real vector bundle $\mathcal{E}$ as a complex line bundle, the Euler number is the same as the Chern number, which is $-1$ . Since $s$ gives a unitary section of $\mathcal{E}$ , we get a contradiction from [TaZw23, Lemma 8.9]. ∎

References

[APY19] J. Ahn, S. Park, and B. Yang. Failure of Nielsen-Ninomiya theorem and fragile topology in two-dimensional systems with space-time inversion symmetry: application to twisted bilayer graphene at magic angle. Physical Review X 9.2 (2019): 021013.
[Am*20] H. Ammari, B. Fitzpatrick, E. O. Hiltunen, H. Lee, and S. Yu, Honeycomb-lattice Minnaert bubbles, SIAM Journal on Mathematical Analysis, 52(6),(2020) 5441-5466.
[Be*22] S. Becker, M. Embree, J. Wittsten and M. Zworski, Mathematics of magic angles in a model of twisted bilayer graphene, to appear in Probability and Mathematical Physics.
[Be*23] S. Becker, T. Humbert, J. Wittsten, and M. Yang, Chiral limit of twisted trilayer graphene, arXiv:2308.10859.
[BHZ22-1] S. Becker, T. Humbert, and M. Zworski, Integrability in chiral model of magic angles, Communications in Mathematical Physics 403.2 (2023): 1153-1169.
[BHZ22-2] S. Becker, T. Humbert, and M. Zworski, Fine structure of flat bands in a chiral model of magic angles, arXiv:2208.01628.
[BHZ23] S. Becker, T. Humbert, and M. Zworski, Degenerate flat bands in twisted bilayer graphene, arXiv:arXiv:2306.02909.
[BeZw23-1] S. Becker and M. Zworski, Dirac points for twisted bilayer graphene with in-plane magnetic field, arXiv:arXiv:2303.00743
[BeZw23-2] S. Becker and M. Zworski, From the chiral model of TBG to the Bistritzer–MacDonald model, Journal of Mathematical Physics 65.6 (2024).
[BeCo18] G. Berkolaiko and A. Comech, Symmetry and Dirac points in graphene spectrum, Journal of Spectral Theory 8.3 (2018): 1099-1147.
[BiMa11] R. Bistritzer and A. MacDonald, Moiré bands in twisted double-layer graphene. PNAS, 108, 12233–12237, 2011.
[Ca23] E. Cancès, L. Garrigue, and D. Gontier, Simple derivation of moiré-scale continuous models for twisted bilayer graphene. Phys. Rev. B, 107, 155403, 2023.
[Ca19] S. Carr, S. Fang, Z. Zhu, and E. Kaxiras, Exact continuum model for low-energy electronic states of twisted bilayer graphene. Phys. Rev. Res., 1, 013001, 2019.
[Ca18] S. Carr, D. Massatt, S. B. Torrisi, P. Cazeaux, M. Luskin, and E. Kaxiras, Relaxation and domain formation in incommensurate two-dimensional heterostructure. Phys. Rev. B, 98, 224102, 2018.
[Ca24] J. Cazalis. Dirac cones for a mean-field model of graphene. Pure and Applied Analysis 6.1 (2024): 129-185.
[CaWe21] M. Cassier, and M. I. Weinstein. High contrast elliptic operators in honeycomb structures, Multiscale Modeling & Simulation 19.4 (2021): 1784-1856.
[Ca20] P. Cazeaux, M. Luskin, and D. Massatt, Energy Minimization of Two Dimensional Incommensurate Heterostructure. Arch. Rat. Mech. Anal., 235, 1289-1325, 2020.
[ChWe24] J. Chaban, and M. I. Weinstein, Instability of quadratic band degeneracies and the emergence of Dirac points, arXiv:arXiv:2404.05886
[Dr21] A. Drouot, Ubiquity of conical points in topological insulators, Journal de l’École polytechnique—Mathématiques 8 (2021): 507-532.
[DuNo80] B.A. Dubrovin and S.P. Novikov, Ground states in a periodic field. Magnetic Bloch functions and vector bundles. Soviet Math. Dokl. 22, 1, 240–244, 1980.
[DyZw19] S. Dyatlov and M. Zworski, Mathematical theory of scattering resonances. Grad. Stud. Math. 200. Amer. Math. Soc., 2019.
[FLW16] C. L. Fefferman, J. P. Lee-Thorp, and M. I. Weinstein. Honeycomb Schrödinger operators in the strong binding regime, Communications on Pure and Applied Mathematics 71.6 (2018): 1178-1270.
[FW12] C. Fefferman and M. I. Weinstein, Honeycomb lattice potentials and Dirac points, Journal of the American Mathematical Society, Volume 25, Number 4, 1169–1220, 2012.
[GZ23] J. Galkowski, and M. Zworski, An abstract formulation of the flat band condition. arXiv:2307.04896.
[Gr09] V. Grushin. Multiparameter perturbation theory of Fredholm operators applied to Bloch functions, Math. Notes, 86(5-6):767–774, 2009.
[Ka55] T. Kato, Notes on Projections and Perturbation Theory, Technical Report No. 9, Univ. Calif., Berkeley, 1955.
[Ka13] T. Kato, Perturbation theory for linear operators, Springer Science & Business Media, 2013.
[Ko*24] T. Kong, D. Liu, M. Luskin, A. B. Watson, Modeling of Electronic Dynamics in Twisted Bilayer Graphene, SIAM Journal on Applied Mathematics, 84, 3, 1011-1038, 2024.
[Ku23] P. Kuchment. Analytic and algebraic properties of dispersion relations (Bloch varieties) and Fermi surfaces. What is known and unknown. Journal of Mathematical Physics 64.11 (2023).
[KrPa02] S. Krantz and H. Parks. A primer of real analytic functions. Springer Science & Business Media, 2002.
[KuPo07] P. Kuchment and O. Post. On the spectra of carbon nano-structures, Comm. Math. Phys., 275(3):805–826, 2007.
[KuPa08] K. Kurdyka and L. Paunescu, Hyperbolic polynomials and multiparameter real analytic perturbation theory, Duke Math. J. 141 (1), 123–149, 2008.
[Le16] M. Lee. Dirac cones for point scatterers on a honeycomb lattice, SIAM Journal on Mathematical Analysis 48.2 (2016): 1459-1488.
[LWZ19] J. P. Lee-Thorp, M. I. Weinstein, and Y. Zhu, Elliptic operators with honeycomb symmetry: Dirac points, edge states and applications to photonic graphene, Archive for Rational Mechanics and Analysis 232 (2019): 1-63.
[LLZ23] W. Li, J. Lin, and H. Zhang, Dirac points for the honeycomb lattice with impenetrable obstacles, SIAM Journal on Applied Mathematics 83.4 (2023): 1546-1571.
[Ma23] D. Massatt, S. Carr, and M. Luskin, Electronic Observables for Relaxed Bilayer Two-Dimensional Heterostructures in Momentum Space, Multiscale Modeling & Simulation 21(4), 1344-1378, 2023.
[P07] G. Panati, Triviality of Bloch and Bloch–Dirac Bundles Ann. Henri Poincaré 8, 995–1011, 2007.
[Q*24] X. Quan, A. B. Watson, D. Massatt, Construction and Accuracy of Electronic Continuum Models of Incommensurate Bilayer 2D Materials, arXiv:arXiv:2406.15712
[S77] B. Simon, Notes on infinite determinants of Hilbert space operators, Advances in Mathematics Volume 24, Issue 3, June 1977, Pages 244-273.
[SlWe58] J. C. Slonczewski and P. R. Weiss. Band structure of graphite, Phys. Rev., 109(2):272–279, 1958.
[SjZw07] J. Sjöstrand and M. Zworski, Elementary linear algebra for advanced spectral problems, Annales de l’Institut Fourier, Volume 57 (2007) no. 7, pp. 2095-2141. doi : 10.5802/aif.2328.
[TaZw23] Z. Tao and M. Zworski, PDE methods in condensed matter physics, Lecture Notes, 2023,
https://math.berkeley.edu/~zworski/Notes_279.pdf.
[TKV19] G. Tarnopolsky, A.J. Kruchkov and A. Vishwanath, Origin of magic angles in twisted bilayer graphene, Phys. Rev. Lett. 122, 106405, 2019.
[Wa47] P. R. Wallace. The band theory of graphite, Phys. Rev.,71:622–634, 1947.
[Wa*22] A. B. Watson, T. Kong, A. H. MacDonald, and M. Luskin, Bistritzer-MacDonald dynamics in twisted bilayer graphene, J. Math. Phys. 64(2023), 031502.
[WaLu21] A. Watson and M. Luskin, Existence of the first magic angle for the chiral model of bilayer graphene, J. Math. Phys. 62(2021), 091502.
[Zw12] M. Zworski. Semiclassical analysis, Vol. 138. American Mathematical Society, 2022.
[Zw23] M. Zworski, M. Yang, and Z. Tao, Mathematical results on the chiral model of twisted bilayer graphene (with an appendix by Mengxuan Yang and Zhongkai Tao), Journal of Spectral Theory, to appear.

Dirac cones and magic angles in the Bistritzer–MacDonald TBG Hamiltonian

Abstract.

1. Introduction

Definition 1.1 (Simple Dirac cone).

Remark 1.

Theorem 1 (Generic existence of Dirac cones).

Theorem 2 (Persistence of chiral magic angles).

Remark 2.

Theorem 3.

Remark 3.

Structure of the paper

Acknowledgements

2. Symmetries of the Hamiltonian and protected states

2.1. BM Hamiltonian

Remark 4.

Remark 5.

2.2. Symmetries revisited

2.3. Additional symmetries

Proposition 2.1.

Proposition 2.2.

3. Existence of Dirac cones

Proposition 3.1.

Proof.

Remark 6.

Remark 7.

Proposition 3.2.

Proof.

Remark 8.

Remark 9.

Proposition 3.3.

Proof.

Example 1.

Definition 3.4.

Proposition 3.5.

Proof.

Proof of Theorem 1.

Corollary 3.6.

Proof.

4. Magic angles: from chiral limit to the BM Hamiltonian

Lemma 4.1.

Proof.

Proof of Theorem 2.

5. Magic angles in the BM Hamiltonian

5.1. Winding number and Euler number

Proposition 5.1.

Proof.

5.2. Topology of band with 𝒫​𝒯\mathcal{PT} symmetry

Proposition 5.2.

5.3. Wannier basis

Proposition 5.3.

Proof.

References

5.2. Topology of band with $\mathcal{PT}$ symmetry