Self-similar occurrence of massless Dirac particles in graphene under a magnetic field

Intricate interplay between the periodicity of the lattice structure and that of the cyclotron motion gives rise to a well-known self-similar fractal structure of the energy eigenvalue, known as the Hofstadter butterfly, for an electron moving in lattice under magnetic field. Connected with the $n=0$ Landau level, the central band of the Hofstadter butterfly is especially interesting in the honeycomb lattice. While the entire Hofstadter butterfly can be in principle obtained by solving Harper's equations numerically, the weak-field limit, most relevant for experiment, is intractable owing to the fact that the size of the Hamiltonian matrix, which needs to be diagonalized, diverges. In this paper, we develop an effective Hamiltonian method that can be used to provide an accurate analytic description of the central Hofstadter band in the weak-field regime. One of the most important discoveries obtained in this work is that massless Dirac particles always exist inside the central Hofstadter band no matter how small the magnetic flux may become. In other words, with its bandwidth broadened by the lattice effect, the $n=0$ Landau level contains massless Dirac particles within itself. In fact, by carefully analyzing the self-similar recursive pattern of the central Hofstadter band, we conclude that massless Dirac particles should occur under arbitrary magnetic field. As a corollary, the central Hofstadter band also contains a self-similar structure of recursive Landau levels associated with such massless Dirac particles. To assess the experimental feasibility of observing massless Dirac particles inside the central Hofstadter band, we compute the width of the central Hofstadter band as a function of magnetic field in the weak-field regime.

Figure 1

Schematic diagram for the gauge used in this work. The yellow parallelogram depicts a magnetic unit cell (MUC). Magnetic unit cells are denoted by the MUC index $\alpha$ along the $y$ direction. Different carbon sites within the same magnetic unit cell are distinguished by the dimer index $n$ and the A/B sublattice index. Note that horizontally connected A and B carbon sites share the same dimer index. Red arrows indicate the directions of the paths, along which nonzero phases are gained via the Peierls substitution. The value of the nonzero Peierls phase is written near each arrow while all the other phases are zero. This gauge is called the optimal gauge.

Figure 2

Hofstadter butterfly showing the energy eigenvalue, $E/t_0$, as a function of the magnetic flux per unit cell in units of magnetic flux quantum, $\varphi /{\varphi }_0$. Here, $t_0$ is the hopping amplitude in the absence of external magnetic field.

Figure 3

Contour plots for the energy dispersion at various flux values of $\varphi /{\varphi }_0=p/q$, with $p=1$ and $q$ increasing from 1 to 6 between panels (a)–(f). In the figure, the energy dispersions are normalized by their respective half bandwidth. The positions of the zero-energy momenta are denoted by little x marks and the magnetic Brillouin zones are enclosed by red solid lines. The above energy dispersions are computed by solving either the original Harper's equations or the effective Hamiltonian method explained in Sec. 4, both of which produce essentially the identical results. It is interesting to note that the effective Hamiltonian method works well even for $p/q=1$ owing to the mirror structure of the Hofstadter butterfly, which maps the region near $p/q=1$ to the weak-field counterpart.

Figure 4

Wave-function profiles for the zero-energy mode as a function of dimer index, $n$, at three different flux values: (a) $\varphi /{\varphi }_0=p/q=1/5$, (b) $1/20$, and (c) $1/100$. As one can see, at moderate flux values, say, $p/q=1/5$ and $1/20$, there is a sizable asymmetry around the maximum position. The asymmetry is seen more clearly in contrast to the Gaussian wave packet (red dashed lines) which is the exact energy eigenstate in the continuum, or weak-field, limit. It is important to note that, while the Gaussian wave packet provides a poor representation of the exact results (solid lines) at moderate flux values, a new analytic expansion method using the Clausen function (open circles) works very well for a wide range of flux values.

Figure 5

Schematic diagram for the construction of basis wave functions in the case of $p/q=3/152$. Basis wave functions for the central Hofstadter band can be constructed in three steps: (i) for a given momentum, $\mathbf{k}=(k_x,k_y)$, one generates a trial basis wave function according to the zero-energy formula in Eqs. (8) and (9), (ii) then, slices the so-obtained trial wave function into equally spaced $p$ pieces such that each piece contains exactly one local maximum, and (iii) finally, normalizes the $p$ piece-wise basis wave functions separately for each of the sublatttices A and B. Refer the text for details.

Figure 6

Overlap integral between the energy eigenstates obtained from the effective Hamiltonian method and the exact counterparts from the original Harper's equations. Stars indicate the averaged value of the overlap integral over all crystal momenta within the magnetic Brillouin zone, while circles denote individual results for different momenta. It is important to note that the overlap integral approaches unity very rapidly as $\varphi /{\varphi }_0=p/q$ decreases.

Figure 7

(a) Comparison between the exact energy dispersions obtained from the original Harper's equations (black dashed lines) and that from the effective Hamiltonian method (red solid line) for various $\varphi /{\varphi }_0=1/q$. Note that, with proper energy and momentum rescaling, all energy dispersions obtained from the effective Hamiltonian at different $\varphi /{\varphi }_0=1/q$ collapse into a single curve. In the figure, the energy dispersions are normalized by their respective half bandwidth $W$, and the momentum is expressed in units of $1/qa$. The inset shows the path in the magnetic Brillouin zone, along which the momentum is scanned. Note that the scanning path is chosen such that it passes through the Dirac points. (b) Comparison between the exact half bandwidth and that obtained from the effective Hamiltonian method as a function of $\varphi /{\varphi }_0=1/q$.

Figure 8

Comparison between the energy dispersions obtained from the original Harper's equations (black lines) and those from the effective Hamiltonian method (red lines) at various $\varphi /{\varphi }_0=p/q$. For clarity, the energy dispersions from the effective Hamiltonian method are plotted only within the window of $-\pi /qa\le \Delta k_y\le \pi /qa$. Note that the momentum is scanned along the same path as in Fig. 7a.

Figure 9

A sequence of zoomed views for the Hofstadter butterfly in graphene showing various self-similar recursive patterns. Note that a fan of narrow energy bands are emanated from each single-band boundary flux (SBF), ${\varphi }_{\mathrm{SBF}}$, which, as indicated by blue guiding curves, scale as $\mathrm{sgn}\left(n\right)\sqrt{\left|n\right(\varphi -{\varphi }_{\mathrm{SBF}}\left)\right|}$ with $n$ being an integer. This scaling behavior is a signature of the formation of recursive Landau levels associated with self-similarly occurring massless Dirac particles.