Existence of the first magic angle for the chiral model of bilayer graphene

Watson, Alexander B.; Luskin, Mitchell

doi:10.1063/5.0054122

We consider the chiral model of twisted bilayer graphene introduced by Tarnopolsky, Kruchkov, and Vishwanath (TKV). TKV proved that for inverse twist angles α such that the effective Fermi velocity at the moiré K point vanishes, the chiral model has a perfectly flat band at zero energy over the whole Brillouin zone. By a formal expansion, TKV found that the Fermi velocity vanishes at α ≈ 0.586. In this work, we give a proof that the Fermi velocity vanishes for at least one α between 0.57 and 0.61 by rigorously justifying TKV’s formal expansion of the Fermi velocity over a sufficiently large interval of α values. The idea of the proof is to project the TKV Hamiltonian onto a finite-dimensional subspace and then expand the Fermi velocity in terms of explicitly computable linear combinations of modes in the subspace while controlling the error. The proof relies on two propositions whose proofs are computer-assisted, i.e., numerical computation together with worst-case estimates on the accumulation of round-off error, which show that round-off error cannot possibly change the conclusion of computation. The propositions give a bound below on the spectral gap of the projected Hamiltonian, an Hermitian 80 × 80 matrix whose spectrum is symmetric about 0, and verify that two real 18-th order polynomials, which approximate the numerator of the Fermi velocity, take values with a definite sign when evaluated at specific values of α. Together with TKV’s work, our result proves the existence of at least one perfectly flat band of the chiral model.

I. INTRODUCTION

A. Outline

Twisted bilayer graphene (TBG) is formed by stacking one layer of graphene on top of another in such a way that the Bravais lattices of the layers are twisted relative to each other. For generic twist angles, the lattices will be incommensurate, so the resulting structure will not be periodic. Bistritzer and MacDonald (BM)¹ introduced an approximate model (BM model) for the electronic states of TBG, which is periodic over the scale of the bilayer moiré pattern, where the twist angle enters as a parameter. Using this model, BM showed that the Fermi velocity, the velocity of electrons at the Fermi level, vanishes at particular twist angles known as “magic angles.” The largest of these angles, known as the first magic angle, is at θ ≈ 1.1°. Numerical computations on the BM model show the stronger result that at magic angles, the Bloch band of the BM model at zero energy is approximately flat over the whole Brillouin zone.^1,2 The flatness of the zero energy Bloch band is thought to be a critical ingredient for recently observed superconductivity of TBG,³ although the precise mechanism for superconductivity in TBG is not yet settled.

Aiming at a simplified model that explains the nearly flat band of TBG, Tarnopolsky, Kruchkov, and Vishwanath (TKV)⁴ introduced a simplification of the BM model, which has an additional “chiral” symmetry, known as the chiral model. TKV showed analytically that at magic angles (of the chiral model, still defined by vanishing of the Fermi velocity), the chiral model has exactly flat bands over the whole Brillouin zone. Using a formal perturbation theory (for the chiral model, the natural parameter is the reciprocal of the twist angle up to a constant), TKV derived approximate values for the magic angles of the chiral model. It is worth noting that the first magic angles of the chiral model and the BM model are nearby, but the higher magic angles are not very close.

Becker et al. introduced a spectral characterization of magic angles of the TKV model where the role of a non-normal operator is emphasized [the operator D^α appearing in (1)]. Using this characterization, they numerically computed precise values for the magic angles of the TKV model [see the discussion below (14)].⁵ In the same work, they also proved that the lowest band of the TKV model becomes exponentially close to flat even away from magic angles as the natural small parameter tends to zero. The same authors also investigated flat bands of the TKV model with more general interlayer coupling potentials and the spectrum of other special cases of the BM model.⁶

In this work, we study the chiral model introduced by TKV and consider the problem of rigorously proving the existence of the first magic angle. We do this by justifying the formal perturbation theory of TKV to make a rigorous expansion of the Fermi velocity to high enough order, and over a large enough parameter range, so that we can prove the existence of a zero. By numerically verifying that the resulting expansion attains a negative value and proving that the result continues to hold when the effect of round-off error is included (Proposition II.2), we obtain existence of the magic angle (Theorem II.2).

The proof of validity of the expansion is challenging because the reciprocal of the twist angle at the zero of the Fermi velocity is large relative to the spectral gap of the unperturbed Hamiltonian, which means that the magic angle falls outside of the interval of twist angles where the perturbation series for the Fermi velocity is obviously convergent. To overcome this difficulty, we start by representing the chiral model Hamiltonian in a basis that takes full advantage of model symmetries. Then, using a rigorous bound on the high frequency components of the error, we reduce the error analysis to analysis of the eigenvalues of the chiral model projected onto finitely many low frequencies. The final stage of the error analysis (Theorem II.1) is to prove a proposition about the eigenvalues of the projected chiral model by a numerical computation that we prove continues to hold when the accumulation of round-off error is considered (Proposition IV.5). We discuss the limitations of our methods and, in particular, whether our methods might be generalized to the more general settings considered by Becker et al.^5,6 in Remarks II.1–II.3.

B. Code availability

We have made code for the numerical computations used in our proofs available at github.com/abwats/magic_angle. We give references to specific scripts in the text.

II. STATEMENT OF RESULTS

A. Tarnopolsky–Kruchkov–Vishwanath’s chiral model

The chiral model, such as the Bistritzer–MacDonald model (B–M model) from which it is derived, is a formal continuum approximation to the atomistic tight-binding model of twisted bilayer graphene. The BM and chiral models aim to capture physics over the length-scale of the bilayer moiré pattern, which is, for small twist angles, much longer than the length-scale of the individual graphene layer lattices. Crucially, even when the graphene layers are incommensurate so that the bilayer is aperiodic on the atomistic scale, the chiral model and BM model are periodic (up to phases) with respect to the moiré lattice so that they can be analyzed via Bloch theory.

We define the moiré lattice to be the Bravais lattice

Λ = \{m_{1} a_{1} + m_{2} a_{2} : (m_{1}, m_{2}) \in Z^{2}\}

generated by the moiré lattice vectors

a_{1} = \frac{2 π}{3} (\begin{matrix} \sqrt{3}, 1 \end{matrix}), a_{2} = \frac{2 π}{3} (\begin{matrix} - \sqrt{3}, 1 \end{matrix}),

and we denote a fundamental cell of the moiré lattice by Ω. The moiré reciprocal lattice is the Bravais lattice

Λ^{*} = \{n_{1} b_{1} + n_{2} b_{2} : (n_{1}, n_{2}) \in Z^{2}\}

generated by the moiré reciprocal lattice vectors defined by a_i · b_j = 2πδ_ij, given explicitly by

b_{1} = \frac{1}{2} (\begin{matrix} \sqrt{3}, 3 \end{matrix}), b_{2} = \frac{1}{2} (\begin{matrix} - \sqrt{3}, 3 \end{matrix}) .

We define $q_{1} = (0, - 1)$ ⁠, which is the (re-scaled) difference of the K points (Dirac points) of each layer, and

q_{1} = (\begin{matrix} 0, - 1 \end{matrix}), q_{2} = q_{1} + b_{1} = \frac{1}{2} (\begin{matrix} \sqrt{3}, 1 \end{matrix}), q_{3} = q_{1} + b_{2} = \frac{1}{2} (\begin{matrix} - \sqrt{3}, 1 \end{matrix}) .

We write Ω^* for a fundamental cell of the moiré reciprocal lattice and refer to such a cell as the Brillouin zone.

Let $ϕ ≔ \frac{2 π}{3}$ ⁠. Tarnopolsky–Kruchkov–Vishwanath’s chiral Hamiltonian is defined as

H^{α} = (\begin{matrix} 0 & D^{α †} \\ D^{α} & 0 \end{matrix}), D^{α} = (\begin{matrix} - 2 i \bar{\partial} & α U (r) \\ α U (- r) & - 2 i \bar{\partial} \end{matrix}),

(1)

where $\bar{\partial} = \frac{1}{2} (\partial_{x} + i \partial_{y})$ ⁠, $U (r) = e^{- i q_{1} \cdot r} + e^{i ϕ} e^{- i q_{2} \cdot r} + e^{- i ϕ} e^{- i q_{3} \cdot r}$ ⁠, † denotes the adjoint (Hermitian transpose), and α is a real parameter, which we will take to be positive α ≥ 0 throughout [see (3)]. The chiral Hamiltonian H^α is an unbounded operator on $H = L^{2} (R^{2}; C^{4})$ with domain $H^{1} (R^{2}; C^{4})$ ⁠. We will write functions in $H$ as

ψ (r) = (\begin{matrix} ψ_{1}^{A} (r), ψ_{2}^{A} (r), ψ_{1}^{B} (r), ψ_{2}^{B} (r) \end{matrix}),

(2)

where $| ψ_{τ}^{σ} (r) |^{2}$ represents the electron density near to the K point (in momentum space) on sublattice σ and on layer τ. The diagonal terms of D^α arise from Taylor expanding the single layer graphene dispersion relation about the K point of each layer, while the off-diagonal terms of D^α couple the A and B sublattices of layers 1 and 2. The chiral model is identical to the BM model except that inter-layer coupling between sublattices of the same type is turned off in the chiral model. The precise form of the interlayer coupling potential U can be derived under quite general assumptions on the real space interlayer hopping.^1,7 The parameter α is, up to unimportant constants, the ratio

α \sim \frac{interlayer hopping strength between A and B sublattices}{twist angle} .

(3)

Although the limit α → 0 can be thought of as the limit of vanishing interlayer hopping strength at fixed twist, it is physically more interesting to view the limit as modeling increasing twist angle at a fixed interlayer hopping strength.

B. Rigorous justification of TKV’s formal expansion of the Fermi velocity and proof of existence of first magic angle

Bistritzer and MacDonald studied the effective Fermi velocity of electrons in twisted bilayer graphene modeled by the BM model and computed values of the twist angle such that the Fermi velocity vanishes, which they called “magic angles.” One can similarly define an effective Fermi velocity for the chiral model and refer to values of α such that the Fermi velocity vanishes as “magic angles” [although technically α is related to the reciprocal of the twist angle (3)].

TKV proved the remarkable result that, at magic angles, the chiral model has a perfectly flat Bloch band at zero energy. Let $L_{K}^{2}$ denote the L² space on a single moiré cell Ω with moiré K point Bloch boundary conditions. The starting point of TKV’s proof is an expression for the Fermi velocity as a function of α, v(α), as a functional of one of the Bloch eigenfunctions, $ψ^{α} \in L_{K}^{2}$ ⁠, of H^α,

v (α) ≔ \frac{| ⟨(ψ^{α *} (- r)| ψ^{α} (r)⟩ |}{| ⟨(ψ^{α}| ψ^{α}⟩ |},

(4)

where $⟨(.| .⟩$ denotes the $L_{K}^{2}$ inner product. We give precise definitions of $L_{K}^{2}$ ⁠, ψ^α, and v(α) in Definition III.2, Proposition III.6, and Definition III.3, respectively. We give a systematic formal derivation of why (4) is the effective Fermi velocity at the moiré K point in Appendix A. To complete the proof, TKV showed that zeros of v(α) imply zeros of ψ^α at special “stacking points” of Ω and that such zeros of ψ^α allow for Bloch eigenfunctions with zero energy to be constructed for all k in the moiré Brillouin zone.

To derive approximate values for magic angles, TKV computed a formal perturbation series approximation of ψ^α,

ψ^{α} = Ψ^{0} + α Ψ^{1} + \dots,

(5)

and then substituted this expression into the functional for v(α) to obtain an expansion of v(α) in powers of α,

v (α) = \frac{1 - 3 α^{2} + α^{4} - \frac{111}{49} α^{6} + \frac{143}{294} α^{8} + \dots}{1 + 3 α^{2} + 2 α^{4} + \frac{6}{7} α^{6} + \frac{107}{98} α^{8} + \dots} .

(6)

By setting v(α) = 0, one obtains an approximation for the smallest magic angle: α ≈ 0.586.

Although TKV proved that flat bands occur at magic angles, they did not prove the existence of magic angles, and hence, they did not prove the existence of flat bands. The contribution of the present work is to prove rigorous estimates on the error in approximation (5), which are sufficiently high order and precise that, once substituted into (4), they suffice to rigorously prove the existence of a zero of v(α) and, hence, via TKV’s proof, the existence of at least one perfectly flat band.

It turns out to be relatively straightforward to prove that series (5) and (6) are uniformly convergent and to derive precise bounds on the error in truncating the series for $| α | < \frac{1}{3}$ ⁠; see Proposition IV.3. The basic challenge, then, is to derive similar error bounds for α over an interval, which includes the expected location of the first magic angle, at $\approx \frac{1}{\sqrt{3}}$ ⁠. The first main theorem we will prove, roughly stated, is the following. See Theorem IV.1 for the more precise statement. The theorem relies on the existence of a spectral gap for an 80 × 80 Hermitian matrix, which requires numerical computation for its proof; see Proposition IV.5.

Theorem II.1.

The K point Bloch function ψ^α satisfies

ψ^{α} = \sum_{n = 0}^{8} α^{n} Ψ^{n} + η^{α},

(7)

where

η^{α} ⊥ \sum_{n = 0}^{8} α^{n} Ψ^{n}

with respect to the

L_{K}^{2}

inner product and

‖ η^{α} ‖_{L_{K}^{2}} \leq \frac{3 α^{9}}{15 - 20 α} f o r a l l 0 \leq α \leq \frac{7}{10} .

(8)

The functions Ψⁿ for 0 ≤ n ≤ 8 are derived recursively: see Appendix C. We stop at eighth order in the expansion because this is the minimal order such that we can guarantee the existence of a zero of v(α), but the functions Ψⁿ are well-defined by a recursive procedure for arbitrary positive integers n; see Proposition IV.1.

Substituting (7) into the functional for Fermi velocity (4) and using $η^{α} ⊥ \sum_{n = 1}^{8} α^{n} Ψ^{n}$ ⁠, we find that

v (α) = \frac{v_{N} (α)}{v_{D} (α)},

where

\begin{aligned} v_{N} (α) ≔ & ⟨(\sum_{n = 0}^{8} α^{n} Ψ^{n^{*}} (- r)| \sum_{n = 0}^{8} α^{n} Ψ^{n} (r)⟩ \\ + ⟨η^{α^{*}}, (- r), |\sum_{n = 0}^{8}), α^{n}, Ψ^{n}, (r)⟩ + ⟨(\sum_{n = 0}^{8} Ψ^{n^{*}} (- r)| η^{α} (r)⟩ \\ + ⟨(η^{α^{*}} (- r)| η^{α} (r)⟩ \end{aligned}

(9)

and

v_{D} (α) ≔ ⟨(\sum_{n = 0}^{8} α^{n} Ψ^{n}| \sum_{n = 0}^{8} α^{n} Ψ^{n}⟩ + ⟨(η^{α}| η^{α}⟩,

where $⟨(.| .⟩$ denotes the $L_{K}^{2}$ inner product and η^α satisfies (8). The following is a straightforward calculation.

Proposition II.1.

The following identities hold:

\begin{aligned} ⟨(\sum_{n = 0}^{8} α^{n} Ψ^{n^{*}} (- r)| \sum_{n = 0}^{8} α^{n} Ψ^{n} (r)⟩ \\ = 1 - 3 α^{2} + α^{4} - \frac{111}{49} α^{6} + \frac{143}{294} α^{8} - \frac{7 536 933}{11 957 764} α^{10} \\ + \frac{45 98 172 331}{47 460 365 316} α^{12} - \frac{30 028 809 212 865 451}{520 327 364 608 478 700} α^{14} + \frac{49 750 141 858 992 227}{12 487 856 750 603 488 800} α^{16}, \end{aligned}

(10)

\begin{aligned} ⟨(\sum_{n = 0}^{8} α^{n} Ψ^{n}| \sum_{n = 0}^{8} α^{n} Ψ^{n}⟩ \\ = 1 + 3 α^{2} + 2 α^{4} + \frac{6}{7} α^{6} + \frac{107}{98} α^{8} + \frac{5119}{48 412} α^{10} \\ + \frac{62 026 511}{356 844 852} α^{12} + \frac{355 691 470 247}{113 410 497 953 025} α^{14} + \frac{2 481 663 780 475 871}{337 509 641 908 202 400} α^{16} . \end{aligned}

(11)

We prove Proposition II.1 in Appendix F. Naïvely, expansions (10) and (11) approximate the formal infinite series expansions of $⟨(\sum_{n = 0}^{\infty} α^{n} Ψ^{n^{*}} (- r)| \sum_{n = 0}^{\infty} α^{n} Ψ^{n} (r)⟩$ and $⟨(\sum_{n = 0}^{\infty} α^{n} Ψ^{n}| \sum_{n = 0}^{\infty} α^{n} Ψ^{n}⟩$ up to terms of order α⁹. We prove in Proposition F.2 that because of some simplifications, expansions (10) and (11) agree with the infinite series up to terms of order α¹⁰.

We are now in a position to state and prove our second result. This result also relies on a proposition, which requires numerical computation for its proof: that one real 18-th order polynomial in α attains a negative value and another attains a positive value, when evaluated at specific values of α; see Proposition II.2.

Theorem II.2.

There exist positive numbers α_min and α_max with 0.57 < α_min < α_max < 0.61 such that the Fermi velocity v(α) defined by (4) has a zero α^* satisfying α_min ≤ α^* ≤ α_max.

Proof.

Equation () and Cauchy and Schwarz imply that

|v_{N} (α) - ⟨(\sum_{n = 0}^{8} α^{n} Ψ^{n^{*}} (- r)| \sum_{n = 0}^{8} α^{n} Ψ^{n} (r)⟩| \leq 2 ‖ η^{α} ‖ \sum_{n = 0}^{8} α^{n} ‖Ψ^{n}‖ + ‖ η^{α} ‖^{2} .

Using Theorem II.1 and Proposition C.3, we see that

v_{N} (α)

is bounded above by the polynomial

\begin{aligned} 1 - 3 α^{2} + α^{4} - \frac{111}{49} α^{6} + \frac{143}{294} α^{8} - \frac{7 536 933}{11 957 764} α^{10} \\ + \frac{4 598 172 331}{47 460 365 316} α^{12} - \frac{30 028 809 212 865 451}{520 327 364 608 478 700} α^{14} + \frac{49 750 141 858 992 227}{12 487 856 750 603 488 800} α^{16} \\ + E (α), \end{aligned}

(12)

where

\begin{aligned} E (α) ≔ & \frac{6 α^{9}}{15 - 20 α} (1 + \sqrt{3} α + \sqrt{2} α^{2} + \frac{\sqrt{14}}{7} α^{3} + \frac{\sqrt{258}}{42} α^{4} + \frac{\sqrt{1 968 837}}{3458} α^{5}) \\ (+ \frac{\sqrt{106 525 799}}{31 122} α^{6} + \frac{2 \sqrt{2 129 312 323 981 473}}{624 696 345} α^{7} + \frac{\sqrt{183 643 119 755 214 454}}{4 997 570 760} α^{8}) \\ + \frac{9 α^{18}}{{(15 - 20 α)}^{2}}, \end{aligned}

where we use Proposition C.3 to calculate the term in brackets for all

0 \leq α \leq \frac{7}{10}

⁠. On the other hand, v(α) is bounded below for all

0 \leq α \leq \frac{7}{10}

by the polynomial

\begin{aligned} 1 - 3 α^{2} + α^{4} - \frac{111}{49} α^{6} + \frac{143}{294} α^{8} - \frac{7 536 933}{11 957 764} α^{10} \\ + \frac{4 598 172 331}{47 460 365 316} α^{12} - \frac{30 028 809 212 865 451}{520 327 364 608 478 700} α^{14} + \frac{49 750 141 858 992 227}{12 487 856 750 603 488 800} α^{16} \\ - E (α) . \end{aligned}

(13)

We now claim the following.

Proposition II.2.

Expression (12), or equivalently the 18-th order polynomial obtained by multiplying (12) by (15 − 20α)², is negative at α = 0.61. Similarly, expression (13) is positive at α = 0.57.

Proposition II.2 obviously implies by continuity (12), (13), and $v_{N} (α)$ ⁠, each has at least one zero in the interval 0.57 < α < 0.61. We denote the largest zero of (12) in the interval by α_max and the smallest zero of (13) in the interval by α_min. Since the zeroes of $v_{N} (α)$ must lie between those of (12) and (13), we are done.□

Proof of Proposition II.2

(computer-assisted). We will first prove that (12) attains a negative value at 0.61 and then explain the modifications necessary to prove that (13) is positive at 0.57. Evaluating using the double-precision floating-point arithmetic, we find that at α = 0.61, (12) attains the negative value −0.020 263 (five significant figures; this value was computed by running the script compute_expansion_symbolically.py in the Github repo). It is straightforward to bound the numerical error, which accumulates when evaluating an 18-th order polynomial using the floating-point arithmetic. Even the simplest exact bound, which does not account for error cancellation (see, e.g., Eq. (8) of Oliver),⁸ yields an upper bound on the possible accumulated round-off error in the evaluation of an nth order polynomial $\sum_{j = 0}^{n} p_{j} α^{j}$ ⁠, for α ∈ [−1, 1], as $(n + 1) [e^{(2 n + 1) ϵ} - 1] \sup_{0 \leq j \leq n} | p_{j} |$ ⁠, where ϵ is “machine epsilon”: roughly speaking, the maximum possible round-off error generated in a single arithmetic operation. Bounding the maximum coefficient in (12) by 1000, taking n = 18, and bounding ϵ by 3 × 10⁻¹⁶ (which was easily attained working in Python on our machine), the maximum possible numerical error in the evaluation of (12) is ≈10⁻¹¹, which is much smaller than 0.020 263. We conclude that the first claim of Proposition II.2 must hold. Regarding the second, evaluating at α = 0.57, we find that (13) equals 0.029 138 (5sf). The same argument as before now shows that accumulated round-off error in the evaluation cannot possibly change the sign of (13) at α = 0.57.□

We do not attempt to rigorously estimate α_min and α_max precisely in this work, but numerically computing roots of polynomials (12) and (13) suggests α_min ≈ 0.576 83 (5sf) and α_max ≈ 0.601 77 (5sf), respectively, where 5sf is an abbreviation for five significant figure. Numerical computation of the first zero of $⟨(\sum_{n = 0}^{8} α^{n} Ψ^{n *} (- r)| \sum_{n = 0}^{8} α^{n} Ψ^{n} (r)⟩$ gives 0.585 97 (5sf); see Fig. 1 (the zero values were computed by running the script compute_expansion_symbolically.py in the Github repo).

FIG. 1.

View large Download slide

Left: plot of the numerator $v_{N} (α)$ of the Fermi velocity approximated by eighth order TKV expansion (6) (orange) and of eighth order expansions with worst-case (12) (blue) and best-case (13) (green) errors. Right: detail showing computed roots of these functions near to $α = \frac{1}{\sqrt{3}}$ ⁠. Numerically computing the zeroes of each curve yields α = 0.585 97 (5sf), α = 0.601 77 (5sf), and α = 0.576 83 (5sf), respectively. The values of α (0.57 and 0.61) where we evaluated expressions (13) (green line) and (12) (blue line) to prove that $v_{N} (α)$ has a zero between 0.57 and 0.61 are shown with the black crosses.

Using Proposition C.1 and the package Sympy⁹ for symbolic computation, we can compute the formal expansion of v(α) up to arbitrarily high order in α. In particular, we find the higher-order terms in expansion (6) to be

v (α) = \frac{1 - 3 α^{2} + α^{4} - \frac{111}{49} α^{6} + \frac{143}{294} α^{8} - \frac{10 227 257}{11 957 764} α^{10} + \frac{6 881 137 015}{47 460 365 316} α^{12} - \frac{130 055 941 435 858 531}{520 327 364 608 478 700} α^{14} + \dots}{1 + 3 α^{2} + 2 α^{4} + \frac{6}{7} α^{6} + \frac{107}{98} α^{8} + \frac{16 011}{48 412} α^{10} + \frac{134 058 653}{356 844 852} α^{12} + \frac{26 407 145 691 649}{226 820 995 906 050} α^{14} + \dots} .

(14)

Truncating the numerator after order α⁴⁰ and setting the numerator equal to zero yield α = 0.585 663 558 389 56 (14sf) for the first zero of the Fermi velocity (to compute this value, run compute_expansion_symbolically.py in the Github repo with N = 40). This is consistent with the numerical computation of Becker et al.,⁵ who found α = 0.585 663 558 389 55 truncated (not rounded) to 14 digits, by diagonalizing a non-normal but compact operator whose reciprocal eigenvalues correspond to magic angles. Note that we do not attempt to rigorously justify series (14) to such large values of α and to such high order in this work; see Remark II.4.

Remark II.1

(Higher magic angles). The chiral model has been conjectured to have infinitely many magic angles,⁵ but it is not straightforward to extend our methods to prove the existence of such higher magic angles. The problem is that calculating the perturbation series centered at α = 0 requires diagonalizing the unperturbed operator H⁰. In principle, it might be possible to calculate the perturbation series to higher order in order to get an accurate approximation of the Fermi velocity near to the higher magic angles. However, this would require significantly more calculation compared with the present work, and we have no guarantee that the error can be made small enough to prove the existence of another zero in that case.

Remark II.2

(More general interlayer hopping potentials). The chiral model (1) is an approximation to the full Hamiltonian of the twisted bilayer, even in the chiral limit where coupling between sublattices of the same type is turned off because the interlayer hopping potential U only allows for hopping between nearest neighbors in the momentum lattice (see Fig. 6 ). More general interlayer hopping potentials have been studied by Becker et al.⁶ In principle, such models should be amenable to the analysis of this work, but longer-range hopping would lead to much more involved calculations, and the construction of the finite-dimensional subspace Ξ of Proposition IV.4 would require more care: the fact that we can choose Ξ so that $‖ P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ = 1$ depends on H¹ only coupling nearest neighbors in the momentum lattice. Locality of hopping in the momentum space lattice has been exploited for efficient computation of density of states¹⁰ of twisted bilayers.

Remark II.3

(Generalization to BM model). Parts of our analysis should also apply to the full Bistritzer–MacDonald model. Specifically, one could study perturbation series for Bloch functions near to zero energy in powers of the inter-layer hopping strength, derive an equivalent expression for the Fermi velocity in terms of that series, and then study the zeroes of that series. However, there are various complications because of the lack of “chiral” symmetry. First, there is no reason for the continuation of the zero eigenvalue of the unperturbed operator to remain at zero. Second, the expression for the Fermi velocity in terms of the associated eigenfunction could be more complicated. Since zeros of the BM model Fermi velocity do not imply the existence of flat bands for that model, we do not consider these complications in this work.

Remark II.4

(Expanding to higher order). Our methods could, in principle, be continued to justify the expansion of the Fermi velocity to arbitrarily high order and potentially over larger intervals of α values. However, these extensions are not immediate: pushing the expansion to higher order or to a larger interval of α values would require a larger set Ξ in Lemma IV.1, and Proposition IV.5 would have to be re-proved for the new set Ξ. Note that the essential difficulty is justifying the perturbation series for large α: the series are easily justified to all orders for $| α | < \frac{1}{3}$ ⁠; see Proposition IV.3.

C. Structure of paper

We review the symmetries, Bloch theory, and symmetry-protected zero modes of TKV’s chiral model in Sec. III. We prove Theorem II.1 in Sec. IV, postponing most details of the proofs to Appendixes A– F. In Appendix A, we show why (4) corresponds to the effective Fermi velocity at the moiré K point. In Appendix B, we construct an orthonormal basis, which we refer to as the chiral basis, which allows for efficient computation and analysis of TKV’s formal expansion. We re-derive TKV’s formal expansions in Appendix C. We give details of the Proof of Theorem II.1 in Appendixes D and E. We prove Proposition II.1 in Appendix F. In the supplementary material, we list the basis functions of the subspace onto which we project the TKV Hamiltonian, give the explicit forms of the higher-order corrections in expansion (7), and present a derivation of the TKV Hamiltonian from the Bistritzer–MacDonald model.

III. SYMMETRIES, BLOCH THEORY, AND ZERO MODES OF TKV’s CHIRAL MODEL

A. Symmetries of the TKV model

In this section, we review the symmetries of the TKV model for the reader’s convenience and to fix notation. Becker et al.⁵ gave a group theoretical account of these symmetries, and further reviews can be found in the physics literature.^11–13 Recall that $ϕ = \frac{2 π}{3}$ ⁠, and let R_ϕ denote the matrix, which rotates vectors counter-clockwise by ϕ, i.e.,

R_{ϕ} = \frac{1}{2} (\begin{matrix} - 1 & - \sqrt{3} \\ \sqrt{3} & - 1 \end{matrix}) .

We define the following.

Definition III.1.

For any v ∈ Λ, we define a phase-shifted translation operator acting on functions

f \in H

by

τ_{v} f ≔ d i a g (1, e^{i q_{1} \cdot v}, 1, e^{i q_{1} \cdot v}) {\tilde{τ}}_{v} f, {\tilde{τ}}_{v} f (r) = f (r + v) .

(15)

We define a phase-shifted version of the operator, which rotates functions

f \in H

clockwise by ϕ by

R f ≔ d i a g (1, 1, e^{- i ϕ}, e^{- i ϕ}) \tilde{R} f, \tilde{R} f (r) = f (R_{ϕ} r) .

(16)

For any

f \in H

⁠, we finally define the “chiral” symmetry operator

S f ≔ d i a g (1, 1, - 1, - 1) f .

(17)

We then have the following.

Proposition III.1.

Operators (15) and (16) are symmetries in the sense that

[H^{α}, τ_{v}] = H^{α} τ_{v} - τ_{v} H^{α} = 0

(18)

for all moiré lattice vectors v ∈ Λ,

[H^{α}, R] = H^{α} R - R H^{α} = 0,

and operator (17) is a “chiral” symmetry in the sense that

\{H^{α}, S\} = H^{α} S + S H^{α} = 0 .

(19)

Proof.

The first claim is a direct calculation using the facts that for any v ∈ Λ,

{\tilde{τ}}_{- v} U (r) {\tilde{τ}}_{v} = e^{- i q_{1} \cdot v} U (r), {\tilde{τ}}_{- v} \bar{\partial} {\tilde{τ}}_{v} = \bar{\partial} .

The second claim is a direct calculation using the facts that

{\tilde{R}}^{- 1} U (r) \tilde{R} = e^{- i ϕ} U (r), {\tilde{R}}^{- 1} \bar{\partial} \tilde{R} = e^{- i ϕ} \bar{\partial} .

The final claim is trivial to check.□

The “chiral” symmetry (19) implies that the spectrum of H^α is symmetric about zero because

H^{α} ψ = E ψ \Leftrightarrow H^{α} S ψ = - E S ψ .

The same calculation implies that zero modes of H^α can always be chosen without loss of generality to be eigenfunctions of $S$ ⁠.

B. Bloch theory for the TKV Hamiltonian

We now want to reduce the eigenvalue problem for H^α using the symmetries just introduced. The symmetry (18) means that eigenfunctions of H^α can be chosen without loss of generality to be simultaneous eigenfunctions of τ_v for all v ∈ Λ. It therefore suffices to seek solutions of

H^{α} ψ = E ψ

for r in a fundamental cell $Ω ≔ R^{2} / Λ$ of the moiré lattice in the symmetry-restricted spaces

L_{k}^{2} ≔ \{f \in L^{2} (Ω; C^{4}) : f (r + v) = e^{i k \cdot v} d i a g (1, e^{i q_{1} \cdot v}, 1, e^{i q_{1} \cdot v}) f (r) \forall v \in Λ\},

(20)

where k is known as the quasimomentum. Since $L_{k + w}^{2} = L_{k}^{2}$ for any w ∈ Λ^*, it suffices to restrict attention to k in a fundamental cell of Λ^*, which we denote $Ω^{*} ≔ R^{2} / Λ^{*}$ and refer to as the Brillouin zone. We also define symmetry-restricted Sobolev spaces $H_{k}^{s}$ for each k ∈ Ω^* and positive integer s by

H_{k}^{s} ≔ \{f \in H^{s} (Ω; C^{4}) : f (r + v) = e^{i k \cdot v} d i a g (1, e^{i q_{1} \cdot v}, 1, e^{i q_{1} \cdot v}) f (r) \forall v \in Λ\} .

We claim the following.

Proposition III.2.

For each fixed k ∈ Ω^* and α ≥ 0, H^α, defined on the domain $H_{k}^{1}$ , extends to an unbounded self-adjoint elliptic operator $L_{k}^{2} \to L_{k}^{2}$ with compact resolvent. In a complex neighborhood of every α ≥ 0, the family H^α is a holomorphic family of type (A) in the sense of Kato.¹⁴

Proof.

Ellipticity is immediate since the principal symbol of H^α is invertible. Self-adjointness is clear using the Fourier transform when α = 0 and for α ≠ 0 because αH¹ is a bounded symmetric perturbation of H⁰ (see, e.g., Theorem 1.4 of Cycon et al.¹⁵). Elliptic regularity implies that the resolvent maps $L_{k}^{2} \to H_{k}^{1}$ ⁠, and compactness of the resolvent then follows by Rellich’s theorem (see, e.g., Proposition 3.4 of Taylor¹⁶). The family H^α is holomorphic of type (A) since the domain of H^α is independent of α, and H^αf is holomorphic for every $f \in H_{k}^{1}$ (see Chap. 7 of Kato¹⁴).□

We now claim the following.

Proposition III.3.

Let $f \in L_{k}^{2}$ . Then, $R f \in L_{R_{ϕ}^{*} k}^{2}$ .

Proof.

By definition, for any v ∈ Λ,

R f (r + v) = d i a g (1, 1, e^{- i ϕ}, e^{- i ϕ}) f (R_{ϕ} r + R_{ϕ} v) .

By the definition of

L_{k}^{2}

⁠, we have

R f (r + v) = e^{i (R_{ϕ}^{*} k) \cdot v} d i a g (1, e^{i (R_{ϕ}^{*} q_{1}) \cdot v}, 1, e^{i (R_{ϕ}^{*} q_{1}) \cdot v}) R f (r) .

The conclusion now follows from

R_{ϕ}^{*} q_{1} = q_{1} + b_{2}

and b₂ · v = 0 mod 2π for all v ∈ Λ.□

In particular, whenever $R_{ϕ}^{*} k = k$ mod Λ^*, we have $R L_{k}^{2} = L_{k}^{2}$ ⁠. Regarding such k, the following is a simple calculation.

Proposition III.4.

The moiré K and K′ points k = 0 and k = −q₁ and the moiré Γ point k = q₁ + b₁ satisfy $R_{ϕ}^{*} k = k$ mod Λ^*.

The moiré K, K′, and Γ points are shown in Fig. 2. Note that the moiré K, K′, and Γ points should not be confused with the single layer K, K′, and Γ points. The moiré K point corresponds to the K point of layer 1, while the moiré K′ point corresponds to the K point of layer 2. Interactions with the K′ points of layers 1 and 2 are formally small for small twist angles and are hence ignored.

FIG. 2.

View large Download slide

Diagram showing locations of moiré K (blue), K′ (red), and Γ (black) points within the moiré Brillouin zone (orange).

In this work, we will be particularly interested in Bloch functions at the moiré K and K′ points. We therefore define the following.

Definition III.2.

L_{K}^{2} ≔ L_{0}^{2}, L_{K^{'}}^{2} ≔ L_{- q_{1}}^{2} .

Let ω = e^iϕ. Since the spaces $L_{K}^{2}$ and $L_{K^{'}}^{2}$ are invariant under $R$ ⁠, they can be divided up into invariant subspaces corresponding to the eigenvalues of $R$ ⁠,

L_{K}^{2} = L_{K, 1}^{2} \oplus L_{K, ω}^{2} \oplus L_{K, ω^{*}}^{2}, L_{K^{'}}^{2} = L_{K^{'}, 1}^{2} \oplus L_{K^{'}, ω}^{2} \oplus L_{K^{'}, ω^{*}}^{2},

where

L_{K, σ}^{2} ≔ \{f \in L_{K}^{2} : R f = σ f\} σ = 1, ω, ω^{*}

and $L_{K^{'}, σ}^{2}, σ = 1, ω, ω^{*}$ ⁠, are defined similarly.

The following, which is trivial to prove, will be important for studying the zero modes of H^α.

Proposition III.5.

The operator $S$ commutes with τ_v and $R$ and hence maps the $L_{K, σ}^{2}$ and $L_{K^{'}, σ}^{2}$ spaces to themselves for σ = 1, ω, ω^*.

Since $S$ has eigenvalues ±1, we can define the spaces

L_{K, σ, \pm 1}^{2} = \{f \in L_{K, σ}^{2} : S f = \pm f\} σ = 1, ω, ω^{*}

and spaces $L_{K^{'}, σ, \pm 1}^{2}, σ = 1, ω, ω^{*}$ similarly.

Remark III.1.

Note that +1 eigenspaces of $S$ correspond to wave-functions, which vanish in their third and fourth entries, which correspond, through (2), to wave-functions supported only on A sites of the layers. Similarly, −1 eigenspaces of $S$ correspond to wave-functions, which vanish in their first and second entries, which are supported only on B sites of the layers.

C. Zero modes of the chiral model

We now want to investigate zero modes of H^α in detail. When α = 0, there are exactly four zero modes given by e_j, j = 1, 2, 3, 4, where e_j equals 1 in its jth entry and 0 in its other entries. It is easy to check that

e_{1} \in L_{K, 1}^{2}, e_{2} \in L_{K^{'}, 1}^{2}, e_{3} \in L_{K, ω^{*}}^{2}, e_{4} \in L_{K^{'}, ω^{*}}^{2},

(21)

and hence, 0 is a simple eigenvalue of H^α when restricted to each of these subspaces. Recall that zero modes can always be chosen as eigenfunctions of $S$ ⁠, and indeed, we have

e_{1} \in L_{K, 1, 1}^{2}, e_{2} \in L_{K^{'}, 1, 1}^{2}, e_{3} \in L_{K, ω^{*}, - 1}^{2}, e_{4} \in L_{K^{'}, ω^{*}, - 1}^{2} .

(22)

We now claim that these zero modes persist for all α. This was already established by TKV,⁴ and the following proposition is also similar to Proposition 3.1 of Becker et al.⁵ We re-state it using our notation and give the proof for completeness.

Proposition III.6.

There exist smooth functions ψ^α with ‖ψ^α‖ = 1 in each of the spaces $L_{K, 1, 1}^{2}$ , $L_{K^{'}, 1, 1}^{2}$ , $L_{K, ω^{*}, - 1}^{2}$ , and $L_{K^{'}, ω^{*}, - 1}^{2}$ such that ψ⁰ is as in (21), α ↦ ψ^α is real-analytic, and H^αψ^α = 0 for all α. The dimension of ker H^α restricted to each of the spaces $L_{K, 1}^{2}$ , $L_{K^{'}, 1}^{2}$ , $L_{K, ω^{*}}^{2}$ , and $L_{K, ω^{*}}^{2}$ is always odd-dimensional.

Proof.

Since $S$ preserves each of the spaces $L_{K, 1}^{2}$ ⁠, $L_{K^{'}, 1}^{2}$ ⁠, $L_{K, ω^{*}}^{2}$ ⁠, and $L_{K, ω^{*}}^{2}$ and anti-commutes with H^α, the spectrum of H^α restricted to each space must be symmetric about 0 for all α. Since H^α restricted to each space has compact resolvent and H^α is a holomorphic family of type (A), the spectrum of H^α consists of finitely degenerate isolated eigenvalues depending real-analytically on α, with associated eigenfunctions also real-analytically depending on α (although the real-analytic choice of eigenfunction at an eigenvalue crossing may not respect ordering); see Theorem 3.9 of Chap. 7 of Kato.¹⁴ The null space of H^α in each of the spaces is one-dimensional at α = 0 by explicit calculation, with the zero modes given by (up to non-zero constants) (21). For small α > 0, real-analyticity and the chiral symmetry force the null space to remain simple and it is clear how to define ψ^α. For large α > 0, the non-zero eigenvalues of H^α may cross 0 at isolated values of α, and in this case, we define ψ^α to be the real-analytic continuation of the zero mode through the crossings. Note that real-analyticity prevents non-zero eigenvalues from equaling zero except at isolated points so that the real-analytic continuation of the zero mode through the crossing must indeed be a zero mode. At crossings, the null space must be odd-dimensional in order to preserve symmetry of the spectrum of H^α about 0. It remains to check that if ψ⁰ is in, say, $L_{K, 1, 1}^{2}$ ⁠, then ψ^α must remain in $L_{K, 1, 1}^{2}$ for all α > 0. However, this must hold because the $S$ -eigenvalue of ψ^α cannot change abruptly while preserving real-analyticity. Smoothness of ψ^α follows from elliptic regularity.□

In this work, we will restrict attention to the moiré K point and, especially, the family $ψ^{α} \in L_{K, 1, 1}^{2}$ ⁠. We expect that our analysis would go through with only minor modifications if we considered instead the moiré K′ point. The zero modes in $L_{K, 1, 1}^{2}$ and $L_{K, ω^{*}, - 1}^{2}$ are related by the following symmetry.

Proposition III.7.

Let $ψ_{1}^{α}$ and $ψ_{- 1}^{α}$ denote the zero modes of H^α in the spaces $L_{K, 1, 1}^{2}$ and $L_{K, ω^{*}, - 1}^{2}$ ⁠, respectively. Then, $ψ_{1}^{α} = {(Φ^{α}, 0)}^{⊤}$ ⁠, where $Φ^{α} \in L^{2} (Ω; C^{2})$ , $Φ^{α} (r + v) = d i a g (1, e^{i q_{1} \cdot v}) Φ^{α} (r)$ for all v ∈ Λ, Φ^α(R_ϕr) = Φ^α(r). Up to gauge transformations $ψ_{- 1}^{α} \mapsto e^{i ϕ (α)} ψ_{- 1}^{α}$ ⁠, which preserve real-analyticity of $ψ_{- 1}^{α}$ , we have $ψ_{- 1}^{α} (r) = {(0, Φ^{α^{*}} (- r))}^{⊤}$ .

Proof.

Since $S ψ_{1}^{α} = ψ_{1}^{α}$ ⁠, the last two entries of $ψ_{1}^{α}$ must vanish, so we can write $ψ_{1}^{α} = {(Φ^{α}, 0)}^{⊤}$ ⁠. That Φ^α satisfies the stated symmetries follows immediately from $ψ_{1}^{α} \in L_{K, 1}^{2}$ ⁠. It is straightforward to check using the definitions of $R$ and τ_v that ${(0, Φ^{α^{*}} (- r))}^{⊤} \in L_{K, ω^{*}, - 1}^{2}$ ⁠. To see that ${(0, Φ^{α^{*}} (- r))}^{⊤}$ is a zero mode, note that Φ^α satisfies D^αΦ^α = 0, which implies that $D^{α †} Φ^{α^{*}} (- r) = 0$ by a simple manipulation. To see that $ψ_{- 1}^{α} (r) = {(0, Φ^{α^{*}} (- r))}^{⊤}$ (up to real-analytic gauge transformations) for all α, note first that this clearly holds for α = 0 [the zero modes are explicit (22)]. For α > 0, the identity must continue to hold by uniqueness (up to real-analytic gauge transformations) of the real-analytic continuation of $ψ_{- 1}^{α}$ starting from α = 0 and continuing first along the non-zero interval where $ψ_{1}^{α}$ is non-degenerate in $L_{K, 1, 1}^{2}$ and then through eigenvalue crossings as in the Proof of Proposition III.6.□

In Appendix A, we use Proposition III.7 to derive the effective Dirac operator with the α-dependent Fermi velocity, which controls the Bloch band structure in a neighborhood of the moiré K point. The Fermi velocity of the effective Dirac operator is given by the following. Note that we drop the subscript +1 when referring to the zero mode of H^α in $L_{K, 1, 1}^{2}$ since the zero mode of H^α in $L_{K, ω^{*}, - 1}^{2}$ plays no further role.

Definition III.3.

Let

ψ^{α} \in L_{K, 1, 1}^{2}

be as in Proposition III.6. Then, we define

v (α) ≔ \frac{| ⟨(ψ^{α *} (- r)| ψ^{α} (r)⟩ |}{| ⟨(ψ^{α}| ψ^{α}⟩ |},

(23)

where

⟨(.| .⟩

denotes the

L_{K}^{2}

inner product.

IV. RIGOROUS JUSTIFICATION OF TKV’s EXPANSION OF THE FERMI VELOCITY

A. Alternative formulation of TKV’s expansion

We now turn to approximating the zero mode $ψ^{α} \in L_{K, 1, 1}^{2}$ by a series expansion in powers of α. We write H^α = H⁰ + αH¹ and formally expand ψ^α as a series

ψ^{α} = Ψ^{0} + α Ψ^{1} + \dots,

(24)

where H⁰Ψ⁰ = 0 and

H^{0} Ψ^{n} = - H^{1} Ψ^{n - 1}

(25)

for all n ≥ 1. To solve H⁰Ψ⁰ = 0, we take Ψ⁰ = e₁. We prove the following in Appendix C.

Proposition IV.1.

Let P denote the projection operator in

L_{K, 1}^{2}

onto e₁ and P^⊥ = I − P. The sequence of Eq. (25) has a unique solution such that

Ψ^{n} \in L_{K, 1, 1}^{2}

for all n ≥ 0 and PΨⁿ = 0 for all n ≥ 1, given by Ψ⁰ = e₁ and

Ψ^{n} = - P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1} Ψ^{n - 1}

(26)

for each n ≥ 1.

Expansion (24) appears different from the series studied by TKV since we work only with the self-adjoint operators H⁰, H¹, and H^α rather than the non-self-adjoint operator D^α [defined in (1)]. Since functions in $L_{K, 1, 1}^{2}$ vanish in their last two components, there is no practical difference. However, working with only self-adjoint operators allows us to use the spectral theorem, which greatly simplifies the error analysis. We compute the first eight terms in expansion (24) in Proposition C.2 after developing some necessary machinery in Appendix B.

B. Rigorous error estimates for the expansion of the moiré K point Bloch function

In this section, we explain the essential challenge in proving error estimates for series (24) and explain how we overcome this challenge. Our goal is to prove the following.

Theorem IV.1.

Let

ψ^{α} \in L_{K, 1, 1}^{2}

be as in Proposition III.6. Then,

ψ^{α} = \sum_{n = 1}^{8} α^{n} Ψ^{n} + η^{α},

where

η^{α} ⊥ \sum_{n = 1}^{8} α^{n} Ψ^{n}

with respect to the

L_{K}^{2}

inner product and

‖ η^{α} ‖_{L_{K, 1}^{2}} \leq \frac{3 α^{9}}{15 - 20 α} f o r a l l 0 \leq α \leq \frac{7}{10} .

Proposition IV.1 guarantees that series (24) is well-defined up to arbitrarily many terms. A straightforward bound on the growth of terms in the series comes from the following proposition.

Proposition IV.2.

The spectrum of H⁰ in

L_{K, 1}^{2}

is

σ_{L_{K, 1}^{2}} (H^{0}) = {\pm | G |, \pm | q_{1} + G | : G \in Λ^{*}},

and hence,

‖ P^{⊥} {(H^{0})}^{- 1} P^{⊥} ‖_{L_{K, 1}^{2} \to L_{K, 1}^{2}} = 1 .

(27)

We also have

‖ H^{1} ‖_{L_{K, 1}^{2} \to L_{K, 1}^{2}} = 3 .

(28)

Proof.

This proposition is a combination of Propositions B.2, B.4, and B.7, proved in Appendix B.□

Proposition IV.2 implies that $‖ P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1} ‖_{L_{K, 1}^{2} \to L_{K, 1}^{2}} \leq 3$ ⁠, which implies the following.

Proposition IV.3.

The formal series (24) converges to ψ^α in $L_{K}^{2}$ , with an explicit error rate, for all $| α | < \frac{1}{3}$ . The formal series for the Fermi velocity v(α) obtained by substituting the series expansion of ψ^α into (23) converges for the same range of α, also with an explicit error rate.

Proof.

For any non-negative integer N, let $ψ^{N, α} ≔ \sum_{n = 0}^{N} α^{n} Ψ^{n}$ ⁠, where Ψⁿ are as in (26). Since Ψ⁰ ⊥ Ψⁿ for all n ≥ 1 and ‖Ψ⁰‖ = 1, we have that ‖ψ^N,α‖ ≥ 1 for all N. Let $ϕ^{N, α} ≔ \frac{ψ^{N, α}}{‖ ψ^{N, α} ‖}$ ⁠; then, we can decompose ϕ^N,α = cψ^α + η^α for some constant c and where η^α ⊥ ψ^α. Applying H^α to both sides, we have that $H^{α} ϕ^{N, α} = \frac{α^{N + 1} H^{1} Ψ^{N}}{‖ ψ^{N, α} ‖} = H^{α} η^{α}$ ⁠. Now, fix α ≥ 0 such that $| α | < \frac{1}{3}$ ⁠. Then, α‖H¹‖ < 1, and hence, the first non-zero eigenvalue of H^α is bounded away from 0 by 1 − 3α (recall that the first non-zero eigenvalues of H⁰ are ±1). Since η^α ⊥ ψ^α, where ψ^α spans the eigenspace of the zero eigenvalue of H^α, we have that $‖ η^{α} ‖ \leq \frac{| α^{N + 1} | ‖ H^{1} Ψ^{N} ‖}{| 1 - 3 α | ‖ ψ^{N, α} ‖}$ ⁠. Using the bound ‖Ψ^N‖ ≤ (3α)^N and the bound below on ‖ψ^N,α‖, we have that $‖ η^{α} ‖ \leq \frac{{(3 α)}^{N + 1}}{| 1 - 3 α |}$ ⁠, which clearly → 0 as N → ∞, so that lim_N→∞ϕ^N,α = ψ^α (up to a non-zero constant). Now, consider v(α) defined by (23). Assuming WLOG that ‖ψ^α‖ = 1 and substituting ψ^α = ϕ^N,α + η^α, we find immediately, using Cauchy–Schwarz, that $|v (α) - ⟨(ϕ^{N, *} (- r)| ϕ^{N} (r)⟩| \leq 2 ‖ η^{α} ‖ + ‖ η^{α} ‖^{2}$ ⁠. In terms of ψ^N, we have $|v (α) - \frac{⟨(ψ^{N, *} (- r)| ψ^{N} (r)⟩}{⟨(ψ^{N}| ψ^{N}⟩}| \leq 2 ‖ η^{α} ‖ + ‖ η^{α} ‖^{2}$ ⁠.□

Proposition IV.3 shows that for $| α | < \frac{1}{3}$ ⁠, series (24) converges to ψ^α and can be used to compute the Fermi velocity. However, this restriction is too strong to prove that the Fermi velocity has a zero, which occurs at the larger value $α \approx \frac{1}{\sqrt{3}}$ ⁠. Of course, Proposition IV.2 establishes only the most pessimistic possible bound on the expansion functions Ψⁿ, and this bound appears to be far from sharp from explicit calculation of each Ψⁿ; see Proposition C.3. We briefly discuss a possible route to a tighter bound in Remark C.2 but do not otherwise pursue this approach in this work.

We now explain how to obtain error estimates over a large enough range of α values to prove v(α) has a zero. We seek a solution of H^αψ^α = 0 in $L_{K, 1, 1}^{2}$ with the form

ψ^{α} = ψ^{N, α} + η^{α}, ψ^{N, α} ≔ \sum_{n = 0}^{N} α^{n} Ψ^{n} .

(29)

For arbitrary α, let Q^α denote the projection in $L_{K, 1}^{2}$ onto ψ^N,α, and Q^α,⊥ ≔ I − Q^α (note that Q⁰ = P). Note that Q^α depends on N, but we suppress this to avoid clutter. We assume WLOG that Q^αη^α(r) = 0. It follows that η^α satisfies

Q^{α, ⊥} H^{α} Q^{α, ⊥} η^{α} = - α^{N + 1} Q^{α, ⊥} H^{1} Ψ^{N} .

To obtain a bound on η^α in L²(Ω), we require a lower bound on the operator $Q^{α, ⊥} H^{α} Q^{α, ⊥} : Q^{α, ⊥} L_{K, 1}^{2} \to Q^{α, ⊥} L_{K, 1}^{2}$ ⁠. The following lemma gives a lower bound on this operator in terms of a lower bound on the projection of this operator onto the finite-dimensional subspace of $L_{K, 1}^{2}$ corresponding to a finite subset of the eigenfunctions of H⁰. The importance of this result is that since H¹ only couples finitely many modes of H⁰, for fixedN, by taking the subset sufficiently large, we can always arrange that ψ^N,α lies in this subspace.

Lemma IV.1.

Let P_Ξ denote the projection onto a subset Ξ of the eigenfunctions of H⁰ in

L_{K, 1}^{2}

, and let μ ≥ 0 be maximal such that

‖ P_{Ξ}^{⊥} H^{0} P_{Ξ}^{⊥} f ‖ \geq μ ‖ f ‖ \forall f \in H_{K, 1}^{1}, P_{Ξ}^{⊥} ≔ I - P_{Ξ},

(30)

(with this notation, the operator P introduced in Proposition IV.1 corresponds to P_Ξ with Ξ being the set {e₁} and μ = 1). Suppose that Q^αP_Ξ = P_ΞQ^α = Q^α, i.e., that ψ^N,α lies in ranP_Ξ. Define g^α by

\begin{aligned} g^{α} & ≔ \min \{| E | : \begin{matrix} E i s a n e i g e n v a l u e o f t h e m a t r i x Q^{α, ⊥} P_{Ξ} H^{α} P_{Ξ} Q^{α, ⊥} \\ a c t i n g Q^{α, ⊥} P_{Ξ} L_{K, 1}^{2} \to Q^{α, ⊥} P_{Ξ} L_{K, 1}^{2} \end{matrix}\} . \end{aligned}

We note that P_ΞQ^α,⊥ is the projection onto the subspace of

P_{Ξ} L_{K, 1}^{2}

orthogonal to ψ^N,α. As long as

3 α \leq μ a n d α ‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ < \min (g^{α}, μ - 3 α),

then

‖ Q^{α, ⊥} H^{α} Q^{α, ⊥} η^{α} ‖ \geq (\min (g^{α}, μ - 3 α) - α ‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖) ‖ Q^{α, ⊥} η^{α} ‖ .

(31)

Note that g^α would be identically zero if not for the restriction that the matrix acts on $Q^{α, ⊥} P_{Ξ} L_{K, 1}^{2}$ since otherwise ψ^N,α would be an eigenfunction with eigenvalue zero for all α. As it is, g⁰ = 1 and α ↦ g^α is real-analytic so that g^α must be positive for a non-zero interval of positive α values.

Proof.

Using Q^αP_Ξ = P_ΞQ^α, we have

P_{Ξ}^{⊥} Q^{α, ⊥} = Q^{α, ⊥} P_{Ξ}^{⊥} = P_{Ξ}^{⊥}

⁠, and hence,

\begin{aligned} ‖ Q^{α, ⊥} H^{α} Q^{α, ⊥} η^{α} ‖ = ‖ Q^{α, ⊥} (P_{Ξ} + P_{Ξ}^{⊥}) H^{α} (P_{Ξ} + P_{Ξ}^{⊥}) Q^{α, ⊥} η^{α} ‖ \\ = ‖ Q^{α, ⊥} P_{Ξ} H^{α} P_{Ξ} Q^{α, ⊥} η^{α} + α Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} η^{α} + α P_{Ξ}^{⊥} H^{α} P_{Ξ} Q^{α, ⊥} η^{α} + P_{Ξ}^{⊥} H^{α} P_{Ξ}^{⊥} η^{α} ‖ . \end{aligned}

By the reverse triangle inequality,

\begin{align} ‖ Q^{α, ⊥} H^{α} Q^{α, ⊥} η^{α} ‖ \\ \geq |‖ Q^{α, ⊥} P_{Ξ} H^{α} P_{Ξ} Q^{α, ⊥} η^{α} + P_{Ξ}^{⊥} H^{α} P_{Ξ}^{⊥} η^{α} ‖ - α ‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} η^{α} + P_{Ξ}^{⊥} H^{α} P_{Ξ} Q^{α, ⊥} η^{α} ‖| . \end{align}

(32)

We want to bound the second term above and the first term below. We start with the second term

\begin{aligned} ‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} η^{α} + P_{Ξ}^{⊥} H^{α} P_{Ξ} Q^{α, ⊥} η^{α} ‖^{2} \\ = ‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} η^{α} ‖^{2} + ‖ P_{Ξ}^{⊥} H^{1} P_{Ξ} Q^{α, ⊥} η^{α} ‖^{2} \\ \leq ‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖^{2} (‖ P_{Ξ}^{⊥} η^{α} ‖^{2} + ‖ P_{Ξ} Q^{α, ⊥} η^{α} ‖^{2}) \\ = ‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖^{2} ‖ Q^{α, ⊥} η^{α} ‖^{2}, \end{aligned}

where we use Pythagoras’ theorem,

P_{Ξ}^{⊥} H^{1} P_{Ξ} Q^{α, ⊥} η^{α} = P_{Ξ}^{⊥} H^{1} P_{Ξ} Q^{α, ⊥} P_{Ξ} Q^{α, ⊥} η^{α}

since P_ΞQ^α,⊥ is a projection, and

‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ = ‖ P_{Ξ}^{⊥} H^{1} P_{Ξ} Q^{α, ⊥} ‖

⁠. Hence, we can bound

‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} η^{α} + P_{Ξ}^{⊥} H^{α} P_{Ξ} Q^{α, ⊥} η^{α} ‖ \leq ‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ ‖ Q^{α, ⊥} η^{α} ‖ .

(33)

For the first term, first note that using Proposition IV.2 and the spectral theorem,

\begin{aligned} ‖ Q^{α, ⊥} P_{Ξ}^{⊥} H^{α} P_{Ξ}^{⊥} Q^{α, ⊥} η^{α} ‖ & \geq | ‖ Q^{α, ⊥} P_{Ξ}^{⊥} H^{0} P_{Ξ}^{⊥} Q^{α, ⊥} η^{α} ‖ - α ‖ Q^{α, ⊥} P_{Ξ}^{⊥} H^{1} P_{Ξ}^{⊥} Q^{α, ⊥} η^{α} ‖ | \\ \geq (μ - 3 α) ‖ P_{Ξ}^{⊥} Q^{α, ⊥} η^{α} ‖ \end{aligned}

as long as μ ≥ 3α. We now estimate

\begin{aligned} ‖ Q^{α, ⊥} P_{Ξ} H^{α} P_{Ξ} Q^{α, ⊥} η^{α} + P_{Ξ}^{⊥} H^{α} P_{Ξ}^{⊥} η^{α} ‖^{2} \\ = ‖ Q^{α, ⊥} P_{Ξ} H^{α} P_{Ξ} Q^{α, ⊥} η^{α} ‖^{2} + ‖ P_{Ξ}^{⊥} H^{α} P_{Ξ}^{⊥} η^{α} ‖^{2} \\ \geq {(g^{α})}^{2} ‖ Q^{α, ⊥} P_{Ξ} η^{α} ‖^{2} + {(μ - 3 α)}^{2} ‖ P_{Ξ}^{⊥} η^{α} ‖^{2} \\ \geq \min ({(g^{α})}^{2}, {(μ - 3 α)}^{2}) (‖ Q^{α, ⊥} P_{Ξ} η^{α} ‖^{2} + ‖ P_{Ξ}^{⊥} η^{α} ‖^{2}) \\ = \min ({(g^{α})}^{2}, {(μ - 3 α)}^{2}) ‖ Q^{α, ⊥} η^{α} ‖^{2} . \end{aligned}

It follows that as long as 3α ≤ μ,

‖ Q^{α, ⊥} P_{Ξ} H^{α} P_{Ξ} Q^{α, ⊥} η^{α} + P_{Ξ}^{⊥} H^{α} P_{Ξ}^{⊥} η^{α} ‖ \geq \min (g^{α}, μ - 3 α) ‖ Q^{α, ⊥} η^{α} ‖ .

(34)

The conclusion now holds as long as 3α ≤ μ and

α ‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ \leq \min (g^{α}, μ - 3 α)

upon substituting () and () into ().□

For Lemma IV.1 to be useful, we must check that it is possible to choose Ξ so that the bound (31) is non-trivial, i.e., so that the constant is positive. We will prove the following in Appendix D.

Proposition IV.4.

There exists a subset Ξ of the eigenfunctions of H⁰ such that the following holds:

The maximal μ such that (30) holds is μ = 7.
ψ^8,α defined by (29) lies in ran P_Ξ.
$‖ P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ = 1$ ⁠, and hence, $‖ Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ \leq 1$ .

The set Ξ constructed in Proposition IV.4 is the set of $L_{K, 1}^{2}$ -eigenfunctions of H⁰ with eigenvalues with magnitude $\leq 4 \sqrt{3}$ ⁠, augmented with two extra basis functions to ensure that $‖ P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ = 1$ ⁠. Including all $L_{K, 1}^{2}$ -eigenfunctions of H⁰ with eigenvalue magnitudes up to and including $4 \sqrt{3}$ ensure that ψ^8,α lies in ran P_Ξ.

We now require the following.

Proposition IV.5.

Let Ξ be as in Proposition IV.4. Then, $g^{α} \geq \frac{3}{4}$ for all $0 \leq α \leq \frac{7}{10}$ .

Proof

(computer-assisted). Consider $H_{Ξ}^{α} ≔ Q^{α, ⊥} P_{Ξ} H^{α} P_{Ξ} Q^{α, ⊥}$ acting on $P_{Ξ} L_{K, 1}^{2}$ ⁠. Assuming that α is restricted to an interval such that the zero eigenspace of $H_{Ξ}^{α}$ is simple, then, using orthogonality of eigenvectors corresponding to different eigenvalues and the fact that Q^α is the spectral projection onto the unique zero mode of $H_{Ξ}^{α}$ ⁠, $H_{Ξ}^{α}$ has the same non-zero eigenvalues as the matrix Q^α,⊥P_ΞH^αP_ΞQ^α,⊥ acting on $Q^{α, ⊥} P_{Ξ} L_{K, 1}^{2}$ ⁠. The matrix $H_{Ξ}^{α}$ is an 81 × 81 matrix whose spectrum is symmetric about 0 because of the chiral symmetry. When α = 0, the spectrum is explicit: 0 is a simple eigenvalue, and the smallest non-zero eigenvalues are ±1, both also simple. Proposition IV.5 is proved if we can prove that the first positive eigenvalue of $H_{Ξ}^{α}$ is bounded away from zero by $\frac{3}{4}$ for all $0 \leq α \leq \frac{7}{10}$ ⁠. Note that if this holds, the zero eigenspace of $H_{Ξ}^{α}$ must be simple for all $0 \leq α \leq \frac{7}{10}$ ⁠, and hence, our basic assumption is justified. The strategy of the proof is as follows.

Define a grid $G ≔ \{\frac{7 n}{10 N} : n \in {0, 1, \dots, N}\}$ ⁠, where N is a positive integer taken sufficiently large that the grid spacing $h ≔ \frac{7}{10 N} < \frac{1}{388 831}$ (the number 388 831 comes from Proposition IV.7).
Numerically compute the eigenvalues of $H_{Ξ}^{α}$ for $α \in G$ ⁠. We find that the numerically computed first positive eigenvalues of these matrices are uniformly bounded below by $\frac{8}{10} > \frac{3}{4}$ ⁠.
Perform a backward error analysis that fully accounts for round-off error in the numerical computation in order to prove that the exact first positive eigenvalues of the matrices $H_{Ξ}^{α}$ must also be bounded below by $\frac{8}{10}$ at each $α \in G$ ⁠.
Use perturbation theory to bound the exact first positive eigenvalue of $H_{Ξ}^{α}$ below by $\frac{3}{4}$ over the whole interval of α values between 0 and $\frac{7}{10}$ ⁠.

When discussing round-off error due to working in the floating-point arithmetic, we will denote “machine epsilon” by ϵ. The significance of this number is that we will assume that all complex numbers a can be represented by floating-point numbers $\tilde{a}$ such that $| a - \tilde{a} | \leq ϵ a$ ⁠. We will also make the standard assumption about creation of round-off error in the floating-point arithmetic operations: if $\tilde{a}$ and $\tilde{b}$ are floating-point complex numbers and if ${(\tilde{a} O \tilde{b})}_{c o m p}$ and $\tilde{a} O \tilde{b}$ represent the numerically computed value and exact value of an arithmetic operation on the numbers $\tilde{a}$ and $\tilde{b}$ ⁠, then ${(\tilde{a} O \tilde{b})}_{c o m p} = \tilde{a} O \tilde{b} + e$ ⁠, where $| e | \leq (\tilde{a} O \tilde{b}) ϵ$ ⁠. In Python, this is indeed the case, for all reasonably sized (such that stack overflow does not occur) complex numbers, with ϵ = 2.220 44 × 10⁻¹⁶ (5sf). We now present the main points of parts 2–4 of the strategy, postponing proofs of intermediate lemmas to Appendix E.

For part 2 of the strategy, for each $α \in G$ ⁠, we let ${\tilde{H}}_{Ξ}^{α}$ denote $H_{Ξ}^{α}$ (which is known exactly) evaluated as floating-point numbers. We generate numerically computed eigenpairs ${\tilde{λ}}_{j}, {\tilde{v}}_{j}$ for 1 ≤ j ≤ 81 for each ${\tilde{H}}_{Ξ}^{α}$ using numpy’s Hermitian eigensolver eigh. We find that the smallest first positive eigenvalue of ${\tilde{H}}_{Ξ}^{α}$ for $α \in G$ is 0.814 719 126 144 543 6 (computed using compute_PHalphaP_enclosures.py in the Github repo). Note that the difference between this number and $\frac{8}{10}$ is bounded below by 0.01.

The main tool for part 3 of the strategy is the following theorem.

Theorem IV.2.

Let m and n denote positive integers with m ≤ n. Let A be a Hermitian n × n matrix, and let

{v_{j}}_{1 \leq j \leq m}

be orthonormal n-vectors satisfying (A − λ_jI)v_j = r_j for scalars λ_j and n-vectors r_j for each 1 ≤ j ≤ m. Then, there are m eigenvalues

{α_{j}}_{1 \leq j \leq m}

of A, which can be put into one-to-one correspondence with λ_j’s such that

| λ_{j} - α_{j} | \leq 2 m \sup_{1 \leq i \leq m} ‖ r_{i} ‖_{2} f o r a l l 1 \leq j \leq m .

Proof.

See Appendix E 1.□

Naïvely, one would hope to be able to calculate enclosure intervals for every eigenvalue of $H_{Ξ}^{α}$ and, in particular, a lower bound on the first positive eigenvalue of $H_{Ξ}^{α}$ by directly applying Theorem IV.2 with $A = H_{Ξ}^{α}$ ⁠, m = 81, and λ_j and v_j given for each 1 ≤ j ≤ 81 by the approximate eigenpairs ${\tilde{λ}}_{j}, {\tilde{v}}_{j}$ computed in part 2. However, we cannot directly apply the theorem because ${{\tilde{v}}_{j}}_{1 \leq j \leq 81}$ are not exactly orthonormal because of round-off error. Hence, we will prove the existence of an exactly orthonormal set ${{\hat{v}}_{j}}_{1 \leq j \leq 81}$ close to the set ${{\tilde{v}}_{j}}_{1 \leq j \leq 81}$ and apply Theorem IV.2 to the set ${{\hat{v}}_{j}}_{1 \leq j \leq 81}$ (with the same ${\tilde{λ}}_{j}$ ⁠) instead. Note that to carry out this strategy, we must bound the residuals ${\hat{r}}_{j} ≔ (H_{Ξ}^{α} - {\tilde{λ}}_{j}) {\hat{v}}_{j}$ ⁠. The result we need to implement this strategy is the following. Note that the result requires numerical computation of inner products and residuals, and we account for round-off error in these computations.

Theorem IV.3.

Let m and n be positive integers with m ≤ n. Let A be an n × n Hermitian matrix, let

\tilde{A}

denote A evaluated in floating-point numbers, and let

{\tilde{v}}_{j}, {\tilde{λ}}_{j}

for 1 ≤ j ≤ m be a set of n-dimensional vectors and real numbers, respectively. Let

{⟨({\tilde{v}}_{i}| {\tilde{v}}_{j}⟩}_{c o m p}

denote their numerically computed inner products, and let

{\tilde{r}}_{j, c o m p} ≔ {[(\tilde{A} - {\tilde{λ}}_{j} I) {\tilde{v}}_{j}]}_{c o m p}

denote their numerically computed residuals. Let ϵ denote machine epsilon, and assume nϵ < 0.01. Let μ be

μ ≔ (1.01) n^{2} ϵ {(\sup_{1 \leq i \leq m} ‖ {\tilde{v}}_{i} ‖_{\infty})}^{2} + \sup_{1 \leq i \leq m} | {⟨({\tilde{v}}_{i}| {\tilde{v}}_{i}⟩}_{c o m p} - 1 | + \sup_{\begin{array}{c} i \neq j \\ 1 \leq i, j \leq m \end{array}} | {⟨({\tilde{v}}_{i}| {\tilde{v}}_{j}⟩}_{c o m p} | .

Then, as long as

m μ < \frac{1}{2}

, there is an orthonormal set of n-vectors

{{\hat{v}}_{j}}_{1 \leq j \leq m}

whose residuals

{\hat{r}}_{j} ≔ (A - {\tilde{λ}}_{j} I) {\hat{v}}_{j}

satisfy the bound

\begin{align} \sup_{1 \leq j \leq m} ‖ {\hat{r}}_{j} ‖_{2} \leq & 2^{- 1 / 2} n (‖ A ‖_{2} + \sup_{1 \leq j \leq m} | {\tilde{λ}}_{j} |) μ + n^{1 / 2} \sup_{1 \leq j \leq m} {‖{\tilde{r}}_{j, c o m p}‖}_{\infty} \\ + (1.01) n^{5 / 2} ϵ (‖ \tilde{A} ‖_{m a x} + \sup_{1 \leq j \leq m} | {\tilde{λ}}_{j} |) \sup_{1 \leq j \leq m} ‖ {\tilde{v}}_{j} ‖_{\infty} + n ϵ ‖ A ‖_{m a x} \sup_{1 \leq j \leq m} ‖ {\tilde{v}}_{j} ‖_{\infty}, \end{align}

(35)

where ‖A‖_max denotes the largest of the absolute values of the elements of the matrix A.

Proof.

See Appendix E 2.□

Numerical computation (using the script compute_PHalphaP_enclosures.py in the Github repo) shows that the maximum of $\sup_{1 \leq i \leq m} | {⟨({\tilde{v}}_{i}| {\tilde{v}}_{i}⟩}_{c o m p} - 1 |$ and $\sup_{\begin{array}{c} i \neq j \\ 1 \leq i, j \leq m \end{array}} | {⟨({\tilde{v}}_{i}| {\tilde{v}}_{j}⟩}_{c o m p} |$ over $α \in G$ is bounded by 7 × 10⁻¹⁵. Hence, we can apply Theorem IV.3 with $A = H_{Ξ}^{α}$ and ${\tilde{λ}}_{j}, {\tilde{v}}_{j}$ given by the numerically computed eigenpairs of ${\tilde{H}}_{Ξ}^{α}$ to obtain orthonormal sets ${{\hat{v}}_{j}}_{1 \leq j \leq 81}$ whose residuals with respect to $H_{Ξ}^{α}$ satisfy (35). The following is straightforward.

Proposition IV.6.

\sup_{0 \leq α \leq \frac{7}{10}} ‖ H_{Ξ}^{α} ‖_{2} \leq 10, \sup_{0 \leq α \leq \frac{7}{10}} ‖ H_{Ξ}^{α} ‖_{m a x} \leq 7 .

Proof.

The first estimate follows from ‖P_ΞH⁰P_Ξ‖ ≤ 7 and ‖H¹‖ ≤ 3. The second estimate follows immediately from writing the matrix $H_{Ξ}^{α}$ in the chiral basis.□

We can now apply Theorem IV.2 with $A = H_{Ξ}^{α}$ and λ_j, v_j given by the numerically computed ${\tilde{λ}}_{j}$ from part 2 and ${\hat{v}}_{j}$ coming from Theorem IV.3, in order to derive rigorous enclosure intervals for every eigenvalue of $H_{Ξ}^{α}$ ⁠. We find that (using the script compute_PHalphaP_enclosures.py in the Github repo) the suprema over $α \in G$ of $\sup_{1 \leq j \leq m} ‖ {\tilde{v}}_{j} ‖_{\infty}$ ⁠, $\sup_{1 \leq j \leq m} ‖ {\tilde{r}}_{j, c o m p} ‖_{\infty}, ‖ {\tilde{H}}_{Ξ}^{α} ‖_{\max}$ ⁠, and $\sup_{1 \leq j \leq m} | {\tilde{λ}}_{j} |$ are bounded by 1, 5 × 10⁻¹⁴, 7, and 8, respectively. It is then easy to see that 2 × 81 times the right-hand side of (35) is much smaller than 0.01 and is hence smaller than the distance between the minimum over $α \in G$ of the numerically computed first positive eigenvalues of ${\tilde{H}}_{Ξ}^{α}$ and $\frac{8}{10}$ ⁠. We can therefore conclude that the first positive eigenvalues of $H_{Ξ}^{α}$ are bounded below by $\frac{8}{10}$ at every $α \in G$ ⁠.

The main tool for part 4 of the strategy is the following.

Theorem IV.4.

Let A^α be an n × n Hermitian matrix real-analytically depending on a real parameter α. Denote the ordered eigenvalues of A^α by

λ_{j}^{α}

for 1 ≤ j ≤ n. Then, for any α and α₀,

| λ_{j}^{α} - λ_{j}^{α_{0}} | \leq | α - α_{0} | \sup_{β \in [α_{0}, α]} ‖ \partial_{β} A^{β} ‖_{2} f o r a l l 1 \leq j \leq n .

Proof.

See Appendix E 3.□

We would like to apply Theorem IV.4 to bound the variation of eigenvalues of $H_{Ξ}^{α}$ ⁠. To this end, we require the following proposition, which bounds the derivative of $H_{Ξ}^{α}$ with respect to α over the interval $0 \leq α \leq \frac{7}{10}$ ⁠.

Proposition IV.7.

\sup_{0 \leq α \leq \frac{7}{10}} {‖\partial_{α} H_{Ξ}^{α}‖}_{2} \leq 38 883 .

Proof.

See Appendix E 4.□

Proposition IV.7 combined with Theorem IV.4 explains the choice of distance $h = \frac{1}{388 831}$ between grid points. Assuming that the first positive eigenvalue of $H_{Ξ}^{α}$ is bounded below by $\frac{8}{10}$ at grid points between 0 and $\frac{7}{10}$ separated by h, we see that as long as

\frac{38 883 h}{2} < \frac{8}{10} - \frac{3}{4} ⟺ h < \frac{1}{388 830},

Proposition IV.7 and Theorem IV.4 guarantee that the first eigenvalue of $H_{Ξ}^{α}$ must be greater than $\frac{3}{4}$ over the whole interval $0 \leq α \leq \frac{7}{10}$ ⁠.□

Remark IV.1.

In the Proof of Proposition IV.5, we bound the round-off error, which can occur in our numerical computations in order to draw rigorous conclusions. A common approach to this is the interval arithmetic; see Rump¹⁷ and references therein. Our approach applies directly to the present problem and is just as rigorous.

The results of a computation of the eigenvalues of $H_{Ξ}^{α}$ are shown in Fig. 3.

FIG. 3.

View large Download slide

Plot of numerically computed eigenvalues of the 81 × 81 matrix $H_{Ξ}^{α}$ acting on $P_{Ξ} L_{K, 1}^{2}$ (blue lines), showing that the first non-zero eigenvalues are bounded away from 0 by $\frac{3}{4}$ (red lines) when α is less than $\frac{7}{10}$ (black line). The zero eigenvalue corresponds to the subspace spanned by ψ^8,α, and the non-zero eigenvalues equal those of the 80 × 80 matrix Q^α,⊥P_ΞH^αP_ΞQ^α,⊥ acting on $Q^{α, ⊥} P_{Ξ} L_{K, 1}^{2}$ since non-zero eigenvectors v of Q^α,⊥P_ΞH^αP_ΞQ^α,⊥ must be orthogonal to ψ^8,α by orthogonality of eigenvectors corresponding to different eigenvalues.

Assuming Proposition IV.4 and Proposition IV.5, the bound (31) becomes, for all $0 \leq α \leq \frac{7}{10}$ ⁠,

‖ Q^{α, ⊥} H^{α} Q^{α, ⊥} η^{α} ‖ \geq (\frac{3}{4} - α) ‖ Q^{α, ⊥} η^{α} ‖ .

We now assume the following, proved in Appendix C.

Proposition IV.8.

$‖ H^{1} Ψ^{8} ‖ \leq \frac{3}{20}$ .

We can now give the Proof of Theorem IV.1.

Proof of Theorem IV.1.

The proof follows immediately from Lemma IV.1, Proposition IV.4, Assumption IV.5, and Proposition IV.8.□

SUPPLEMENTARY MATERIAL

In the supplementary material, we list the chiral basis functions, which span the space Ξ, list terms Ψ⁵ − Ψ⁸ in the formal expansion of ψ^α, and derive the TKV Hamiltonian from the Bistritzer–MacDonald model.

ACKNOWLEDGMENTS

M.L. and A.B.W. were supported, in part, by ARO MURI Award No. W911NF-14-0247 and NSF DMREF Award No. 1922165.

DATA AVAILABILITY

The data that support the findings of this study are available within the article and its supplementary material.

APPENDIX A: DERIVATION OF EXPRESSION FOR FERMI VELOCITY IN TERMS OF $L_{K, 1, 1}^{2}$ ZERO MODE OF H^α

The Bloch eigenvalue problem for the TKV Hamiltonian at quasi-momentum k is

H^{α} ψ_{k}^{α} = E_{k} ψ_{k}^{α},

where H^α is as in (1) and

ψ_{k}^{α} (r + v) = e^{i k \cdot v} d i a g (1, e^{i q_{1} \cdot v}, 1, e^{i q_{1} \cdot v}) ψ_{k}^{α} (r) \forall v \in Λ .

By Propositions III.6 and III.7, 0 is a twofold (at least) degenerate eigenvalue at the moiré K point k = 0, with associated eigenfunctions $ψ_{\pm 1}^{α}$ as in Proposition III.7. In what follows, we assume that 0 is exactly twofold degenerate so that $ψ_{\pm 1}^{α}$ form a basis of the degenerate eigenspace. This assumption is clearly true for small α but could, in principle, be violated for α > 0.

Introducing $χ_{k}^{α} ≔ e^{- i k \cdot r} ψ_{k}^{α}$ ⁠, we derive the equivalent Bloch eigenvalue problem with k-independent boundary conditions

H_{k}^{α} χ_{k}^{α} = E_{k} χ_{k}^{α},

(A1)

where

H_{k}^{α} ≔ (\begin{matrix} 0 & D_{k}^{α †} \\ D_{k}^{α} & 0 \end{matrix}), D_{k}^{α} = (\begin{matrix} D_{x} + k_{x} + i (D_{y} + k_{y}) & α U (r) \\ α U (- r) & D_{x} + k_{x} + i (D_{y} + k_{y}) \end{matrix}),

where D_x,y ≔ −i∂_x,y and

χ_{k}^{α} (r + v) = d i a g (1, e^{i q_{1} \cdot v}, 1, e^{i q_{1} \cdot v}) χ_{k}^{α} (r) \forall v \in Λ .

Clearly, $ψ_{\pm 1}^{α}$ remain a basis of the zero eigenspace for the problem (A1) at k = 0.

Differentiating the operator $D_{k}^{α}$ ⁠, we find $\partial_{k_{x}} D_{k}^{α} = I_{2}$ and $\partial_{k_{y}} D_{k}^{α} = i I_{2}$ ⁠, where I₂ denotes the 2 × 2 identity matrix, so that

\partial_{k_{x}} H_{k}^{α} = (\begin{matrix} 0 & I_{2} \\ I_{2} & 0 \end{matrix}), \partial_{k_{y}} H_{k}^{α} = (\begin{matrix} 0 & - i I_{2} \\ i I_{2} & 0 \end{matrix}) .

(A2)

By degenerate perturbation theory,¹⁸ for small k, we have that eigenfunctions $χ_{k}^{α}$ of (A1) are given by

χ_{k}^{α} \approx \sum_{σ = \pm 1} c_{σ, k} ψ_{σ}^{α},

where the coefficients c_σ,k and associated eigenvalues E_k ≈ ϵ_k are found by solving the matrix eigenvalue problem

(\begin{matrix} \frac{⟨(ψ_{1}^{α}| k \cdot \nabla_{k} H_{0}^{α} ψ_{1}^{α}⟩}{⟨(ψ_{1}^{α}| ψ_{1}^{α}⟩} & \frac{⟨(ψ_{1}^{α}| k \cdot \nabla_{k} H_{0}^{α} ψ_{- 1}^{α}⟩}{⟨(ψ_{1}^{α}| ψ_{1}^{α}⟩} \\ \frac{⟨(ψ_{- 1}^{α}| k \cdot \nabla_{k} H_{0}^{α} ψ_{1}^{α}⟩}{⟨(ψ_{- 1}^{α}| ψ_{- 1}^{α}⟩} & \frac{⟨(ψ_{- 1}^{α}| k \cdot \nabla_{k} H_{0}^{α} ψ_{- 1}^{α}⟩}{⟨(ψ_{- 1}^{α}| ψ_{- 1}^{α}⟩} \end{matrix}) (\begin{matrix} c_{+ 1, k} \\ c_{- 1, k} \end{matrix}) = ϵ_{k} (\begin{matrix} c_{+ 1, k} \\ c_{- 1, k} \end{matrix}) .

(A3)

Using (A2) and the explicit forms of $ψ_{\pm 1}^{α}$ given by Proposition III.7, we find that the matrix on the left-hand side of (A3) can be simplified to

(\begin{matrix} 0 & λ (α) (k_{x} - i k_{y}) \\ λ^{*} (α) (k_{x} + i k_{y}) & 0 \end{matrix}), λ (α) ≔ \frac{⟨(ψ_{1}^{α} (r)| ψ_{1}^{α^{*}} (- r)⟩}{⟨(ψ_{1}^{α}| ψ_{1}^{α}⟩} .

It follows that, for small k, we have E_k ≈ ±v(α)|k|, where v(α) = |λ(α)| is as in (23).

APPENDIX B: THE CHIRAL BASIS OF $L_{K, 1}^{2}$ AND ACTION OF H⁰ AND H¹ WITH RESPECT TO THIS BASIS

1. The spectrum and eigenfunctions of H⁰ in $L_{K}^{2}$

The first task is to understand the spectrum and eigenfunctions of H⁰ in $L_{K}^{2}$ ⁠. In the Appendix B 2, we will discuss the spectrum and eigenfunctions of H⁰ in $L_{K, 1}^{2}$ ⁠. Recall that

H^{0} = (\begin{matrix} 0 & D^{0 †} \\ D^{0} & 0 \end{matrix}), D^{0} = (\begin{matrix} - 2 i \bar{\partial} & 0 \\ 0 & - 2 i \bar{\partial} \end{matrix}),

where $\bar{\partial} = \frac{1}{2} (\partial_{x} + i \partial_{y})$ ⁠. To describe the eigenfunctions of H⁰ in $L_{K}^{2}$ ⁠, we introduce some notation. Let $v = (\begin{matrix} v_{1}, v_{2} \end{matrix})$ be a vector in $R^{2}$ ⁠. Then, we will write

z_{v} = v_{1} + i v_{2}, {\hat{z}}_{v} = \frac{v_{1} + i v_{2}}{| v |} .

Finally, let V denote the area of the moiré cell Ω.

Proposition B.1.

The zero eigenspace of H⁰ in

L_{K}^{2}

is spanned by

χ_{\pm}^{0} = \frac{1}{\sqrt{2 V}} (\begin{matrix} 1, 0, \pm 1, 0 \end{matrix}) .

For all G ≠ 0 in the reciprocal lattice, then

χ_{\pm}^{G} (r) = \frac{1}{\sqrt{2 V}} (\begin{matrix} 1, 0, \pm {\hat{z}}_{G}, 0 \end{matrix}) e^{i G \cdot r}

are eigenfunctions with eigenvalues ±|G|. For all G in the reciprocal lattice,

χ_{\pm}^{q_{1} + G} (r) = \frac{1}{\sqrt{2 V}} (\begin{matrix} 0, 1, 0, \pm {\hat{z}}_{G + q_{1}} \end{matrix}) e^{i (q_{1} + G) \cdot r}

are eigenfunctions with eigenvalues ±|q₁ + G|. The operator H⁰ has no other eigenfunctions in

L_{K}^{2}

other than linear combinations of these, and hence, the spectrum of H⁰ in

L_{K}^{2}

is

σ_{L_{K}^{2}} (H^{0}) = \{\pm | G |, \pm | q_{1} + G | : G \in Λ^{*}\} .

Proof.

The proof is a straightforward calculation taking into account the $L_{K}^{2}$ boundary conditions given by (20) with k = 0. For example, e₂ and e₄ are zero eigenfunctions of H⁰ but in $L_{K^{'}}^{2}$ ⁠, not $L_{K}^{2}$ ⁠.□

Note that (as it must be because of the chiral symmetry) the spectrum is symmetric about 0 and all of the eigenfunctions with negative eigenvalues are given by applying $S$ to the eigenfunctions with positive eigenvalues.

The union of the lattices Λ^* and Λ^* + q₁ has the form of a honeycomb lattice in momentum space, where the lattice Λ^* corresponds to “A” sites and Λ^* + q₁ corresponds to “B” sites (or vice versa); see Fig. 4.

FIG. 4.

View large Download slide

Diagram showing A (blue) and B (red) sites of the momentum space lattice. Each site corresponds to two $L_{K}^{2}$ -eigenvalues of H⁰, given by ± the distance between the site and the origin (black). The lattice vectors b₁ and b₂ are shown as well as the A site nearest neighbor vectors q₁, q₂, and q₃.

2. The spectrum and eigenfunctions of H⁰ in $L_{K, 1}^{2}$

We now discuss the spectrum of H⁰ in $L_{K, 1}^{2}$ ⁠.

Proposition B.2.

The zero eigenspace of H⁰ in

L_{K, 1}^{2}

is spanned by

χ^{\tilde{0}} ≔ \frac{1}{\sqrt{V}} e_{1} .

For all G ≠ 0 in the reciprocal lattice Λ^*,

χ_{\pm}^{\tilde{G}} ≔ \frac{1}{\sqrt{3}} \sum_{k = 0}^{2} R^{k} χ_{\pm}^{G} = \frac{1}{\sqrt{3}} \sum_{k = 0}^{2} χ_{\pm}^{{(R_{ϕ}^{*})}^{k} G}

are eigenfunctions of H⁰ in

L_{K, 1}^{2}

with associated eigenvalues ±|G|. For all G in the reciprocal lattice Λ^*,

χ^{\pm \tilde{G + q_{1}}} = \frac{1}{\sqrt{3}} \sum_{k = 0}^{2} R^{k} χ_{\pm}^{G + q_{1}} = \frac{1}{\sqrt{3}} \sum_{k = 0}^{2} χ_{\pm}^{{(R_{ϕ}^{*})}^{k} (G + q_{1})}

are eigenfunctions of H⁰ in

L_{K, 1}^{2}

with associated eigenvalues ±|q₁ + G|. The operator H⁰ has no other eigenfunctions in

L_{K, 1}^{2}

other than linear combinations of these, and hence, the spectrum of H⁰ in

L_{K, 1}^{2}

is

σ_{L_{K, 1}^{2}} (H^{0}) = \{\pm | G |, \pm | q_{1} + G | : G \in Λ^{*}\} .

Proof.

The proof is another straightforward calculation starting from Proposition B.1.□

For an illustration of the support of the $L_{K, 1}^{2}$ -eigenfunctions of H⁰ on the momentum space lattice, see Fig. 5. It is important to note that the notation introduced in Proposition B.2 is not one-to-one because, for example,

χ^{\pm \tilde{G}} = χ^{\pm \tilde{R_{ϕ}^{*} G}} = χ^{\pm \tilde{{(R_{ϕ}^{*})}^{2} G}}

for any G ≠ 0 in Λ^*.

FIG. 5.

View large Download slide

Diagram showing support of $L_{K, 1}^{2}$ -eigenfunctions of H⁰ superposed on the momentum space lattice. Each eigenfunction is given by superposing an $L_{K}^{2}$ -eigenfunction of H⁰ with its rotations by $\frac{2 π}{3}$ and $\frac{4 π}{3}$ ⁠. The support of the eigenfunctions $χ^{\pm \tilde{q_{1}}}$ with eigenvalues ±1 is shown with the black crosses, while the support of the eigenfunctions $χ^{\pm \tilde{b_{1}}}$ with eigenvalues $\pm \sqrt{3}$ is shown with the black circles.

3. The chiral basis of $L_{K, 1}^{2}$

Recall that zero modes of H^α can be assumed to be eigenfunctions of the chiral symmetry operator $S$ ⁠. It follows that the most convenient basis for our purposes is not be the spectral basis just introduced but the basis of $L_{K, 1}^{2}$ consisting of eigenfunctions of $S$ ⁠. We call this basis the chiral basis.

Definition B.1.

The chiral basis of

L_{K, 1}^{2}

is defined as the union of the functions

χ^{\tilde{0}} = \frac{1}{\sqrt{V}} e_{1},

χ^{\tilde{G}, \pm 1} ≔ \frac{1}{\sqrt{2}} (χ^{\tilde{G}} \pm χ^{- \tilde{G}}), G \in Λ^{*} \ {0},

and

χ^{\tilde{q_{1} + G}, \pm 1} ≔ \frac{1}{\sqrt{2}} (χ^{\tilde{q_{1} + G}} \pm χ^{- \tilde{q_{1} + G}}), G \in Λ^{*} .

The following is straightforward.

Proposition B.3.

The chiral basis is an orthonormal basis of $L_{K, 1}^{2}$ . The modes $χ^{\tilde{0}}$ , $χ^{\tilde{G}, 1}$ , and $χ^{\tilde{q_{1} + G}, 1}$ are +1 eigenfunctions of $S$ , while the modes $χ^{\tilde{G}, - 1}$ and $χ^{\tilde{q_{1} + G}, - 1}$ are −1 eigenfunctions of $S$ .

Written out, chiral basis functions have a very simple form. We have

χ^{\tilde{0}} = \frac{1}{\sqrt{V}} e_{1},

(B1)

and for all G ∈ Λ^*\{0},

χ^{\tilde{G}, 1} (r) = \frac{1}{\sqrt{3 V}} e_{1} \sum_{k = 0}^{2} e^{i ({(R_{ϕ}^{*})}^{k} G) \cdot r}, χ^{\tilde{G}, - 1} (r) = \frac{1}{\sqrt{3 V}} {\hat{z}}_{G} e_{3} \sum_{k = 0}^{2} e^{- i k ϕ} e^{i ({(R_{ϕ}^{*})}^{k} G) \cdot r},

(B2)

and for all G ∈ Λ^*,

\begin{aligned} χ^{\tilde{G + q_{1}}, 1} (r) & = \frac{1}{\sqrt{3 V}} e_{2} \sum_{k = 0}^{2} e^{i ({(R_{ϕ}^{*})}^{k} (q_{1} + G)) \cdot r}, \\ χ^{\tilde{G + q_{1}}, - 1} (r) & = \frac{1}{\sqrt{3 V}} {\hat{z}}_{G + q_{1}} e_{4} \sum_{k = 0}^{2} e^{- i k ϕ} e^{i ({(R_{ϕ}^{*})}^{k} (q_{1} + G)) \cdot r} . \end{aligned}

(B3)

We use the chiral basis to divide up $L_{K, 1}^{2}$ as follows.

Definition B.2.

We define spaces $L_{K, 1, \pm 1}^{2}$ to be the spans of the ±1 eigenfunctions of $S$ in $L_{K, 1}^{2}$ , respectively.

Clearly, we have

L_{K, 1}^{2} = L_{K, 1, 1}^{2} \oplus L_{K, 1, - 1}^{2} .

We can divide up the chiral basis more finely as follows.

Definition B.3.

We define

\begin{aligned} L_{K, 1, 1, A}^{2} & ≔ \{χ^{\tilde{0}}\} \cup \{χ^{\tilde{G}, 1} : G \in Λ^{*} \ {0}\}, \\ L_{K, 1, 1, B}^{2} & ≔ \{χ^{\tilde{G + q_{1}}, 1} : G \in Λ^{*}\}, \\ L_{K, 1, - 1, A}^{2} & ≔ \{χ^{\tilde{G}, - 1} : G \in Λ^{*} \ {0}\}, \\ L_{K, 1, - 1, B}^{2} & ≔ \{χ^{\tilde{G + q_{1}}, - 1} : G \in Λ^{*}\} . \end{aligned}

Remark B.1.

Note that notations A and B in Definition B.3 refer to A and B sites of the momentum space lattice, not to the A and B sites of the real space lattice. Recalling Remark III.1 and comparing (B2) and (B3) with (2), we see that $L_{K, 1, 1, A}^{2}$ corresponds to wave-functions supported on A sites of layer 1, $L_{K, 1, 1, B}^{2}$ corresponds to wave-functions supported on A sites of layer 2, $L_{K, 1, - 1, A}^{2}$ corresponds to wave-functions supported on B sites of layer 1, and $L_{K, 1, - 1, B}^{2}$ corresponds to wave-functions supported on B sites of layer 2.

Clearly, we have

L_{K, 1}^{2} = L_{K, 1, 1, A}^{2} \oplus L_{K, 1, 1, B}^{2} \oplus L_{K, 1, - 1, A}^{2} \oplus L_{K, 1, - 1, B}^{2} .

The following propositions are straightforward to prove. For the first claim, note that ${S, H^{0}} = 0$ ⁠.

Proposition B.4.

Operator H⁰ maps

L_{K, 1, \pm 1, σ}^{2} \to L_{K, 1, \mp 1, σ}^{2}

for σ = A, B. The action of H⁰ on chiral basis functions is as follows:

H^{0} χ^{\tilde{0}} = 0

for all G ∈ Λ^* with G ≠ 0,

H^{0} χ^{\tilde{G}, \pm 1} = | G | χ^{\tilde{G}, \mp 1},

and for all G ∈ Λ^*,

H^{0} χ^{\tilde{q_{1} + G}, \pm 1} = | q_{1} + G | χ^{\tilde{q_{1} + G}, \mp 1} .

Proposition B.5.

Let P denote the projection operator onto

χ^{\tilde{0}}

in

L_{K, 1}^{2}

, and P^⊥ = 1 − P. Then, the operator

P^{⊥} {(H^{0})}^{- 1} P^{⊥}

maps

L_{K, 1, \pm 1, σ}^{2} \to L_{K, 1, \mp 1, σ}^{2}

for σ = A, B, and

P^{⊥} {(H^{0})}^{- 1} P^{⊥} χ^{\tilde{G}, \pm 1} = \frac{1}{| G |} χ^{\tilde{G}, \mp 1}

for all G ∈ Λ^* with G ≠ 0, and

P^{⊥} {(H^{0})}^{- 1} P^{⊥} χ^{\tilde{q_{1} + G}, \pm 1} = \frac{1}{| q_{1} + G |} χ^{\tilde{q_{1} + G}, \mp 1}

for all G ∈ Λ^*.

In the Appendixes B 4 and 5, we will study the action of the operator H¹ on $L_{K, 1}^{2}$ with respect to the chiral basis.

4. The spectrum of H¹ in $L_{K}^{2}$ and $L_{K, 1}^{2}$

Recall that

H^{1} = (\begin{matrix} 0 & D^{1 †} \\ D^{1} & 0 \end{matrix}), D^{1} = (\begin{matrix} 0 & U (r) \\ U (- r) & 0 \end{matrix}),

where $U (r) = e^{- i q_{1} \cdot r} + e^{i ϕ} e^{- i q_{2} \cdot r} + e^{- i ϕ} e^{- i q_{3} \cdot r}$ ⁠. We claim the following.

Proposition B.6.

For each r₀ ∈ Ω, ±|U(r₀)| and ±|U(−r₀)| are eigenvalues of

H^{1} : L_{K}^{2} \to L_{K}^{2}

. For r₀ such that U(r₀) ≠ 0, the ±|U(r₀)| eigenvectors are

(\begin{matrix} 0, 1, \pm \frac{U (r_{0})}{| U (r_{0}) |}, 0 \end{matrix}) δ (r - r_{0}) .

For r₀ such that U(−r₀) ≠ 0, the ±|U(−r₀)| eigenvectors are

(\begin{matrix} 1, 0, 0, \pm \frac{U (- r_{0})}{| U (- r_{0}) |} \end{matrix}) δ (r - r_{0}) .

When U(r₀) = 0, zero is a degenerate eigenvalue with associated eigenfunctions e₂δ(r − r₀) and e₃δ(r − r₀). When U(−r₀) = 0, zero is a degenerate eigenvalue with associated eigenfunctions e₁δ(r − r₀) and e₄δ(r − r₀). Finally,

σ_{L_{K}^{2}} (H^{1}) = [- 3, 3] .

(B4)

Proof.

We prove only () since the other assertions are clear. The triangle inequality yields the obvious bound

| U (r_{0}) | \leq 3

so that the

L_{K}^{2}

spectrum of H¹ must be contained in the interval [−3, 3]. To see that the spectrum actually equals [−3, 3], note that if

r_{0} = (\begin{matrix} \frac{4 π}{3 \sqrt{3}}, 0 \end{matrix})

⁠; then,

q_{1} \cdot r_{0} = 0, (q_{1} + b_{1}) \cdot r_{0} = \frac{1}{2} (\begin{matrix} \sqrt{3}, 1 \end{matrix}) \cdot r_{0} = \frac{2 π}{3}, (q_{1} + b_{2}) \cdot r_{0} = \frac{1}{2} (\begin{matrix} - \sqrt{3}, 1 \end{matrix}) \cdot r_{0} = - \frac{2 π}{3},

and hence, U(r₀) = 3. On the other hand, when r₀ = 0, we have U(r₀) = 0 so that the spectrum of H¹ in

L_{K}^{2}

equals [−3, 3].□

By taking linear combinations of rotated copies of the H¹ eigenfunctions, just as we did with the H⁰ eigenfunctions, it is straightforward to prove an analogous result to Proposition B.6 in $L_{K, 1}^{2}$ ⁠. We record only the following.

Proposition B.7.

σ_{L_{K, 1}^{2}} (H^{1}) = [- 3, 3] .

5. The action of H¹ on $L_{K, 1}^{2}$ with respect to the chiral basis

We now want to study the action of H¹ on $L_{K, 1}^{2}$ with respect to the chiral basis. We will prove two propositions, which parallel Proposition B.4.

Proposition B.8.

The operator H¹ maps

L_{K, 1, 1, A}^{2} \to L_{K, 1, - 1, B}^{2}

and

L_{K, 1, 1, B}^{2} \to L_{K, 1, - 1, A}^{2}

. The action of H¹ on chiral basis functions is as follows:

H^{1} χ^{\tilde{0}} = \sqrt{3} \bar{{\hat{z}}_{q_{1}}} χ^{\tilde{q_{1}}, - 1}

(B5)

and

H^{1} χ^{\tilde{q_{1}}, 1} = e^{i ϕ} \bar{{\hat{z}}_{q_{1} - q_{2}}} χ^{\tilde{q_{1} - q_{2}}, - 1} + e^{- i ϕ} \bar{{\hat{z}}_{q_{1} - q_{3}}} χ^{\tilde{q_{1} - q_{3}}, - 1} .

(B6)

For all G ∈ Λ^*\{0},

H^{1} χ^{\tilde{G}, 1} = \bar{{\hat{z}}_{G + q_{1}}} χ^{\tilde{G + q_{1}}, - 1} + e^{i ϕ} \bar{{\hat{z}}_{G + q_{2}}} χ^{\tilde{G + q_{2}}, - 1} + e^{- i ϕ} \bar{{\hat{z}}_{G + q_{3}}} χ^{\tilde{G + q_{3}}, - 1} .

(B7)

For all G ∈ Λ^*\{0},

H^{1} χ^{\tilde{G + q_{1}}, 1} = \bar{{\hat{z}}_{G}} χ^{\tilde{G}, - 1} + e^{i ϕ} \bar{{\hat{z}}_{G + q_{1} - q_{2}}} χ^{\tilde{G + q_{1} - q_{2}}, - 1} + e^{- i ϕ} \bar{{\hat{z}}_{G + q_{1} - q_{3}}} χ^{\tilde{G + q_{1} - q_{3}}, - 1} .

(B8)

Note that H¹ exchanges chirality (⁠ $S$ eigenvalue) and the A and B momentum space sublattices, while H⁰ only exchanges chirality. Proposition B.8 has a simple interpretation in terms of nearest neighbor hopping in the momentum space lattice; see Figs. 6 and 7.

FIG. 6.

View large Download slide

Illustration of the action of H¹ in $L_{K, 1}^{2}$ as hopping in the momentum space lattice described by Eq. (B7) (left, starting at b₁) and (B8) (right, starting at q₁ + b₁ − b₂). The origin is marked by a black dot.

FIG. 7.

View large Download slide

Illustration of the action of H¹ as hopping in the momentum space lattice described by Eq. (B5) (left, starting at 0) and (B6) (right, starting at q₁). Although it appears that the hopping in these cases does not respect $\frac{2 π}{3}$ rotation symmetry, this is an artifact of working with chiral basis functions, which individually respect the rotation symmetry; see (B9).

Remark B.2.

At first glance, (B5) and (B6) appear different from (B7) and (B8) because they appear to violate

\frac{2 π}{3}

rotation symmetry. However, this is not the case since every chiral basis function individually respects this symmetry. For example, using

χ^{\tilde{q_{1}}, - 1} = χ^{\tilde{q_{2}}, - 1} = χ^{\tilde{q_{3}}, - 1}

and

\bar{{\hat{z}}_{q_{1}}} = e^{i ϕ} \bar{{\hat{z}}_{q_{2}}} = e^{- i ϕ} \bar{{\hat{z}}_{q_{3}}}

, we can re-write (B5) in a way that manifestly respects the

\frac{2 π}{3}

rotation symmetry as

H^{1} χ^{\tilde{0}} = \frac{1}{\sqrt{3}} (\bar{{\hat{z}}_{q_{1}}} χ^{\tilde{q_{1}}, - 1} + e^{i ϕ} \bar{{\hat{z}}_{q_{2}}} χ^{\tilde{q_{2}}, - 1} + e^{- i ϕ} \bar{{\hat{z}}_{q_{3}}} χ^{\tilde{q_{3}}, - 1}) .

(B9)

Equation (B6) can also be written in a manifestly rotationally invariant way, but the expression is long, and hence, we omit it. Note that (B6) cannot have a term proportional to

χ^{\tilde{0}}

since

χ^{\tilde{0}} \in L_{K, 1, 1}^{2}

and H¹ maps

L_{K, 1, 1}^{2} \to L_{K, 1, - 1}^{2}

.

Proof of Proposition B.8.

We will prove (); the proofs of the other identities are similar and hence omitted. We have

H^{1} χ^{\tilde{G}, 1} = \frac{1}{\sqrt{3 V}} (e^{i q_{1} \cdot r} + e^{i ϕ} e^{i (q_{1} + b_{1}) \cdot r} + e^{- i ϕ} e^{i (q_{1} + b_{2}) \cdot r}) (e^{i G \cdot r} + e^{i (R_{ϕ}^{*} G) \cdot r} + e^{i ({(R_{ϕ}^{*})}^{2} G) \cdot r}) e_{4} .

Multiplying out, we have

\begin{aligned} \frac{1}{\sqrt{3 V}} (e^{i (q_{1} + G) \cdot r} + e^{i ϕ} e^{i (q_{1} + G + b_{1}) \cdot r} + e^{- i ϕ} e^{i (q_{1} + G + b_{2}) \cdot r}) \\ + e^{i (q_{1} + (R_{ϕ}^{*} G)) \cdot r} + e^{i ϕ} e^{i (q_{1} + (R_{ϕ}^{*} G) + b_{1}) \cdot r} + e^{- i ϕ} e^{i (q_{1} + (R_{ϕ}^{*} G) + b_{2}) \cdot r} \\ (+ e^{i (q_{1} + ({(R_{ϕ}^{*})}^{2} G)) \cdot r} + e^{i ϕ} e^{i (q_{1} + ({(R_{ϕ}^{*})}^{2} G) + b_{1}) \cdot r} + e^{- i ϕ} e^{i (q_{1} + ({(R_{ϕ}^{*})}^{2} G) + b_{2}) \cdot r}) \\ = & \frac{1}{\sqrt{3 V}} (e^{i (q_{1} + G) \cdot r} + e^{i ϕ} e^{i (q_{1} + G + b_{1}) \cdot r} + e^{- i ϕ} e^{i (q_{1} + G + b_{2}) \cdot r}) \\ e^{i (R_{ϕ}^{*} (q_{1} + G + b_{1})) \cdot r} + e^{i ϕ} e^{i (R_{ϕ}^{*} (q_{1} + G + b_{2})) \cdot r} + e^{- i ϕ} e^{i (R_{ϕ}^{*} (q_{1} + G)) \cdot r} \\ (+ e^{i ({(R_{ϕ}^{*})}^{2} (q_{1} + G + b_{2})) \cdot r} + e^{i ϕ} e^{i ({(R_{ϕ}^{*})}^{2} (q_{1} + G)) \cdot r} + e^{- i ϕ} e^{i () {(R_{ϕ}^{*})}^{2} (q_{1} + G + b_{1}) \cdot r}) \\ = & \frac{1}{\sqrt{3 V}} (e^{i (q_{1} + G) \cdot r} + e^{- i ϕ} e^{i (R_{ϕ}^{*} (q_{1} + G)) \cdot r} + e^{i ϕ} e^{i ({(R_{ϕ}^{*})}^{2} (q_{1} + G)) \cdot r}) \\ + \frac{1}{\sqrt{3 V}} e^{i ϕ} (e^{i (q_{1} + G + b_{1}) \cdot r} + e^{- i ϕ} e^{i (R_{ϕ}^{*} (q_{1} + G + b_{1})) \cdot r} + e^{i ϕ} e^{i ({(R_{ϕ}^{*})}^{2} (q_{1} + G + b_{1})) \cdot r}) \\ + \frac{1}{\sqrt{3 V}} e^{- i ϕ} (e^{i (q_{1} + G + b_{2}) \cdot r} + e^{- i ϕ} e^{i (R_{ϕ}^{*} (q_{1} + G + b_{2})) \cdot r} + e^{i ϕ} e^{i ({(R_{ϕ}^{*})}^{2} (q_{1} + G + b_{2})) \cdot r}), \end{aligned}

from which () follows.□

Proposition B.9.

The operator H¹ maps

L_{K, 1, - 1, A}^{2} \to L_{K, 1, 1, B}^{2}

and

L_{K, 1, - 1, B}^{2} \to L_{K, 1, 1, A}^{2}

. The action of H¹ on chiral basis functions is as follows:

H^{1} χ^{\tilde{q_{1}}, - 1} = {\hat{z}}_{q_{1}} (\sqrt{3} χ^{\tilde{0}} + e^{- i ϕ} χ^{\tilde{q_{1} - q_{2}}, 1} + e^{i ϕ} χ^{\tilde{q_{1} - q_{3}}, 1}) .

For all G ∈ Λ^*\{0},

H^{1} χ^{\tilde{G}, - 1} = {\hat{z}}_{G} (χ^{\tilde{G + q_{1}}, 1} + e^{- i ϕ} χ^{\tilde{G + q_{2}}, 1} + e^{i ϕ} χ^{\tilde{G + q_{3}}, 1}) .

For all G ∈ Λ^*,

H^{1} χ^{\tilde{G + q_{1}}, - 1} = {\hat{z}}_{G + q_{1}} (χ^{\tilde{G}, 1} + e^{- i ϕ} χ^{\tilde{G + q_{1} - q_{2}}, 1} + e^{i ϕ} χ^{\tilde{G + q_{1} - q_{3}}, 1}) .

Proof.

The proof is similar to that of Proposition B.8 and is hence omitted.□

APPENDIX C: FORMAL EXPANSION OF THE ZERO MODE

We now bring to bear the developments of the Appendix B on the asymptotic expansion of the zero mode $ψ^{α} \in L_{K, 1, 1}^{2}$ starting from $Ψ^{0} = e_{1} = χ^{\tilde{0}}$ ⁠. We first give the Proof of Proposition IV.1.

Proof of Proposition IV.1.

We have seen that

χ^{\tilde{0}} \in L_{K, 1, 1}^{2}

⁠. By the calculations of the Appendix B 5,

H^{1} χ^{\tilde{0}} \in L_{K, 1, - 1}^{2}

⁠, which is orthogonal to the null space of H⁰. The general solution of H⁰Ψ¹ = −H¹Ψ⁰ is

Ψ^{1} = - P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1} Ψ^{0} + C Ψ^{0},

where C is an arbitrary constant, which is in

L_{K, 1, 1}^{2}

by Proposition B.4. To ensure that Ψ¹ is orthogonal to Ψ⁰, we take C = 0. It is clear that this procedure can be repeated to derive an expansion to all orders satisfying the conditions of Proposition IV.1.□

Our goal is to calculate $Ψ^{n} \in L_{K, 1, 1}^{2}$ satisfying the conditions of Proposition B.4 up to n = 8. This amounts to calculating, for n = 1 to n = 8,

Ψ^{n} = - P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1} Ψ^{n - 1} .

We do this algorithmically by repeated application of the following proposition, which combines Proposition B.8 and Proposition B.5.

Proposition C.1.

The operator

- P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1}

maps

L_{K, 1, 1, A}^{2} \to L_{K, 1, 1, B}^{2}

and

L_{K, 1, 1, B}^{2} \to L_{K, 1, 1, A}^{2}

. Its action on chiral basis functions is as follows:

- P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1} χ^{\tilde{0}} = - \sqrt{3} \bar{{\hat{z}}_{q_{1}}} χ^{\tilde{q_{1}}, 1}

(C1)

and

- P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1} χ^{\tilde{q_{1}}, 1} = - \frac{e^{i ϕ} \bar{{\hat{z}}_{q_{1} - q_{2}}}}{| q_{1} - q_{2} |} χ^{\tilde{q_{1} - q_{2}}, 1} - \frac{e^{- i ϕ} \bar{{\hat{z}}_{q_{1} - q_{3}}}}{| q_{1} - q_{3} |} χ^{\tilde{q_{1} - q_{3}}, 1} .

(C2)

For all G ∈ Λ^*\{0},

\begin{aligned} - P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1} χ^{\tilde{G}, 1} = \\ - \frac{\bar{{\hat{z}}_{G + q_{1}}}}{| G + q_{1} |} χ^{\tilde{G + q_{1}}, 1} - \frac{e^{i ϕ} \bar{{\hat{z}}_{G + q_{2}}}}{| G + q_{2} |} χ^{\tilde{G + q_{2}}, 1} - \frac{e^{- i ϕ} \bar{{\hat{z}}_{G + q_{3}}}}{| G + q_{3} |} χ^{\tilde{G + q_{3}}, 1} . \end{aligned}

(C3)

For all G ∈ Λ^*\{0},

- P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1} χ^{\tilde{G + q_{1}}, 1} = - \frac{\bar{{\hat{z}}_{G}}}{| G |} χ^{\tilde{G}, 1} - \frac{e^{i ϕ} \bar{{\hat{z}}_{G + q_{1} - q_{2}}}}{| G + q_{1} - q_{2} |} χ^{\tilde{G + q_{1} - q_{2}}, 1} - \frac{e^{- i ϕ} \bar{{\hat{z}}_{G + q_{1} - q_{3}}}}{| G + q_{1} - q_{3} |} χ^{\tilde{G + q_{1} - q_{3}}, 1} .

We now claim the following.

Proposition C.2.

Let Ψⁿ be the sequence defined by Proposition IV.1. Then,

Ψ^{1} = - \sqrt{3} i χ^{\tilde{q_{1}}, 1},

(C4)

Ψ^{2} = (\frac{\sqrt{3} - i}{2}) χ^{\tilde{- b_{1}}, 1} + (\frac{\sqrt{3} + i}{2}) χ^{\tilde{- b_{2}}, 1},

(C5)

Ψ^{3} = \frac{1}{\sqrt{7}} (\frac{\sqrt{7} - 3 \sqrt{21} i}{14}) χ^{\tilde{q_{1} - b_{2}}, 1} + \frac{1}{\sqrt{7}} (\frac{- \sqrt{7} - 3 \sqrt{21} i}{14}) χ^{\tilde{q_{1} - b_{1}}, 1},

(C6)

\begin{aligned} Ψ^{4} = & \frac{1}{\sqrt{21}} (\frac{- 5 \sqrt{7} + \sqrt{21} i}{14}) χ^{\tilde{- b_{2}}, 1} + \frac{1}{2 \sqrt{21}} (\frac{2 \sqrt{7} + \sqrt{21} i}{7}) χ^{\tilde{- 2 b_{2}}, 1} \\ + \frac{1}{\sqrt{21}} (\frac{- 5 \sqrt{7} - \sqrt{21} i}{14}) χ^{\tilde{- b_{1}}, 1} + \frac{1}{2 \sqrt{21}} (\frac{2 \sqrt{7} - \sqrt{21} i}{7}) χ^{\tilde{- 2 b_{1}}} \\ + \frac{2 \sqrt{3}}{21} χ^{\tilde{- b_{1} - b_{2}}, 1} . \end{aligned}

(C7)

Proof.

Equations () and () follow immediately from () and () and using q₂ = q₁ + b₁ and q₃ = q₁ + b₂. The derivation of Eq. () is more involved, so we give details. Using linearity and applying () twice, we find that

\begin{aligned} - P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1} Ψ^{2} = \\ (\frac{\sqrt{3} - i}{2}) (\frac{\bar{{\hat{z}}_{q_{1} - b_{2}}}}{| q_{1} - b_{2} |} χ^{\tilde{q_{1} - b_{2}}, 1} + \frac{e^{i ϕ} \bar{{\hat{z}}_{q_{1} + b_{1} - b_{2}}}}{| q_{1} + b_{1} - b_{2} |} χ^{\tilde{q_{1} + b_{1} - b_{2}}, 1} + e^{- i ϕ} \bar{{\hat{z}}_{q_{1}}} χ^{\tilde{q_{1}}, 1}) \\ + (\frac{\sqrt{3} + i}{2}) (\frac{\bar{{\hat{z}}_{q_{1} - b_{1}}}}{| q_{1} - b_{1} |} χ^{\tilde{q_{1} - b_{1}}, 1} + e^{i ϕ} \bar{{\hat{z}}_{q_{1}}} χ^{\tilde{q_{1}}, 1} + \frac{e^{- i ϕ} \bar{{\hat{z}}_{q_{1} + b_{2} - b_{1}}}}{| q_{1} + b_{2} - b_{1} |} χ^{\tilde{q_{1} + b_{2} - b_{1}}, 1}) . \end{aligned}

First, the terms proportional to

χ^{\tilde{q_{1}}, 1}

cancel. Next, since R_ϕ(q₁ + b₁ − b₂) = q₁ + b₂ − b₁, we have

χ^{\tilde{q_{1} + b_{1} - b_{2}}, 1} = χ^{\tilde{q_{1} + b_{2} - b_{1}}, 1}

⁠. These terms also cancel, leaving (). The derivation of () (and the higher corrections) is involved but does not depend on any new ideas and is therefore omitted.□

We give the explicit forms of Ψ⁵-Ψ⁸ in the supplementary material.

Remark C.1.

Written out, (C4) and (C5) become

Ψ^{1} = - \sqrt{3} i \frac{1}{\sqrt{3 V}} e_{2} (e^{i q_{1} \cdot r} + e^{i q_{2} \cdot r} + e^{i q_{3} \cdot r})

and

Ψ^{2} = - i e^{i ϕ} \frac{1}{\sqrt{3 V}} e_{1} (e^{i b_{1} \cdot r} + e^{i (b_{2} - b_{1}) \cdot r} + e^{- i b_{2} \cdot r}) + i e^{- i ϕ} \frac{1}{\sqrt{3 V}} e_{1} (e^{i b_{2} \cdot r} + e^{- i b_{1} \cdot r} + e^{i (b_{1} - b_{2}) \cdot r}),

which agree with Eq. (24) of Tarnopolsky et al.⁴ up to an overall factor of

\sqrt{V}

(this factor cancels in the Fermi velocity, so there is no discrepancy).

Using orthonormality of the chiral basis functions, it is straightforward to calculate the norms of each of the Ψⁿ. We have the following.

Proposition C.3.

‖ Ψ^{0} ‖ = 1, ‖ Ψ^{1} ‖ = \sqrt{3}, ‖ Ψ^{2} ‖ = \sqrt{2}, ‖ Ψ^{3} ‖ = \frac{\sqrt{14}}{7}, ‖ Ψ^{4} ‖ = \frac{\sqrt{258}}{42}, ‖ Ψ^{5} ‖ = \frac{\sqrt{1 968 837}}{3458},

‖ Ψ^{6} ‖ = \frac{\sqrt{106 525 799}}{31 122}, ‖ Ψ^{7} ‖ = \frac{2 \sqrt{2 129 312 323 981 473}}{624 696 345}, ‖ Ψ^{8} ‖ = \frac{\sqrt{183 643 119 755 214 454}}{4 997 570 760} .

Remark C.2.

Note that the sequence of norms of the expansion functions grows much slower than the pessimistic bound ‖Ψ^N+1‖ ≤ 3‖Ψ^N‖, N = 0, 1, 2, … guaranteed by Proposition IV.2. The reason is that bounds (27) and (28) are never attained. As N becomes larger, the bound (27) is very pessimistic because Ψ^N is mostly made up of eigenfunctions of H⁰ with eigenvalues strictly larger than 1. The bound (28) is also very pessimistic because it is attained only at delta functions, which can only be approximated with a superposition of a large number of eigenfunctions of H⁰. It seems possible that a sharper bound could be proved starting from these observations, but we do not pursue this in this work.

We finally give the Proof of Proposition IV.8.

Proof of Proposition IV.8.

Explicit computation using Proposition B.8 and orthonormality of the chiral basis functions gives

‖ H^{1} Ψ^{8} ‖ = \frac{\sqrt{4 855 076 200 233 765 642}}{14 992 712 280} \approx 0.147 \leq \frac{3}{20} .

□

APPENDIX D: PROOF OF PROPOSITION IV.4

We choose Ξ as

Ξ ≔ \{\begin{matrix} L_{K, 1}^{2} - eigenfunctions of H^{0} with \\ eigenvalues with magnitude \leq 4 \sqrt{3} \end{matrix}\} ⋃ \{χ^{\tilde{q_{1} - 4 b_{1} + b_{2}}, \pm 1}, χ^{\tilde{q_{1} + b_{1} - 4 b_{2}}, \pm 1}\} .

Part 1 of Proposition IV.4 follows immediately from observing that $χ^{\tilde{q_{1} - 2 b_{1} - 2 b_{2}}, \pm 1}$ is not in Ξ, but |q₁ − 2b₁ − 2b₂| = 7. That μ = 7 is optimal can be seen from Fig. 8.

FIG. 8.

View large Download slide

Illustration of Ξ in the momentum space lattice. The circle has radius $4 \sqrt{3}$ so that every dot within the circle corresponds to two chiral basis vectors included in Ξ. Chiral basis vectors exactly $4 \sqrt{3}$ away from the origin, marked with black dots, are also included in Ξ. We also include in Ξ the chiral basis vectors $\{χ^{\tilde{q_{1} - 4 b_{1} + b_{2}}, \pm 1}, χ^{\tilde{q_{1} + b_{1} - 4 b_{2}}, \pm 1}\}$ ⁠, which correspond to the dots marked with circles, which are distance 7 (NB. $7 > 4 \sqrt{3}$ ⁠) from the origin. We do not include the chiral basis vectors $χ^{\tilde{q_{1} - 2 b_{1} - 2 b_{2}}, \pm 1}$ ⁠, marked with the black crosses, which are also a distance 7 from the origin. The reason for this is so that part 3 of Proposition IV.4 holds. With this choice, every dot in Ξ has at most one nearest neighbor lattice point outside of Ξ. It follows immediately from Propositions B.8 and B.9 (H¹ acts by nearest neighbor hopping in the momentum space lattice) that $‖ P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ = 1$ ⁠. Note that if we chose Ξ to include $χ^{\tilde{q_{1} - 2 b_{1} - 2 b_{2}}, \pm 1}$ ⁠, this would no longer hold because these basis functions would have two nearest neighbors outside Ξ, resulting in the worse bound $‖ P_{Ξ} H^{1} P_{Ξ}^{⊥} ‖ \leq \sqrt{2}$ ⁠.

Part 2 follows from the fact that ψ^8,α depends only on eigenfunctions of H⁰ with eigenvalues with magnitude less than or equal to $4 \sqrt{3}$ ⁠. The largest eigenvalue is $4 \sqrt{3}$ ⁠, coming from dependence of Ψ⁸ on $χ^{\tilde{- 4 b_{2}}, 1}$ ⁠, since $| - 4 b_{2} | = 4 \sqrt{3}$ ⁠.

Part 3 can be seen from Fig. 8.

APPENDIX E: PROOF OF PROPOSITION IV.5

1. Proof of Theorem IV.2

We will prove Theorem IV.2 starting from Theorem 11.5.1 of Parlett,¹⁹ where the proof can be found.

Lemma E.1.

Let Q be an n × m matrix, which satisfies Q^†Q = I_m. Define H = Q^†AQ and R = AQ − QH. Let

{θ_{j}}_{1 \leq j \leq m}

denote the eigenvalues of H (the Ritz values). Then, m of A’s eigenvalues

{α_{j}}_{1 \leq j \leq m}

can be put into one-to-one correspondence with the

{θ_{j}}_{1 \leq j \leq m}

in such a way that

| θ_{j} - α_{j} | \leq ‖ R ‖_{2} 1 \leq j \leq m .

Proof of Theorem IV.2.

Let Q be the matrix whose columns are v₁, …, v_m. Using orthonormality of the v_j, Q^†Q = I_m and

H = Q^{†} A Q = d i a g (λ_{1}, \dots, λ_{m}) + (\begin{matrix} ⟨(v_{1}| r_{1}⟩ & \dots & ⟨(v_{1}| r_{m}⟩ \\ \dots & \dots & \dots \\ ⟨(v_{n}| r_{1}⟩ & \dots & ⟨(v_{n}| r_{n}⟩ \end{matrix}) .

We now prove that the eigenvalues of H, denoted by θ_j, are close to λ_j’s. By the Gershgorin circle theorem, we have

| θ_{i} - (λ_{i} + ⟨(v_{i}| r_{i}⟩) | \leq \sum_{j \neq i}^{m} | ⟨(v_{i}| r_{j}⟩ |,

which implies, using ‖v_j‖₂ = 1,

| θ_{i} - λ_{i} | = | θ_{i} - λ_{i} - ⟨(v_{i}| r_{i}⟩ + ⟨(v_{i}| r_{i}⟩ | \leq \sum_{j = 1}^{m} | ⟨(v_{i}| r_{j}⟩ | \leq m \sup_{1 \leq i \leq m} ‖ r_{i} ‖_{2} .

We can now use Lemma E.1 to bound the difference between λ_js and exact eigenvalues α_j,

| λ_{j} - α_{j} | = | λ_{j} - θ_{j} + θ_{j} - α_{j} | \leq m \sup_{1 \leq i \leq m} ‖ r_{i} ‖_{2} + ‖ R ‖_{2},

where R ≔ AQ − QH = (I − QQ^†)AQ. Since QQ^† projects onto v_j, R simplifies to

R = (I - Q Q^{†}) R^{'}, R^{'} ≔ (\begin{matrix} r_{1} & \dots & r_{m} \end{matrix}) .

Since QQ^† is a projection, so is I − QQ^†, and hence, ‖R‖₂ ≤ ‖R′‖₂. To bound ‖R′‖₂, note that for any v with ‖v‖₂ = 1, we have

‖ R^{'} v ‖_{2} = ⟨(e_{1}| v⟩ r_{1} + \dots + ⟨(e_{m}| v⟩ r_{m} \leq m \sup_{1 \leq i \leq m} ‖ r_{i} ‖_{2},

where e_j denote the standard orthonormal basis vectors. The result now follows.□

2. Proof of Theorem IV.3

Proof of Theorem IV.3.

We start with the following lemma, which guarantees that numerically computed approximately orthonormal sets can be approximated by exactly orthonormal sets.

Lemma E.2.

Let

{\tilde{v}}_{1}, \dots, {\tilde{v}}_{m}

be n-dimensional vectors, let

{⟨({\tilde{v}}_{i}| {\tilde{v}}_{j}⟩}_{c o m p}

for 1 ≤ i, j ≤ m denote their numerically computed inner products, let ϵ denote machine epsilon, and assume that nϵ < 0.01. Define

μ ≔ (1.01) n^{2} ϵ \sup_{1 \leq i \leq m} ‖ {\tilde{v}}_{i} ‖_{\infty}^{2} + \sup_{i} | {⟨({\tilde{v}}_{i}| {\tilde{v}}_{i}⟩}_{c o m p} - 1 | + \sup_{i \neq j} | {⟨({\tilde{v}}_{i}| {\tilde{v}}_{j}⟩}_{c o m p} | .

(E1)

Then, as long as

m μ < \frac{1}{2}

⁠, there is a set of n-dimensional orthonormal vectors

{\hat{v}}_{1}, \dots, {\hat{v}}_{m}

⁠, which satisfy

‖ {\hat{v}}_{j} - {\tilde{v}}_{j} ‖_{2} \leq 2^{- 1 / 2} m μ, 1 \leq j \leq m .

Proof.

Bounding the round-off error in computing inner products in the usual way (see, for example, Chap. 2.7 of Golub and Van Loan²⁰) and assuming that nϵ < 0.01, we have that for each 1 ≤ i, j ≤ m,

⟨({\tilde{v}}_{i}| {\tilde{v}}_{j}⟩ = {⟨({\tilde{v}}_{i}| {\tilde{v}}_{j}⟩}_{c o m p} + e_{i j}

⁠, where

| e_{i j} | \leq (1.01) n ϵ | {\tilde{v}}_{i} |^{⊤} | {\tilde{v}}_{j} | \leq (1.01) n ϵ ‖ {\tilde{v}}_{i} ‖_{2} ‖ {\tilde{v}}_{j} ‖_{2}

⁠. Letting

\tilde{Q}

denote the matrix whose columns are

{\tilde{v}}_{i}

’s, then

{\tilde{Q}}^{†} \tilde{Q} - I_{m} = E

⁠, where, for all i ≠ j,

| E_{i j} | \leq | {⟨({\tilde{v}}_{i}| {\tilde{v}}_{j}⟩}_{c o m p} | + (1.01) n ϵ ‖ {\tilde{v}}_{i} ‖_{2} ‖ {\tilde{v}}_{j} ‖_{2}

⁠, and for all i,

| E_{i i} | \leq | 1 - {⟨({\tilde{v}}_{i}| {\tilde{v}}_{i}⟩}_{c o m p} | + (1.01) n ϵ ‖ {\tilde{v}}_{i} ‖_{2}^{2}

⁠. Paying the price of factors of

\sqrt{n}

to replace ‖·‖₂ norms by ‖·‖_∞ norms, we can obtain a trivial bound on the maximal element of E, denoted as ‖E‖_max, by

‖ E ‖_{\max} \leq (1.01) n^{2} ϵ {(\sup_{i} ‖ {\tilde{v}}_{i} ‖_{\infty})}^{2} + \sup_{i} | {⟨({\tilde{v}}_{i}| {\tilde{v}}_{i}⟩}_{c o m p} - 1 | + \sup_{i, j} | {⟨({\tilde{v}}_{i}| {\tilde{v}}_{j}⟩}_{c o m p} | .

Note that this is nothing but μ in the statement of the theorem. Using the Gershgorin circle theorem, we then have that the eigenvalues λ of

{\tilde{Q}}^{†} \tilde{Q}

satisfy |λ − 1| ≤ m‖E‖_max. We claim that there are exact orthonormal vectors

{\hat{v}}_{j}

near (in the ‖·‖₂-norm) to the

{\tilde{v}}_{j}

⁠. To see this, note that

\hat{Q} ≔ \tilde{Q} {({\tilde{Q}}^{†} \tilde{Q})}^{- 1 / 2}

satisfies

{\hat{Q}}^{†} \hat{Q} = I_{m}

and

‖ \hat{Q} - \tilde{Q} ‖_{2} \leq ‖ \hat{Q} ‖_{2} ‖ {({\tilde{Q}}^{†} \tilde{Q})}^{1 / 2} - I_{m} ‖_{2} = ‖ {({\tilde{Q}}^{†} \tilde{Q})}^{1 / 2} - I_{m} ‖_{2} .

Let λ_max and λ_min denote the maximum and minimum eigenvalues of

{\tilde{Q}}^{†} \tilde{Q}

⁠, respectively. Then,

{‖{({\tilde{Q}}^{†} \tilde{Q})}^{1 / 2} - I_{m}‖}_{2} \leq \max \{| λ_{\min}^{1 / 2} - 1 |, | λ_{\max}^{1 / 2} - 1 |\}

⁠. Since λ_min is bounded below by 1 − m‖E‖_max and λ_max is bounded above by 1 + m‖E‖_max, we have

‖ {({\tilde{Q}}^{†} \tilde{Q})}^{1 / 2} - I_{m} ‖_{2} \leq \max \{| {(1 + m ‖ E ‖_{\max})}^{1 / 2} - 1 |, | {(1 - m ‖ E ‖_{\max})}^{1 / 2} - 1 |\} .

Using Taylor’s theorem, for

| x | < \frac{1}{2}

⁠, we have that |(1 + x)^1/2 − 1| ≤ 2^−1/2|x| and |(1 − x)^1/2 − 1| ≤ 2^−1/2|x|. Since by assumption

m ‖ E ‖_{\max} < \frac{1}{2}

⁠, we conclude

‖ \tilde{Q} - \hat{Q} ‖_{2} \leq ‖ {({\tilde{Q}}^{†} \tilde{Q})}^{1 / 2} - I_{m} ‖_{2} \leq 2^{- 1 / 2} m ‖ E ‖_{\max} .

Letting

{\hat{v}}_{j}

denote the columns of

\hat{Q}

and noting that

‖ {\hat{v}}_{j} - {\tilde{v}}_{j} ‖_{2} \leq ‖ \hat{Q} - \tilde{Q} ‖_{2}

for all 1 ≤ j ≤ m, the result is proved.□

Using Lemma E.2, we have that there exists an exactly orthonormal set ${{\hat{v}}_{j}}_{1 \leq j \leq m}$ nearby to the set ${{\tilde{v}}_{j}}_{1 \leq j \leq m}$ ⁠. We now want to bound the residuals of the ${\hat{v}}_{j}$ in terms of numerically computable quantities. We start with the following easy lemma whose proof is a straightforward manipulation.

Lemma E.3.

Let A be an n × n Hermitian matrix, and suppose that

\hat{r} ≔ (A - \tilde{λ} I) \hat{v}

and

\tilde{r} ≔ (A - \tilde{λ} I) \tilde{v}

. Then,

‖ \hat{r} ‖_{2} \leq (‖ A ‖_{2} + | \tilde{λ} |) ‖ \hat{v} - \tilde{v} ‖_{2} + ‖ \tilde{r} ‖_{2} .

The following lemma quantifies the error in approximating exact residuals by numerically computed values.

Lemma E.4.

Let A be a Hermitian n × n matrix, and let

\tilde{A}

denote the matrix whose entries are those of A evaluated as floating-point numbers. Let

{[(\tilde{A} - \tilde{λ} I) \tilde{v}]}_{c o m p}

denote the numerically computed value of

(\tilde{A} - \tilde{λ} I) \tilde{v}

in the floating-point arithmetic. Then,

\tilde{r} ≔ (A - \tilde{λ} I) \tilde{v}

satisfies

‖ \tilde{r} ‖_{2} \leq n^{1 / 2} {‖{[(\tilde{A} - \tilde{λ} I) \tilde{v}]}_{c o m p}‖}_{\max} + (1.01) n^{5 / 2} ϵ ‖ \tilde{A} - \tilde{λ} I ‖_{\max} ‖ \tilde{v} ‖_{\infty} + n ϵ ‖ A ‖_{\max} ‖ \tilde{v} ‖_{\infty} .

Proof.

For matrices A and B with entries A_ij and B_ij, we will write |A| to denote the matrix with entries |A_ij| for all i, j, and |A| ≤ |B| if |A_ij| ≤ |B_ij| for all i, j. It is straightforward to see that (see Chap. 2.7 of Golub and Van Loan²⁰)

\tilde{A} = A + F

⁠, where |F| < ϵ|A|. In addition,

(\tilde{A} - \tilde{λ} I) \tilde{v} = {[(\tilde{A} - \tilde{λ} I) \tilde{v}]}_{c o m p} + g

⁠, where

| g | \leq (1.01) n ϵ | (\tilde{A} - \tilde{λ} I) | | \tilde{v} |

⁠. Now, note that

(A - \tilde{λ} I) \tilde{v} = (\tilde{A} - \tilde{λ} I) \tilde{v} - F \tilde{v}

so that

‖ \tilde{r} ‖_{2} \leq {‖{[(\tilde{A} - \tilde{λ} I) \tilde{v}]}_{c o m p}‖}_{2} + ‖ g ‖_{2} + ‖ F \tilde{v} ‖_{2} .

Noting that

‖ g ‖_{2} = (1.01) n ϵ ‖ | \tilde{A} - \tilde{λ} I | | \tilde{v} | ‖_{2} \leq (1.01) n^{5 / 2} ϵ ‖ \tilde{A} - \tilde{λ} I ‖_{\max} ‖ \tilde{v} ‖_{\infty},

‖ F \tilde{v} ‖_{2} \leq ϵ ‖ | A | | \tilde{v} | ‖_{2} \leq n ϵ ‖ A ‖_{\max} ‖ \tilde{v} ‖_{\infty},

and

{‖{[(\tilde{A} - \tilde{λ} I) \tilde{v}]}_{c o m p}‖}_{2} \leq n^{1 / 2} {‖{[(\tilde{A} - \tilde{λ} I) \tilde{v}]}_{c o m p}‖}_{\max},

the result is proved.□

We now prove estimate (35). Applying Lemma E.2 to the set ${{\tilde{v}}_{j}}_{1 \leq j \leq m}$ yields an orthonormal set ${{\hat{v}}_{j}}_{1 \leq j \leq m}$ such that $‖ {\hat{v}}_{j} - {\tilde{v}}_{j} ‖_{2} \leq 2^{- 1 / 2} m μ$ ⁠, where μ is as in (E1). By Lemma E.3, we have that

‖ {\hat{r}}_{j} ‖_{2} \leq 2^{- 1 / 2} m (‖ A ‖_{2} + | {\tilde{λ}}_{j} |) μ + ‖ {\tilde{r}}_{j} ‖_{2}, 1 \leq j \leq m .

The estimate now follows easily upon applying Lemma E.4 and taking the sup over j.□

3. Proof of Theorem IV.4

Proof of Theorem IV.4.

The proof is a simple consequence of the min–max characterization of eigenvalues of Hermitian matrices. By min–max (here, U denotes a subspace of

C^{n}

⁠),

| λ_{j}^{α} - λ_{j}^{α_{0}} | = \min_{\dim U = j} \max_{\begin{array}{c} v \in U \\ v \neq 0 \end{array}} |\frac{⟨(v| (A^{α} - A^{α_{0}}) v⟩}{⟨(v| v⟩}| .

On the other hand, for any fixed v ≠ 0, we have by Taylor’s theorem

|\frac{⟨(v| (A^{α} - A^{α_{0}}) v⟩}{⟨(v| v⟩}| \leq | α - α_{0} | \sup_{β \in [α_{0}, α]} |\partial_{β} \frac{⟨(v| (A^{β} - A^{α_{0}}) v⟩}{⟨(v| v⟩}| \leq | α - α_{0} | \sup_{β \in [α_{0}, α]} ‖ \partial_{β} A^{β} ‖_{2},

and the result follows immediately.□

4. Proof of Proposition IV.7

Proof of Proposition IV.7.

Differentiating

H_{Ξ}^{α}

yields

\partial_{α} H_{Ξ}^{α} = (- \partial_{α} Q^{α}) P_{Ξ} H^{α} P_{Ξ} Q^{α, ⊥} + Q^{α, ⊥} P_{Ξ} H^{1} P_{Ξ} Q^{α, ⊥} + Q^{α, ⊥} P_{Ξ} H^{α} P_{Ξ} (- \partial_{α} Q^{α}) .

For α < 1, we have ‖P_ΞH^αP_Ξ‖₂ ≤ 10, and ‖H¹‖₂ ≤ 3. It remains only to estimate ‖∂_αQ^α‖₂. Using the Dirac notation to represent

L_{K, 1}^{2}

-projections, we have

Q^{α} = |Ψ^{(8)}⟩ ⟨Ψ^{(8)}| = |\sum_{m = 0}^{8} α^{m} Ψ^{m}⟩ ⟨\sum_{n = 0}^{8} α^{n} Ψ^{n}| = \sum_{m = 0}^{8} \sum_{n = 0}^{8} α^{m + n} |Ψ^{m}⟩ ⟨Ψ^{n}|

so that

\partial_{α} Q^{α} = \sum_{m = 0}^{8} \sum_{n = 0}^{8} (m + n) α^{m + n - 1} |Ψ^{m}⟩ ⟨Ψ^{n}| .

Using

‖ |Ψ^{m}⟩ ⟨Ψ^{n}| ‖_{2} \leq ‖ Ψ^{m} ‖_{2} ‖ Ψ^{n} ‖_{2}

and

\max_{0 \leq j \leq 8} ‖ Ψ^{j} ‖_{2} \leq \sqrt{3}

by Proposition C.3, we have, for α ≤ 1,

‖ \partial_{α} Q^{α} ‖_{2} \leq 3 \sum_{m = 0}^{8} \sum_{n = 0}^{8} (m + n) = 1944 .

Putting everything together, we conclude

\sup_{0 \leq α \leq \frac{7}{10}} ‖ \partial_{α} H_{Ξ}^{α} ‖_{2} \leq 2 \times 10 \times 1944 + 3 = 38 883 .

□

APPENDIX F: PROOF OF PROPOSITION II.1

We can now prove Proposition II.1. We start by proving (11).

1. Proof of (11)

We now prove (11). It is straightforward to derive

⟨(\sum_{n = 0}^{8} α^{n} Ψ^{n}| \sum_{n = 0}^{8} α^{n} Ψ^{n}⟩ = \sum_{n = 0}^{8} \sum_{j = 0}^{n} ⟨(Ψ^{j}| Ψ^{n - j}⟩ α^{n} + \sum_{n = 0}^{7} \sum_{j = 0}^{n} ⟨(Ψ^{8 - j}| Ψ^{8 - (n - j)}⟩ α^{16 - n} .

(F1)

We now make two observations which simplify the computation. First, recall that the operator $- P^{⊥} {(H^{0})}^{- 1} P^{⊥} H^{1}$ maps $L_{K, 1, 1, A}^{2} \to L_{K, 1, 1, B}^{2}$ and $L_{K, 1, 1, B}^{2} \to L_{K, 1, 1, A}^{2}$ ⁠. It follows that $Ψ^{0} \in L_{K, 1, 1, A}^{2}$ ⁠, $Ψ^{1} \in L_{K, 1, 1, B}^{2}$ ⁠, $Ψ^{2} \in L_{K, 1, 1, A}^{2}$ ⁠, and so on, and hence,

⟨(Ψ^{2 i}| Ψ^{2 j + 1}⟩ = 0 \forall i, j \in {0, 1, 2, \dots} .

It follows that all terms in (F1) with odd powers of α vanish. Second, note that since Ψ⁰ ∈ ranP, while Ψⁿ ∈ ranP^⊥ for all n ≥ 1, we have that

⟨(Ψ^{n}| Ψ^{0}⟩ = ⟨(Ψ^{0}| Ψ^{n}⟩ = 0 \forall n \in {1, 2, \dots} .

Deriving (11) is then just a matter of computation using the properties of the chiral basis. For the leading term, we have

⟨(Ψ^{0}| Ψ^{0}⟩ = ⟨(χ^{\tilde{0}}| χ^{\tilde{0}}⟩ = 1 .

For the α² term, the only non-zero term is

⟨(Ψ^{1}| Ψ^{1}⟩ = ⟨(- \sqrt{3} i χ^{\tilde{q_{1}}, 1}| - \sqrt{3} i χ^{\tilde{q_{1}}, 1}⟩ = 3,

using (C4). For the α⁴ term, the possible non-zero terms are

⟨(Ψ^{3}| Ψ^{1}⟩ + ⟨(Ψ^{2}| Ψ^{2}⟩ + ⟨(Ψ^{1}| Ψ^{3}⟩,

but Ψ³ and Ψ¹ depend on orthogonal chiral basis vectors [see (C4) and (C6)], so we are left with

⟨(Ψ^{2}| Ψ^{2}⟩ = ⟨((\frac{\sqrt{3} - i}{2}) χ^{\tilde{- b_{1}}, 1} + (\frac{\sqrt{3} + i}{2}) χ^{\tilde{- b_{2}}, 1}| (\frac{\sqrt{3} - i}{2}) χ^{\tilde{- b_{1}}, 1} + (\frac{\sqrt{3} + i}{2}) χ^{\tilde{- b_{2}}, 1}⟩ = 2,

using (C5) and orthgonality of $χ^{\tilde{- b_{1}}, 1}$ and $χ^{\tilde{- b_{2}}, 1}$ ⁠. We omit the derivation of the higher terms since the derivations do not require any new ideas.

2. Proof of (10)

It is straightforward to derive

\begin{align} ⟨(\sum_{n = 0}^{8} α^{n} Ψ^{n *} (- r)| \sum_{n = 0}^{8} α^{n} Ψ^{n} (r)⟩ = & \sum_{n = 0}^{8} \sum_{j = 0}^{n} ⟨(Ψ^{j *} (- r)| Ψ^{n - j} (r)⟩ α^{n} \\ + \sum_{n = 0}^{7} \sum_{j = 0}^{n} ⟨(Ψ^{8 - j *} (- r)| Ψ^{8 - (n - j)} (r)⟩ α^{16 - n} . \end{align}

(F2)

We now note the following.

Proposition F.1.

Let χ be a chiral basis function in $L_{K, 1, 1}^{2}$ . Then, χ^*(−r) = χ(r).

Proof.

The proof follows immediately from the explicit forms of the chiral basis functions in $L_{K, 1, 1}^{2}$ given by (B1)–(B2)–(B3) and the observation that for any $k \in R^{2}$ ⁠, ${(e^{i k \cdot (- r)})}^{*} = e^{i k \cdot r}$ ⁠.□

Using Proposition F.1 and the same two observations as in the Appendix F 1, we have that the only non-zero terms in (F2) are those with even powers of α and that other than the leading term, terms involving Ψ⁰ do not contribute. The calculation is then similar to the previous case. For the leading order term, we have

⟨(Ψ^{0 *} (- r)| Ψ^{0} (r)⟩ = ⟨(χ^{\tilde{0}}| χ^{\tilde{0}}⟩ = 1 .

The only non-zero α² term is

⟨(Ψ^{1 *} (- r)| Ψ^{1} (r)⟩ = ⟨(\sqrt{3} i χ^{\tilde{q_{1}}, 1}| - \sqrt{3} i χ^{\tilde{q_{1}}, 1}⟩ = - 3 .

The only non-zero α⁴ term is

\begin{aligned} ⟨(Ψ^{2 *} (- r)| Ψ^{2} (r)⟩ \\ = ⟨((\frac{\sqrt{3} + i}{2}) χ^{\tilde{- b_{1}}, 1} + (\frac{\sqrt{3} - i}{2}) χ^{\tilde{- b_{2}}, 1}| (\frac{\sqrt{3} - i}{2}) χ^{\tilde{- b_{1}}, 1} + (\frac{\sqrt{3} + i}{2}) χ^{\tilde{- b_{2}}, 1}⟩ \\ = {(\frac{\sqrt{3} - i}{2})}^{2} + {(\frac{\sqrt{3} + i}{2})}^{2} = 1 . \end{aligned}

We omit the derivation of the higher terms since the derivations do not require any new ideas.

Proposition IV.1 implies that the series expansion of ψ^α exists up to any order. We can therefore define formal infinite series by

⟨(\sum_{n = 0}^{\infty} α^{n} Ψ^{n *} (- r)| \sum_{n = 0}^{\infty} α^{n} Ψ^{n} (r)⟩,

(F3)

⟨(\sum_{n = 0}^{\infty} α^{n} Ψ^{n}| \sum_{n = 0}^{\infty} α^{n} Ψ^{n}⟩ .

(F4)

We then have the following.

Proposition F.2.

Expansions (10) and (11) approximate the formal series (F3) and (F4) up to terms of order α¹⁰.

Proof.

The series agree exactly without any simplifications up to terms of α⁹. However, because the even and odd terms in the expansion of ψ^α are orthogonal (since they lie in $L_{K, 1, 1, A}^{2}$ and $L_{K, 1, 1, B}^{2}$ ⁠, respectively), all terms with odd powers of α vanish in expansions (F3) and (F4). The series may disagree at order α¹⁰ because the infinite series includes terms arising from inner products of Ψ¹ and Ψ⁹.□

REFERENCES

1.

R.

Bistritzer

and

A. H.

MacDonald

,

Proc. Natl. Acad. Sci. U. S. A.

108

,

12233

(

2011

); arXiv:1009.4203.

https://doi.org/10.1073/pnas.1108174108

Crossref

PubMed

2.

P.

San-Jose

,

J.

González

, and

F.

Guinea

,

Phys. Rev. Lett.

108

,

216802

(

2012

); arXiv:1110.2883.

https://doi.org/10.1103/physrevlett.108.216802

Crossref

PubMed

3.

Y.

Cao

,

V.

Fatemi

,

S.

Fang

,

K.

Watanabe

,

T.

Taniguchi

,

E.

Kaxiras

, and

P.

Jarillo-Herrero

,

Nature

556

,

43

(

2018

).

https://doi.org/10.1038/nature26160

Crossref

PubMed

4.

G.

Tarnopolsky

,

A. J.

Kruchkov

, and

A.

Vishwanath

,

Phys. Rev. Lett.

122

,

106405

(

2019

); arXiv:1808.05250.

https://doi.org/10.1103/PhysRevLett.122.106405

Crossref

PubMed

5.

S.

Becker

,

M.

Embree

,

J.

Wittsten

, and

M.

Zworski

, “

Mathematics of magic angles in a model of twisted bilayer graphene

,” arXiv:2008.08489 [math-ph] (

2020

).

Google Scholar

6.

S.

Becker

,

M.

Embree

,

J.

Wittsten

, and

M.

Zworski

,

Phys. Rev. B

103

,

165113

(

2020

); arXiv:2010.05279.

https://doi.org/10.1103/PhysRevB.103.165113

Crossref

7.

G.

Catarina

,

B.

Amorim

,

E. V.

Castro

,

J. M. V. P.

Lopes

, and

N.

Peres

,

Handbook of Graphene

(

Wiley

,

2019

), pp.

177

–

231

.

Google Scholar

Crossref

8.

J.

Oliver

,

J. Comput. Appl. Math.

5

,

85

(

1979

).

https://doi.org/10.1016/0771-050x(79)90002-0

Crossref

9.

A.

Meurer

,

C. P.

Smith

,

M.

Paprocki

,

O.

Čertík

,

S. B.

Kirpichev

,

M.

Rocklin

,

A.

Kumar

,

S.

Ivanov

,

J. K.

Moore

,

S.

Singh

,

T.

Rathnayake

,

S.

Vig

,

B. E.

Granger

,

R. P.

Muller

,

F.

Bonazzi

,

H.

Gupta

,

S.

Vats

,

F.

Johansson

,

F.

Pedregosa

,

M. J.

Curry

,

A. R.

Terrel

,

Š.

Roučka

,

A.

Saboo

,

I.

Fernando

,

S.

Kulal

,

R.

Cimrman

, and

A.

Scopatz

,

PeerJ Comput. Sci.

3

,

e103

(

2017

).

https://doi.org/10.7717/peerj-cs.103

Crossref

10.

D.

Massatt

,

S.

Carr

,

M.

Luskin

, and

C.

Ortner

,

Multiscale Model. Simul.

16

,

429

(

2018

).

https://doi.org/10.1137/17m1141035

Crossref

11.

L.

Zou

,

H. C.

Po

,

A.

Vishwanath

, and

T.

Senthil

,

Phys. Rev. B

98

,

085435

(

2018

); arXiv:1806.07873.

https://doi.org/10.1103/physrevb.98.085435

Crossref

12.

H. C.

Po

,

L.

Zou

,

A.

Vishwanath

, and

T.

Senthil

,

Phys. Rev. X

8

,

031089

(

2018

); arXiv:1803.09742.

https://doi.org/10.1103/PhysRevX.8.031089

Crossref

13.

H. C.

Po

,

L.

Zou

,

T.

Senthil

, and

A.

Vishwanath

,

Phys. Rev. B

99

,

195455

(

2019

); arXiv:1808.02482.

https://doi.org/10.1103/physrevb.99.195455

Crossref

14.

T.

Kato

,

Perturbation Theory for Linear Operators

(

Springer-Verlag Berlin Heidelberg

,

1995

), Vol. 132.

Google Scholar

Crossref

15.

H. L.

Cycon

,

R. G.

Froese

,

W.

Kirsch

, and

B.

Simon

,

Schrödinger Operators: With Applications to Quantum Mechanics and Global Geometry

(

Springer-Verlag Berlin Heidelberg

,

1987

).

Google Scholar

Crossref

16.

M. E.

Taylor

,

Partial Differential Equations I

(

Springer

,

New York, NY

,

1996

).

Google Scholar

Crossref

17.

S. M.

Rump

,

Acta Numer.

19

,

287

–

449

(

2010

).

https://doi.org/10.1017/s096249291000005x

Crossref

18.

A.

Messiah

,

Quantum Mechanics: Volume II

(

North-Holland Publishing Company

,

Amsterdam

,

1962

).

Google Scholar

19.

B. N.

Parlett

,

The Symmetric Eigenvalue Problem

(

Society for Industrial and Applied Mathematics

,

1998

).

Google Scholar

Crossref

20.

G. H.

Golub

and

C. F.

Van Loan

,

Matrix Computations

, 4th ed. (

The Johns Hopkins University Press

,

2013

).

Google Scholar

Crossref

2021

Author(s)

Existence of the first magic angle for the chiral model of bilayer graphene

I. INTRODUCTION

A. Outline

B. Code availability

II. STATEMENT OF RESULTS

A. Tarnopolsky–Kruchkov–Vishwanath’s chiral model

B. Rigorous justification of TKV’s formal expansion of the Fermi velocity and proof of existence of first magic angle

C. Structure of paper

III. SYMMETRIES, BLOCH THEORY, AND ZERO MODES OF TKV’s CHIRAL MODEL

A. Symmetries of the TKV model

B. Bloch theory for the TKV Hamiltonian

C. Zero modes of the chiral model

IV. RIGOROUS JUSTIFICATION OF TKV’s EXPANSION OF THE FERMI VELOCITY

A. Alternative formulation of TKV’s expansion

B. Rigorous error estimates for the expansion of the moiré K point Bloch function

SUPPLEMENTARY MATERIAL

ACKNOWLEDGMENTS

DATA AVAILABILITY

APPENDIX A: DERIVATION OF EXPRESSION FOR FERMI VELOCITY IN TERMS OF LK,1,12 ZERO MODE OF Hα

APPENDIX B: THE CHIRAL BASIS OF LK,12 AND ACTION OF H0 AND H1 WITH RESPECT TO THIS BASIS

1. The spectrum and eigenfunctions of H0 in LK2

2. The spectrum and eigenfunctions of H0 in LK,12

3. The chiral basis of LK,12

4. The spectrum of H1 in LK2 and LK,12

5. The action of H1 on LK,12 with respect to the chiral basis

APPENDIX C: FORMAL EXPANSION OF THE ZERO MODE

APPENDIX D: PROOF OF PROPOSITION IV.4

APPENDIX E: PROOF OF PROPOSITION IV.5

1. Proof of Theorem IV.2

2. Proof of Theorem IV.3

3. Proof of Theorem IV.4

4. Proof of Proposition IV.7

APPENDIX F: PROOF OF PROPOSITION II.1

1. Proof of (11)

2. Proof of (10)

REFERENCES

Supplementary Material

Citing articles via

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

This Feature Is Available To Subscribers Only

APPENDIX A: DERIVATION OF EXPRESSION FOR FERMI VELOCITY IN TERMS OF $L_{K, 1, 1}^{2}$ ZERO MODE OF H^α

APPENDIX B: THE CHIRAL BASIS OF $L_{K, 1}^{2}$ AND ACTION OF H⁰ AND H¹ WITH RESPECT TO THIS BASIS

1. The spectrum and eigenfunctions of H⁰ in $L_{K}^{2}$

2. The spectrum and eigenfunctions of H⁰ in $L_{K, 1}^{2}$

3. The chiral basis of $L_{K, 1}^{2}$

4. The spectrum of H¹ in $L_{K}^{2}$ and $L_{K, 1}^{2}$

5. The action of H¹ on $L_{K, 1}^{2}$ with respect to the chiral basis