Geometry effects in topologically confined bilayer graphene loops

We investigate the electronic confinement in bilayer graphene by topological loops of different shapes. These loops are created by lateral gates acting via gap inversion on the two graphene sheets. For large-area loops the spectrum is well described by a quantization rule depending only on the loop perimeter. For small sizes, the spectrum depends on the loop shape. We find that zero-energy states exhibit a characteristic pattern that strongly depends on the spatial symmetry. We show this by considering loops of higher to lower symmetry (circle, square, rectangle and irregular polygon). Interestingly, magnetic field causes valley splittings of the states, an asymmetry between energy reversal states, flux periodicities and the emergence of persistent currents.

We investigate the electronic confinement in bilayer graphene by topological loops of different shapes. These loops are created by lateral gates acting via gap inversion on the two graphene sheets. For large-area loops the spectrum is well described by a quantization rule depending only on the loop perimeter. For small sizes, the spectrum depends on the loop shape. We find that zeroenergy states exhibit a characteristic pattern that strongly depends on the spatial symmetry. We show this by considering loops of higher to lower symmetry (circle, square, rectangle and irregular polygon). Interestingly, magnetic field causes valley splittings of the states, an asymmetry between energy reversal states, flux periodicities and the emergence of persistent currents.

I. INTRODUCTION
One of the most fundamental predictions of Quantum Mechanics is the existence of confined or bound states which are eigenstates of the Hamiltonian describing a particular system 1 . It is well known that the Schrödinger equation with a spatially dependent potential yields in many cases bound eigenstates located around the absolute potential minima. This potential binding mechanism around minima is crucial in many natural physical systems like, e.g., electrons in atoms and molecules 2 ; and it is also behind the formation of artificial bound states by potential gating with microelectrodes, as in semiconductor quantum dots 3 . In bilayer graphene (BLG) with Bernal stacking, the material we consider in this work, a similar potential confinement mechanism exists but it is related to the potential difference V a applied to split dual gates on the two graphene sheets; electrons are bound to the regions where V a is lowest. These regions of lowest V a can take different shapes, such as circular dots and rings for instance, and they have been intensively studied both theoretically and experimentally [4][5][6][7][8][9][10][11][12][13][14][15][16] . The technological interest of BLG nanostructures is twofold. First, they could serve as scalable spin qubits with long decoherence times due to the weak spin-orbit interaction in graphene [17][18][19] . Second, in combination with magnetic fields they could be utilized for valleytronic operations thanks to the manipulation of the valley degree of freedom [20][21][22] .
The possibility of an even more robust type of confinement in BLG was suggested in Ref. 23. This is a qualitatively different confinement, of a topological character, emerging near the domain wall that separates two regions where V a changes sign. The confined states are intrinsically 1D-like, with a characteristic transverse decay length, and they propagate along the domain wall with locked relative orientations of momentum and valley pseudospin. In the absence of valley mixing potentials states for the K and K' valleys propagate along the wall, without backscattering, showing opposite chiralities. This property is similar to the helical states within the quantum spin Hall effect [24][25][26] if we replace spin with the valley pseudospin of graphene. Now, let us suppose that the external gates are arranged in such a way that the domain wall closes into itself. A schematic representation is shown in Fig. 1a. Then, a loop forms supporting electronic bound states 27 and the following question naturally arises. What is the effect of different loop shapes onto the energy spectrum of these bound states? Being 1D-like, topologically confined sates in BLG are expected to essentially depend on the loop perimeter unlike the trivial potential confinement, which depends both on the surface and the shape of the confining region. Nevertheless, it is important to notice that the domain-wall states depicted with a red line in Fig. 1a are in reality quasi-1D with a finite spatial width. Therefore, the discrete energy spectrum will be in general a function of the system's geometrical symmetries, the region aspect ratio, etc. Below, we present a detailed investigation for topological loops of different shapes and symmetries. First, we show with the aid of an exactly solvable model that the energy spectrum is uniquely determined by the loop perimeter provided that the system size is large. Yet, for smaller loops finite size effects become important and the spectrum depends on the particular geometrical structure. This dependence is better seen in the sequence of zero energy crossing and anticrossings as the perimeter increases.
In presence of a magnetic field B, the dispersion is affected by the flux piercing the topological loop. This creates an energy splitting for the two valleys and an asymmetry with respect to energy inversion (E → −E) in the spectra as a function of perimeter. The B-dependence of the spectra with a fixed perimeter and shape is characterized by magnetic periodicities of Aharonov-Bohm type. Therefore, the loops host persistent currents at finite fields that, similarly to the B = 0 spectra, exhibit geometry dependent features in the small size limit, showing that currents would serve as an excellent tool to probe topologically confined states in BLG nanostructures.

A. A Quantized-perimeter model (QPM)
An analytic model of topological loops in BLG can be devised utilizing the results of Ref. 23 for a straight domain wall with an abrupt transition +V a to −V a (a kink). The analytic relations between linear momentum p along the kink and energy E of the two topological branches, corresponding to the two valleys (K, K ), were given in Ref. 23 as where the Fermi velocity v F = 660 meV nm and the interlayer coupling t = 380 meV are BLG intrinsic parameters.
In a closed loop of perimeter P we can assume that the tangent momentum p is quantized, such that an integer number of wavelengths must fit into P, where negatives n's represent negative momenta (i.e., opposite propagation) and the perimeter quantization in terms of wavelength λ reads P = |n|λ. We can also add a magnetic field, represented by a vector potential A, and write the circulation integral where A t is the component of the vector potential tangent to the kink. Noticing that P dl = P and P dl A t = Φ, Φ being the magnetic flux across the loop, we find an implicit condition from which one can derive the bound state energy E for a given loop perimeter, magnetic flux and principal quantum number n: where Φ 0 = h/e is the flux quantum. Equation (4) is our quantized-perimeter model (QPM) for bound states in BLG topological loops. As anticipated, when Φ = 0 the relation depends only on the perimeter and is totally independent of the loop geometry. At finite fields, however, the bound-state energies depend on the loop surface through the Aharonov-Bohm flux ratio Φ/Φ 0 . We expect the QPM to be reliable for large enough loops, when different parts of the loop do not interfere. In small loops, the traversal extension of the topological states becomes comparable to the distances inside the loop and the QPM breaks down. Below, we investigate deviations from the QPM using an exact approach for BLG within the continuum limit.

B. A quantum 2D model (Q2DM)
The low-energy Hamiltonian describing the states that are formed in two-dimensional (2D) BLG nanostructures for energies near the Dirac points reads 28,29 where σ x,y,z , τ x,y,z and λ x,y,z are Pauli matrices for the sublattice, valley and layer pseudo spins, respectively. A topological loop forms from a space dependent function V a (x, y) that, as sketched in Fig. 1a, takes the constant values +V a (−V a ) outside (inside) the loop. In our numerical simulations, we consider a smooth spatial transition with a diffusivity 30 s to mimic a realistic experiment 31,32 . Details of how we model the closed loops of different shapes are given in the Appendix. Finally, the field strength is included in the magnetic length l z = /eB. We numerically look for the eigenstates of the Hamiltonian given by Eq. (5) using finite difference discretization of the xy plane in a square grid 33 . The spurious solutions due to Fermion doubling [34][35][36] have been filtered out as in Ref. 30 by coarse graining the wave functions. Our method can handle any symmetry of the loop, in contrast to radial grids which are of smaller dimension and more efficient computationally, but can only describe radially symmetric structures by construction. In all the cases treated below we have checked that good convergence with the grid size is reached.
The Hamiltonian Eq. (5) fulfills a chiral symmetry given by operator C = σ x τ x λ y relating eigenstates of opposite energies CHC = −H with C 2 = 1. In absence of magnetic field it also fulfills time reversal symmetry Θ = iτ y K, with K representing complex conjugation, relating states of the same energy ΘHΘ = −H and with Θ 2 = −1. Figure 2 compares the bound state energies in the analytical (QPM) and numerical (Q2DM) approaches as a function of the loop perimeter when the structure is a circle. Qualitatively, both approaches nicely agree. At zero energy, the QPM predicts a sequence of uniformly spaced branch crossings, with spacing given by the zero-energy wavelength

III. RESULTS
which for V a = 10 meV yields λ 0 160 nm. We then find an excellent correspondence with the zero energy crossings obtained numerically in the Q2DM. Physically, each branch corresponds to an angular momentum and, therefore, can produce crossings with additional branches of different angular momenta. For energies departing from zero, the main difference between panels a and b of Fig.  2 is that the energy branches are more densely spaced in Q2DM than QPM. This difference can be attributed to the small diffusivity s = 12.5 nm used in the QD2M. A smoother 1D kink is known to reduce the spacing of energy branches 30,37 . In Fig. 3 we show the Q2DM energy spectra for different loop shapes as a function of their perimeters. Remarkably, the spectrum for noncircular loops show anticrossings of the branches. The square (b) shows an alternation of zero energy crossings and anticrossings as the perimeter increases. Both the rectangle (c) and the irregular polygon (d) show zero energy anticrossings in all cases, although the size of the anticrossing gap is not uniform. There is a general tendency for all shapes to decrease the anticrossing gaps for large perimeter, making the spectra shape independent in this limit.
The fact that, in addition to the circle, only the square shows crossings for some branches can be understood as a symmetry consequence. Inspired by the QPM, the branch crossing is possible in a square when, simultaneously, two conditions are fufilled P = n 1 λ 0 , e ik0Lx = ±i ⇒ P = (2n 2 + 1)λ 0 , n 1 , n 2 ∈ N , (7) where k 0 = 2π/λ 0 is the mode wavenumber. Only with n 1 odd the crossing is then possible, as indeed obtained in Fig. 3b. The second condition in Eq. (7) corresponds to the requirement that a translation along the perimeter by the length of a side changes the state in a ±i phase, which occurs when the total number of wavelengths to distribute in the four square sides is odd n 1 = 2n 2 + 1.
In a rectangle one should replace L x by L x + L y in this second condition, which then becomes not compatible with the first one and no crossing is therefore allowed. For the irregular shape (d), the lack of any symmetry is responsible for the absence of crossings in the spectrum. We discuss next the results in presence of a magnetic field B. Figure 4 compares the analytic and numeric results for a circular loop when B = 0.5 T. The magnetic field breaks the valley degeneracy of the loop eigenmodes, yielding asymmetric spectra with respect to energy inversion. We remark that Fig. 4 shows the results for a given valley. The results for the opposite valley, not shown, are shifted in the opposite direction restoring the symmetry of the spectrum for both energy and valley inversion (C symmetry). There is again a good qualitative agreement between QPM and Q2DM, the most noticeable difference being the denser bunching of the energy branches within Q2DM, as in the B = 0 cases discussed above.
In Fig. 5 we show the results of the 2D model at finite magnetic field for the different shapes. The field-induded asymmetry is similar when the geometrical structure changes. Branch anticrossings are present for nonsymmetic shapes, although they are not centered at zero energy but shifted vertically. Similarly to B = 0 (cf. Fig. 3), we find that with an increase of the perimeter the shape-dependent features of the spectrum are washed out.
As a function of magnetic field, for a fixed loop dimension and shape, the loop spectrum shows a periodic behavior that is reminiscent of the Aharonov-Bohm effect. Figure 6 shows the result for a square in the QPM and Q2DM, both in good qualitative agreement. The en- ergy branches come in doublets, originating from the two kink branches f 1,2 given by Eq. (1). The results in Fig. 6 are in general asymmetric with respect to energy inversion, in the same way as for the results as a function of perimeter (Fig. 5); the energy inversion symmetry being fulfilled at B = 0 only. For different shapes the results show a similar periodicity, only that the field spacings between branch doublets must be scaled by the enclosed flux Φ.
As with previous figures, the results of Fig. 6 are for a single valley; the energy branches for the reversed valley having the exactly opposite slope. The magnetic field slope of the branches is clearly indicating that the loop sustains persistent currents in finite fields, as discussed already for trivial and topological circular rings in Refs. 27 and 38, respectively. Hence, we investigate the shape dependence of the persistent current in topological loops as a possible probe of the geometric effects on topologically confined states. For definiteness, we focus on the current associated to the lowest positive-energy quasiparticle of the spectrum E 1 from the derivative The persistent current J 1 for the different shapes is shown in Fig. 7 for a fixed perimeter of 500 nm. Due to the aforementioned valley and energy inversion symmetry, the persistent current has an odd symmetry with respect to field inversion, implying a vanishing value at zero field. The current is characterized by a sequence of alternating plateaus of almost constant values. This is an effect of the nearly linear dispersion of the energy levels with magnetic field. Positive and negative plateaus correspond to persistent currents of opposite valleys, circulating with opposed chiralities. Along the sequence from high symmetry (a, circle) to low symmetry (c rectangle and d irregular shapes) the behavior looks qualitatively similar, but the plateau lengths are reduced due to the increased number of field transitions. This is also expected since the flux Φ is reduced for the deformed loops. Therefore, measuring the persistent current yields valuable information on the geometrical structure of the loop.

IV. CONCLUSIONS
We have investigated the confined states in topological loops built in bilayer graphene. We have devised an analytical model, the quantized perimeter model (QPM), based on the infinite straight kink, that provides the wavelength quantization along the loop perimeter. In addition, a full quantum 2D model (Q2DM), valid at low energies, has been numerically solved in order to ascertain the validity of the analytical model. A general good agreement between both models is found. For large sizes, the energy spectra are almost insensitive to the loop shape, as expected from QPM. For small sizes, we have found that Q2DM reflects shape dependence in the emergence of zero-energy crossings for circles, alternating crossings-anticrossings for squares, and only anticrossings for other more irregular structures. The magnetic field introduces energy-inversion asymmetries of the spectrum for a single valley and Aharonov-Bohm periodicities. Shape-dependent anticrossings in small loops are also present with magnetic fields. Finally, we have found that persistent currents are sensitive to shape-dependent features in presence of magnetic fields, which makes them a useful tool to look into bound states in topological sys-tems.
Overall, our results are helpful for tailoring potential structures, finding the optimal geometry for possible applications in quantum computation with topological qubits or for valleytronic devices. Since there could be deviations from the intended symmetry due to fabrication limitations, our study emphasizes the importance of taking into account emerging spatial asymmetries for a careful characterization of the underlying physical system.

Appendix A: Modeling
The Q2DM low energy Hamiltonian, Eq. (5), requires a specific spatial dependence of the asymmetric potential V a (x, y). This appendix contains the details of the modeling of different loop shapes with V a (x, y) smooth profiles parameterized with a diffusivity s. The circle is described by a simple combination of logistic functions in the radial coordinate r = x 2 + y 2 , To describe smooth square and rectangle shapes we first introduce the ranges of the Cartessian coordinates x ∈ [x a , x b ], y ∈ [y a , y b ] representing the loop inner area and then use a similar smooth parametrization in each Cartessian direction. Namely, where F 2D (x, y) = ( F(x; x b , s) − F(x; x a , s) ) × ( F(y; y b , s) − F(y; y a , s) ) .
Finally, for the irregular shape we use the same modeling of Eq. (A3) but with the modification that y b is no longer constant but the following piecewise function of x Notice that we assume centered positions of the loop, with x a,b = ∓L x /2 and y a,b = ∓L y /2. For completeness, the perimeters for the different shapes are P cir = 2πR 0 , P squ = 4L x , P rect = 2(L x + L y ), P irr = 3 2 (L x + L y ) + 1 2 L 2 x + L 2 y .