Moir\'e ordered current loops in the graphene twist bilayer

TL;DR

This paper reveals that twisted bilayer graphene under a magnetic field hosts an ordered array of persistent current loops across the material, driven by moiré patterns and broken time-reversal symmetry.

Contribution

It demonstrates the existence of moiré-scale current loops in twisted bilayer graphene under magnetic fields, a novel intrinsic feature linked to the moiré pattern and broken time-reversal symmetry.

Findings

Current loops form throughout the material at small magnetic fields.

Current loops exist in non-equilibrium transport states, flowing against bias.

The loops are tied to the moiré pattern and broken time-reversal symmetry.

Abstract

While a typical material exhibits field induced currents only at the boundary, a uniform out-of-plane magnetic field applied to two mutually rotated layers of graphene is shown to result in an ordered array of permanent current loops throughout the material. Each current loop consists of an interlayer current flowing through the open AA stacked regions of the moir\'e created by rotation, which then flows back through the neighboring AB regions to form a circuit, with significant current strength even at small fields. Similar moir\'e ordered arrays of current loops are also shown to exist in non-equilibrium transport states, where they manifest as current back flowing against the applied bias in the device. Such current loops thus represent an intrinsic feature of the twist bilayer in conditions of broken time reversal symmetry, and exist both as a low field imprint of the moir\'e lattice on Landau physics, and as measurable moir\'e scale current configurations in transport states.

Figures

c. However, one should note that the width of the ribbon here is only $\sim 4l_{B}$ and a “bulk” Landau level structure is not seen in this geometry.

c. However, one should note that the width of the ribbon here is only $\sim 4l_{B}$ and a “bulk” Landau level structure is not seen in this geometry.

Equations

$$t_{ij}\to t_{ij}\text{ exp}\left\{i\frac{e}{\hslash }{\int }_{{\mathbf{r}}_i}^{{\mathbf{r}}_j}\mathbf{A}\cdot d\mathbf{r}\right\},$$

$$I_{ij}=\frac{2e}{\hslash }\Im \sum _{\gamma {\gamma }^{\prime }}{H}_{j\gamma ,i{\gamma }^{\prime }}{c}_{i{\gamma }^{\prime }}{c}_{j\gamma }^{\star },$$

$$\mathbf{j}\left({\mathbf{R}}_i\right)=\sum _j{I}_{ij}{\hat{\mathbf{r}}}_{ij}$$

$$V\left({\mathbf{R}}_i-{\mathbf{R}}_j\right)=A\text{exp}\left(-B{\left|{\mathbf{R}}_i-{\mathbf{R}}_j\right|}^2\right),$$

Full text

Introduction. For a typical material the magnetic length exceeds by many orders of magnitude the lattice parameter. The material can thus be approximated as homogeneous, leading to the well known non-relativistic ($E_{n}\propto n+1/2$) or Dirac-Weyl ($E_{n}\propto\sqrt{n}$) Landau level spectra. This situation may change in a layered two dimensional system, for which a mutual rotation between the layers can result in a long range stacking modulation known as a moiré lattice. Such materials have a moiré lattice parameter comparable in magnitude to the magnetic length of experimentally achievable magnetic fields and, as a consequence, dramatically different physics in the presence of a magnetic field1, 2, 3. Most notably this can be seen through the fact that the $B$ or $\sqrt{B}$ field dependence of the standard non-relativistic or Dirac-Weyl spectra goes over to a fractal dependence on the field strength (the so called Hofstadter butterfly4).

The graphene twist bilayer5, 6, 7 is perhaps the most studied example of a two dimensional moiré. This system exhibits a remarkable range of electronic structure phenomena, encompassing both graphene like (large angle) and charge confined (small angle) limits8, 9, 10. For small twist angles the moiré is comparable to the magnetic length of achievable (1-5 Tesla) experimental fields, and indeed for a $\theta=0.52^{\circ}$ twist bilayer signatures of the Hofstadter butterfly were recently observed in experiment11.

In this work we show that in addition to displaying the Hofstadter butterfly, the twist bilayer possesses a low field and intermediate angle pre-Hofstadter phase in which although the Landau spectrum has the $E_{n}\propto\sqrt{n}$ Dirac-Weyl form of graphene, the wave functions carry a remarkable imprint of the moiré in the form of an array of permanent current loops, with one current loop per moiré cell. This current loop lattice occurs at significantly weaker fields than those required to observe the Hofstadter butterfly: for example, for a twist angle of $6^{\circ}$ it is found already at a field of the order of 1 Tesla, while a field of $\approx$118 Tesla would be required to observe the Hofstadter butterfly. While such current loops might be difficult to observe in equilibrium conditions – although they will certainly impact intercalation of the bilayer – we show they exist also in transport states, both with and without an external magnetic field. To see this we gate a single layer of graphene and overlay this with a twisted layer. This system exhibits an ordered array of current loops aligned in the bias direction, resulting in remarkable back flowing current against the bias direction in the ungated layer. Such ordered arrays of current loops should thus be observable in transport experiments under ballistic conditions, for instance with scanning probe techniques12, 13, 14.

Results. In order to address physics in the pre-Hofstadter regime, we must first establish a criteria for the critical field strength at which the Hofstadter butterfly is first seen. Such a criteria has been suggested by Moon and Koshino15 to be $l_{B}<D$, with $l_{B}=\sqrt{\hbar/(eB)}$ the magnetic length and $D=a/(2\sin\theta/2)$ the moiré length ($a$ is the lattice constant of graphene). From the semi-classical point of view this is an intuitive result, as once the magnetic length exceeds the underlying moiré periodicity it is effectively averaged out over one cyclotron orbit, and thus one expects that quasiparticles will feel only a structureless average of the moiré. We will therefore first numerically investigate whether this criteria represents a useful marker of the pre-Hofstadter regime.

For this purpose the width of the zero mode Landau level provides a useful guide as (i) at the onset of the Hofstadter butterfly the degeneracy of this mode will be broken, and thus will acquire a finite width and, (ii), the zero mode is the last of the Landau levels to be destroyed by the clustering of van Hove singularities towards the Dirac point as $\theta\to 0$9, 10. We therefore consider four different twist angles in the moiré regime of $\theta<15^{\circ}$ and plot the zero mode width as a function of applied field. As shown in Fig. 1, for all four twist angles the onset of zero mode broadening is always at $l/D\approx 1$, confirming the intuitive criteria proposed in Ref. 15. Note that the zero mode broadening curves are seen to be, to a good approximation, a function only of the ratio $l/D$ suggesting that the lattice constant $a$ no longer plays an important role at the large $D$ found at these twist angles.

The current density of the ground state zero mode wave function, as we now show, turns out to be profoundly sensitive to the moiré lattice even in the pre-Hofstadter regime of $l_{B}>>D$. In Fig. 2 we present the current density obtained by summing over all states in the zero mode of a twist bilayer nanoribbon with $\theta=3.89^{\circ}$. Shown are the current densities for three different fields: $B=3~{}\mathrm{T}$ [panels (a-d)], $B=10~{}\mathrm{T}$ [panels (e-h)], and $B=109~{}\mathrm{T}$ [panels (i-l)]. The band structures and corresponding density of states for these systems can be found in SI; one can note that the moiré induced van Hove singularity is at $\sim\pm 0.3$eV, greater than $\hbar\omega_{c}$ except for very large fields. The first two of these fields have $l_{B}/D$ ratios of 4.14 and 2.23 respectively, and exhibit clear $\sqrt{n}$ Dirac-Weyl Landau spectra (see SI). In contrast, the $B=109~{}\mathrm{T}$ field has $l_{B}/D=0.69$ and a zero mode broadening greater than $\hbar\omega_{c}$ (see SI): this is clearly in the Hofstadter regime. In all cases, however, the current density of the zero mode is qualitatively the same and exists throughout the bulk of the material. Remarkably, therefore, the Landau wave functions exhibit a current loop lattice at magnetic fields for which the Landau spectra shows no trace of the presence of the moiré.

We now characterize more closely the current density. This consists of permanent current loops in the bulk of the material in which an interlayer charge flows through the AA spots, and back through the neighboring AB spots of the lattice in the opposite direction. As may be seen from panels (d), (h), and (l), the interlayer current is almost zero in the intermediate stacking regions with, however, an in-plane current such that a closed loop is formed between the “source” and “sink” regions of interlayer current in each layer. The topology of the current flow thus has the structure of a “convection cell” attached to each moiré spot of the lattice. Furthermore this current loop lattice, irrespective of the strength of the field, has a lattice parameter equal to the moiré lattice parameter; a fact evident from the similar form of the current density seen in all three panels of Fig. 2 even though the magnetic field changes by more than an order of magnitude.

Note that although we have shown the current density for a zero mode filling factor of one, a similar structure of the current density is found for any finite filling factor, although for low filling factors “vacancies” appear in the loop lattice (see SI). This simply reflects the fact that individual eigenstates all exhibit a similar current density, modulated however by an exponential envelope on the length scale of $l_{B}$. Summing these over the zero mode then leads to the current shown in Fig. 2 in which only the presence of the length scale $D$ can be seen.

We now consider the strength of these current loops. As the interlayer current averaged over the bulk of the nanoribbon must obviously be zero, to characterize its strength we instead sum over all currents within a single AA spot of the moiré, within which the interlayer current always flows in the direction of the applied field. The field dependence of this AA current grows to a very good approximation linearly in all cases (see Fig. 3(a)), provided $l_{B}>D$. This arises from the fact that (i) the degeneracy of the zero mode increases linearly with field and, (ii), as noted above, all eigenstates in the zero mode have a similar current density. These net currents are significant in size, of the order of 10-100 $\mu$A, larger by an order of magnitude than the currents typically found in driven nanostructures16, 17, 18, 19, 20. To investigate how the moiré impacts the intrinsic current strength in Fig. 3(b) we present the current density of the AA spots, i.e. the total current through the AA region divided by its area, as a function of twist angle for several fields. Interestingly, the current density increases significantly as the twist angle decreases, with for the 5.9 Tesla field strength an almost order of magnitude increase in the current density found on reducing the twist angle from 12∘ to $2^{\circ}$.

Closed current loops in equilibrium are likely to be difficult to detect experimentally and so we now turn to the question of whether signatures of this unusual ground state exist in transport states, for which we employ a non-equilibrium Green’s function technique (see the Methods section). Our transport setup is shown in Fig. 4(a) and consists of a single layer graphene nanoribbon junction (width $\sim 100$ nm) overlaid by a rotated layer ($\theta=3.89^{\circ}$) not connected to the leads. Current flow in the top layer will therefore arise purely from the interlayer interaction.

In Fig. 4(b-d) we show the local currents in the area of the junction indicated in Fig. 4(a) for three representative low energy scattering states. Each of (b-d) display in the topmost panels the in-plane currents in layers 1 and 2, and in the bottom panels the interlayer current density (left hand panel) with the lattice sites displayed as small black dots, and in the right hand panel the left-to-right transmission function of the junction. In this latter panel the vertical line indicates the energy of the scattering state for which the currents are shown.

We first consider transport in the absence of the out-of-plane magnetic field, Fig. 4(b). We note that the transmission function is lowered by the presence of the overlaid twist layer, in agreement with recent finite temperature Boltzmann calculations of the twist bilayer21. However most interesting is the structure of the electron currents in the overlaid layer which flow against the transport direction, a remarkable situation that arises from the formation of a moiré ordered array of interlayer current loops. Contrasting Figs. 4(b) and 2 one notes that this loop lattice no longer possesses the approximate $C_{3}$ symmetry of the underlying moiré, with the current loops now centered on the AA regions of the lattice. For transport states at higher energies the loop lattice disorders, although strong interlayer currents remain a feature of the transport current. Thus the ordered array of current loops found in the ground state at finite magnetic field persists, in modified form, into low energy transport states in the absence of the magnetic field. This indicates that such current loops are intrinsic to the twist bilayer in situations of broken time reversal symmetry, and is consistent with the independence of the lattice parameter of the equilibrium current loop lattice from the magnetic length. In a low energy picture the twist bilayer can be described as two Dirac-Weyl sub-systems coupled by a non-abelian interlayer moiré potential22, 10. As this potential is known to drive an interlayer precession of quasi-particles22, we conclude that transport is governed by an interplay of the applied bias potential and the moiré potential which, in the low energy transport states, manifests as ordered current loop lattices.

In the presence of an out-of-plane magnetic field similar current loop structures can be observed in transport states, see Fig. 4(c-d), with, interestingly, those that block longitudinal transport showing a similar interlayer current modulation to that found in the ground state and those that transmit a structure closer to the $B=0$ current modulation. In fact, in the blocked state the current flow close to the edge of the gated layer is, due to the current loop array, against the direction of the applied bias, see Fig. 4c. However, one should note that the width of the ribbon here is only $\sim 4l_{B}$ and a “bulk” Landau level structure is not seen in this geometry.

Discussion. We have uncovered imprints of the moiré in the current densities of both equilibrium and non-equilibrium states of the twist bilayer. In equilibrium the zero mode Landau level exhibits an ordered array of permanent current loops throughout the bulk, consisting of a “convection cell” current configuration attached to each moiré spot of the lattice. Such current loop lattices, in modified form, are also found in transport states near the Dirac point, as seen in the low energy scattering states of a graphene nanoribbon junction overlaid by a twisted layer. Remarkably, here the current loop array orders such that current in the second layer – which has a magnitude of $\sim 50$% the current strength in the gated layer – flows against the transport direction. Loop lattice transport states are found for both zero and finite magnetic field, showing the phenomena to be intrinsic to the twist bilayer in situations of broken time reversal symmetry. The physics of the equilibrium loop lattice is therefore underpinned by the interlayer precession property of the moiré potential22, rather than the physics of a confining potential of length scale comparable to the magnetic length which drives bulk currents in 2DEGs at rough semi-conductor interfaces23, 24. While closed current loops in equilibrium are likely to be difficult to observe in experiment, although they may interact in interesting ways with ordered arrays of intercalated impurities for which the moire forms a natural geometry25, the persistence of such current loops in transport states shows that this phenomena should be observable by local probes such as scanning gate microscopy.

Method and computational details. In order to study the twist bilayer in an out-of-plane magnetic field we employ a nanoribbon geometry that allows free choice of the magnetic field strength, as well as allowing for the study of edge currents. For all systems studied we ensure that the nanoribbon width is at least an order of magnitude larger than the magnetic length in order to achieve a good bulk. The large system sizes that small angle twist bilayers nanoribbons then inevitably entail - our nanoribbon unit cells typically contain up to 40,000 carbon atoms - necessitates the use of a semi-empirical tight-binding method. We will use the tight-binding method of Ref. 8 which was deployed in that work for the study of twist bilayer flakes; it consists of the environment dependent method of Tang et al.26, but re-parameterized by performing a least squares fit to the high symmetry eigenvalues from a number of small unit cell few layer graphene systems generated ab-initio; for details we refer the reader to Ref. 8.

We introduce an external magnetic field into our calculations via the standard Peierls substitution27 as an additional phase in the two-center hopping integrals:

[TABLE]

with a Landau gauge $\mathbf{A}=\left(0,Bx,0\right)$ such that $\mathbf{\nabla}\times\mathbf{A}=B\mathbf{\hat{z}}$ (the graphene plane taken to lie in the $xy$ Cartesian plane).

For calculating the electron and hole currents we follow the derivation of T. Todorov28 and determine the current from atom $n^{\prime}$ to atom $n$ as:

[TABLE]

where $H$ is the real space tight-binding Hamiltonian, $c_{i\gamma}$ are the expansion coefficients of the eigenvectors in the tight-binding basis, and $\gamma$ is the orbital index. The current at each site in the lattice is found from

[TABLE]

with $\hat{\mathbf{r}}_{ij}$ the unit vector pointing from site $i$ to site $j$.

The problem of finding a commensuration cell of the twist bilayer leads to a discrete set of unit cells labeled by two integers, which we will refer to as $p$ and $q$ following the notation of Refs. 29, 9, 10. For a nanoribbon geometry, these unit cells evidently define the edge structure which, as might be expected, is always a low index facet. For $p=1$ nanoribbons the edge comprises one zigzag segment per nanoribbon cell edge, with the remainder of the edge having the armchair type termination.

For the transport calculations presented in Fig. 4, we have simplified the tight-binding model to the two-center hopping function

[TABLE]

and constructed the geometry and Hamiltonian using a code specifically written for this purpose. We have then converted them to the format of the TBtrans code 30, which was used to calculate the transmissions and orbital currents using Green’s function formalism. The model parameters have been chosen such that they reproduce the DOSs for AA- and AB-stacked bilayer graphene as given in Ref. 31. We have then made use of the SISL code 32 to calculate the bond currents $I_{ij}$ from the lattice site $\mathbf{R}_{j}$ to $\mathbf{R}_{i}$. From the bond currents, we have calculated a vector field via Eq. (3), which has then been coarse grained in order to obtain the representation shown in Fig. 4.

Data availability. The data that support the findings of this study are available from the corresponding author upon request.

Author contributions. The project was framed by Shallcross and Meded, calculations were performed by Landgraf, Fleischmann, and Weckbecker, using a code written by Weckbecker and Shallcross. All authors contributed to the interpretation of results and writing of the manuscript.

Competing financial interests. The authors declare no competing financial interests.

Acknowledgments. This work was supported by the Collaborative Research Center SFB 953 of the Deutsche Forschungsgemeinschaft (DFG).