Controlling the layer localization of gapless states in bilayer graphene with a gate voltage

Experiments in gated bilayer graphene with stacking domain walls present topological gapless states protected by no-valley mixing. Here we research these states under gate voltages using atomistic models, which allow us to elucidate their origin. We find that the gate potential controls the layer localization of the two states, which switches non-trivially between layers depending on the applied gate voltage magnitude. We also show how these bilayer gapless states arise from bands of single-layer graphene by analyzing the formation of carbon bonds between layers. Based on this analysis we provide a model Hamiltonian with analytical solutions, which explains the layer localization as a function of the ratio between the applied potential and interlayer hopping. Our results open a route for the manipulation of gapless states in electronic devices, analogous to the proposed writing and reading memories in topological insulators.


Introduction
Two-dimensional Dirac materials, like graphene, have attracted remarkable interest for novel nanoelectronic applications due to their reduced dimensionality and extraordinary transport properties [1][2][3][4][5][6]. Recently, great effort has been devoted to seek and develop methods to control the transport of different degrees of freedom in these new materials. Spin [7], valley Gapless states in bilayer graphene have also been studied theoretically in works focusing on their topological features [27][28][29]. As long as the no-valley mixing condition is fulfilled, one can specify a topological quantity, namely, the valley Chern number [27,29]. Its change across the DWs accounts for the number of gapless states that appear under a gate potential. However, experiments and applications based on these conducting states need a clarification of other important physical features, such as spatial localization, where using topology arguments would be complex. A microscopic model is better suited to identify the symmetry and localization of the states. A few relevant works on lattice models specifically study different stacking boundaries for a gate voltage close to the interlayer coupling [27,28]. However, a systematic analysis of the layer localization in these stacking DW states under varying gate voltages is still missing.
In this work we show that the layer localization of gapless states in a gated BLG domain wall is tuned by an externally applied voltage. We demonstrate that at the domain wall the carriers are concentrated in the upper or bottom layer depending on the ratio between gate voltage and interlayer coupling. This dependence allows for switching the localization of topological states between layers by a change of the gate magnitude, implying that an additional degree of freedom, the layer, can play a role in BLG-based devices. Next, to explain the aforementioned behavior, we consider the symmetries of the bands and bonds formed between atoms in the gated system. We model a periodic array of stacking domain walls that, in contrast to the single DW, allows us to work with energy bands and wave functions. In this way, we identify the specific bands of the uncoupled layers which create the gapless states under interlayer coupling.
Controlling the layer localization by an external voltage opens a novel route for the manipulation of gapless states in electronic devices that could be denominated layertronics, in analogy to valley-and spintronics. We propose that layer localization would be another tunable degree of freedom in BLG, in addition to valley and spin.

Gated stacking domain wall: preliminary remarks
We first review the main features of the system. The left panel of figure 1 shows a single domain wall between AB and BA regions in bilayer graphene joined in the zigzag direction. Strained or corrugated graphene presents larger and more realistic domain walls, but for our theoretical analysis we choose an abrupt boundary that allows for an AB/BA stacking change. Previous works [23,24,28,30] have shown the robustness of topological states for smoother boundaries, with their main features preserved. For such a single boundary under a gate there are two topological gapless modes around the K point, as shown in the right panel of figure 1. The plot shows the local density of states (LDOS) resolved in energy E and wave vector k. These calculations employ a Green's function matching method 7 and a p z tightbinding model, with an intralayer hopping parameter γ 0 = −2.7 eV and a single interlayer hopping γ 1 = 0.1γ 0 [12,13]. There is another valley with another couple of gapless states with negative wave vectors with respect to figure 1 due to time-reversal symmetry. Note that the valley separation has motivated the proposal to employ such topologically protected modes for graphene valleytronics [31,32].

Layer localization with gate voltage magnitude
We report on the effect of layer localization exchange of topological states at a single stacking domain wall for small and large gate voltages. We use the same Green's function matching method and graphene model when there is an applied V. The LDOS around the K valley projected in the boundary nodes is shown in figure 2 for different voltages. The localization in the top and bottom layers is presented by the color scale from blue to red 8 . We start with the LDOS for ungated bilayer in panel (a). The LDOS localizes differently at the top and bottom layers because the symmetry between them is broken by the stacking domain wall. For a small gate voltage applied to the bottom layer, e.g. V = 0.1 eV in panel (b), the gap states appear, as it is well known, but they turn out to be separated in the two layers. The state on the left side of the cone is more localized at the bottom layer, while the state on the right is at the top layer. This localization can be explained as a perturbation of the LDOS of the ungated bilayer for small voltages applied. For increasing voltages, around V = γ 1 in panel (c), the states become fully mixed between layers. Next, for a large voltage, e.g. V = 0.5 eV in panel (d), the states are again separated in the top and bottom layer. This time, however, the state on the left side of the cone is more localized at the top layer, while the state on the right is at the bottom layer. Therefore, the topological states at the boundary atoms reverse localization for small and large voltages of the same sign.
Concerning the maximum experimental voltage between layers we have to take into account the maximum field applied by gate materials used in experiments. For instance, the breakdown field in standard SiO 2 substrates would allow a bias in bilayer graphene under dual gating up to 0.6-0.7 eV. These biases are larger than the 0.5 eV bias for which the layer switching is clearly seen, as discussed above in figure 3. As a limiting voltage and for the 7 The LDOS is calculated at the boundary in the zigzag direction. For details see [38,39,40,41]. 8 The LDOS ratio is plotted only when LDOS is greater than some set value. sake of completeness, we show in figure 3 the case for V = 1.2 eV bias, that allows us to comment on the asymmetry of the cones around specific K and K′ points.
One can also ponder the role of screening on layer localization under such voltages. It is noteworthy that DFT calculations [33], which take into account such effects, show that electron Coulomb screening between layers results in a smaller interlayer coupling. We have performed calculations for other smaller γ 1 values keeping either the γ 0 value or the γ 0 /γ 1 ratio fixed, as reported in DFT calculations. Including further electron screening by reducing the interlayer hopping γ 1 would decrease the required gate potentials for top/bottom layer polarization to smaller values. The layer polarization with gate value is preserved even when modifying the screening between layers to more realistic values (See suppl. inf. (stacks.iop.org/TDM/5/025006/mmedia), SI1, figure S1 for details). Additionally, note that the inclusion of an on-site Coulomb term in our model amounts to a higher localization of the electron. We have checked that it does not change our results for reasonable values of U around γ 0 [34,35].

Spatial distribution of TPS in each layer
Now we explore how the localization of topological states extends away from the boundary atoms. We investigate their spatial distributions at the Fermi level, as shown in figure 3. The LDOS of each state is averaged per rectangular zigzag unit cell (four atoms in each layer) and decomposed in the top and bottom layers. The spatial distributions show that localized states have their maxima close to the boundary region and decay in the adjacent unit cells [23]. The weight of the LDOS at the top and bottom layers is interchanged by increasing the voltage, around V ≈ γ 1 , as commented above.
It is also important to compare with the experiment to consider how large the region of layer interchange extends. Because of their gapless character, the states decay away from the domain wall in AB and BA stacking regions. Their decay rates increase with magnitude of the gate voltage. Note that the region where a state exists with a given layer polarization can be limited by its oscillatory character or fast decay (see cases V = 0.1 eV and V = 1.2 eV in figure 3). In an abrupt domain wall, the width of the top/bottom layer polarization extends at least for 6 u.c., i.e. about 4 nm. In the case of a smooth domain wall, such as induced by a tensile, shear or corrugated layer, we can expect two effects. On one hand, the finite transition region may affect the exact spatial distribution of the states and the gate voltage needed for the polarization change may slightly decrease. On the other hand, the interchange region will be extended with respect to the case for a sharp domain wall. Because the transition region acts as a gapless medium [23,28,36], TPS will propagate therein as extended states preserving their high or low localization in each layer [23], and showing the top/bottom layer localization even better resolved.
The layer localization of the gapless states with gate voltage exists also in the case of a smooth domain walls. The gate voltage for which we observe the change of the top/bottom polarization can slightly decrease, e.g. for the case of corrugated domain walls. Note that when the domain wall region is smooth, the interchange region increases. The major finding on top/bottom layer localization of TPS with gate value remains the same when the domain walls are smooth.
Topological states in gated bilayer graphene run across the insulator gap, so that their behavior can also be investigated by varying the chemical potential. Note that up to now we have looked at the topological states at the Fermi energy. However, experiments are usually performed away from the neutrality point, due to doping or to the interaction with different substrates. Figure 4 shows the LDOS distributions of topological states above and below the Fermi level. The left and right topological states are no longer complementary. We find that the total charge is more localized on one of the layers, i.e. layer localization is tuned by doping. This effect seems similar to the above case with high voltages, but with significant differences because: (i) the topological states can be on either top or bottom layer, and (ii) the effect is available at lower gate voltage values. Additionally, the distribution of the topological states in the top or bottom layers can be controlled by the experimental gate voltage applied in doped samples, a finding that has to be taken into account for practical applications.
As far as we know, this effect based on the asymmetry in localization of topological states has not been noticed and not explored before. We believe that the reason for this omission is that for V ≈ γ 1 , as normally set in most calculations, the energy gap saturates at that value and the topological states are equally disrtibuted between the two layers [23,28,37]. The layer distribution for these potential values is fully mixed at the boundary atoms.

Periodic stacking domain walls
To gain physical insight into the remarkable variation of the spatial distribution of these modes, we need to examine the symmetry of the corresponding wave functions. To this end, we study periodic systems, i.e. bilayer superlattices of stacking domain walls along the zigzag direction [38,39], which we label BS-DW 9 . The unit cell and band structures near the K point for different superlattice lengths are collected in figure S2 in supplementary information. It is important to note that for V > 0 and γ 1 = 0 the bands of the constituent graphene layers overlap. For nonzero γ 1 the overlapping bands interact and split yielding the energy gap in the case of pristine bilayer. However, when a stacking domain wall is imposed, two bands still persist in the gap. The treatment of γ 1 as a perturbation allows recognizing the bands of single graphene layers that give rise to the presence of topological bands in the energy gap, to be analyzed below.

Gaps and band crossing points near the Fermi level
We consider two uncoupled graphene layers with a gate potential V applied to the bottom layer and we switch on the interlayer hopping to study its role in the appearance of the topologically protected bands. Without hopping, the energy structure of two pristine graphene layers with a voltage difference applied is gapless, as shown in figure 5(a). We have chosen an 8-atom unit cell to compare more easily to the band structures of stacking DW superlattices. The bands of the gated bottom layer are shifted in energy by an amount V 9 The choice of periodic systems allows us to avoid the interplay between topological states and edge states characteristic for zigzag graphene edge [42]. with respect to the top layer bands. At k = 2/3π the pair that belongs to the top ungated layer crosses at E = 0, while the pair of bands of the bottom gated layer cross at E = V. Note that the bands originating from different layers cross at E = 1 2 V . One pair crosses for k ≲ 2/3π and another pair for k ≳ 2/3π. The bands are labeled B, T, indicating that they belong to the bottom and top layers, respectively, and s, a, due to the symmetric or antisymmetric character of the corresponding wavefunctions.
When the interlayer interaction is switched on in the BLG case ( figure 5(b)), the crossing bands mix and split because they form pairs of bonding and antibonding states and a gap opens. The resulting band structure near the Fermi level has the well-known Mexican hat shape (see the bands plotted in figure 5(b) below the unit cell) [12]. However, when the nodes are connected producing the stacking defects, as in figure 5(c), a pair of states remains in the gap, shown in the band structure depicted under the corresponding unit cell. The resulting structure is a BS-DW with W = 1. A key question that we address below is why these states are gapless and how they arise from the bands of pristine graphene.
We focus on the left crossing point marked with a circle in figure 5(a). The two bands crossing therein are labeled Ts and Ba, due to their symmetry and localization. The bonding and antibonding p z orbital combinations formed due to the interlayer interaction are schematically illustrated in figure 6(a). However, for a domain wall, there are always topological modes crossing the gap. Figure 6(b) graphically shows that it happens because it is not possible to form bonding and antibonding combinations of the top and bottom layer wave functions simultaneously for both pairs of the overlapping nodes from different layers in a DW.

Basis functions for the crossing bands
We next treat the interlayer hopping as a perturbation to the bands of the uncoupled layers. The starting basis is given in terms of the uncoupled bands depicted in figure 5(a). For the 8atom unit cell employed therein, and explicitly labeled in figure 5(b), we have Ψ Ba, k = 1 2 ϕ α + e ikd ϕ β − e ikd ϕ γ − ϕ δ Ψ Bs, k = 1 2 ϕ α − e ikd ϕ β − e ikd ϕ γ + ϕ δ Ψ T a, k = 1 2 e ikd ϕ α′ + ϕ β′ − ϕ γ′ − e ikd ϕ δ′ Ψ T s, k = 1 2 e ikd ϕ α′ − ϕ β′ − ϕ γ′ + e ikd ϕ δ′ (1) where the labels refer to their layer localization and their symmetry, as in figure 5(a); d is the distance between contiguous rows of atoms in the zigzag direction (see figure 6); ϕ μ denotes the p z orbital in the μ atom, with μ running from α to δ in the bottom layer and the same labels with primes in the top layer (see figures 5(b) and (c)). The wave vector k dependence is explicitly indicated as a subscript, but we omit it from now on for the sake of simplicity.
For connected layers, the degeneracy at the crossing points of the uncoupled system k As commented before, in a BS-DW one pair of the bands crossing the gap has positive velocity and the other pair has a negative slope. In the smallest BS-DW with W = 1, the states at the boundaries are unavoidably mixed. For large W, each pair belongs to a different stacking domain wall in the unit cell, being spatially separated, as illustrated in figure 8. The bands for W = 40 are resolved in the two stacking boundaries, DW and DW. When one stacking DW is present in the system, one pair of bands is in the gap with the same velocity, as shown in figure 1. The sign of the velocity is related to the change of stacking, either AB to BA or vice versa. Further, the different localization between layers predicted from the periodic calculations are compared with the case for a single DW, as shown in section 3. In figure 9 we present the wavefunctions corresponding to E f ± δ around the DW. We find that a DW either isolated or in a periodic arrangement behaves the same way.

Continuum Hamiltonian for topological states
We generalize the above discussion to k values away from the band crossing points. The Hamiltonian H of each layer close to the Dirac point is represented in the basis (1), given by Ψ Ba and Ψ Bs for one layer, and Ψ Ta and Ψ Ts for the other layer. For a monolayer, H is a 2 × 2 diagonal matrix with linear k, and −k values in the diagonal elements. For gated bilayer graphene we have to double this matrix [2,37]. Using the basis{Ψ Ba , Ψ Bs , Ψ Ta , Ψ Ts } we obtain the interaction Hamiltonian between layers that introduces the hopping γ 1 in the nondiagonal elements. Employing the expression of these basis vectors in terms of the localized orbitals ϕ μ given in equation (1) and the geometry of the interlayer bonds shown in figure   5(b), the Hamiltonian for bilayer graphene is written where α = 3aγ 0 /2ℏ, and a is the graphene lattice constant. For V = 0 the eigenvalues of this Hamiltonian have a quadratic dependence on k, as it should be in bilayer graphene [37]. Two of them are degenerate at the Dirac point. For V ≠ 0 the gap opens and the low energy bands versus k show the Mexican hat shape. In the gated system the energy gap E g depends on the parameters V and γ 1 . For V < γ 1 , E g ≈ V, while for γ 1 < V, E g ≈ γ 1 . The Hamiltonian (2) gives an appropriate description of the low-energy bands of BLG.
In the same way, i.e. using the basis vectors given in (1), we derive a Hamiltonian for bilayer graphene with periodic stacking domain walls as, The Hamiltonian (3) includes a k dependence, but it does not mix k and −k values inside or between the layers. The eigenvalues remain linear in k, so for V ≠ 0 they cross the energy gap constituting the topologically protected gap states. Its eigenvalues have the following analytical expression: The corresponding eigenvector components for velocity +k are 0, −V ± γ 1 2 + V 2 γ 1 , 0, 1 (5) and for velocity −k V ± γ 1 2 + V 2 γ 1 , 0, 1, 0 . (6) Note that the localization of the states will behave differently in the two gate voltage regimes. In the limit of V ≪ γ 1 the eigenvectors are mixed between layers-in contrast to the single domain wall case, where we observe localization reversal at the boundary. The reason for layer separation was the asymmetry between layers, which is no longer present in BS-DW.
In the limit of V ≫ γ 1 the components of the eigenvectors given by equation (5) approach (0, 0, 0, 1) and those given by equation (6) approach (0, V/γ 1 , 0, 1), and to (0, 0, 1, 0) and (V/γ 1 , 0, 1, 0). This means that each pair of the gap states with the same velocity (+k or −k) has the corresponding wave functions localized in different layers. The wave functions in each pair have the same symmetry, i.e. a or s (see equation (1)). Therefore, the topologically protected gap states look like the single layer bands that cross at the Fermi level. This continuum Hamiltonian with four bands allows us to predict the relevant properties of topological states in stacking domain walls.

Discussion and conclusions
We have investigated the gapless states with topological character that appear in gated bilayer graphene with stacking domain walls. By employing atomistic models, we find that each of the two topological states in a valley is layer-resolved; furthermore, their localization is switched between the top and bottom layer by varying the magnitude, but not the sign, of the gate voltage. Therefore, besides the valley and sublattice degrees of freedom, these states can also be labeled by a layer index.
Our findings also bring forth an atomistic understanding of the origin of gapless states with topological character. We have analyzed the chemical bond formation between gated layers of graphene with domain walls to elucidate the appearance of gapless states. We have also provided a continuum model for the gapless states that correctly describes the swithing in layer localization of the topologically protected states. For a large gate voltage, i.e. above the interlayer coupling, the layer localization presents the standard trend between top and bottom layer for a particular voltage orientation. However, for a small gate voltage, the asymmetry between layers introduced by the stacking domain walls prevails, so the gapless states are slightly perturbed by this voltage, having an opposite layer localization. Any TPS in other materials could show similar layer localization control with the gate voltage, at least if there are two layers and two states involved, although the asymmetry induced by domain walls seems also required. The focus here has been on graphene-based materials. The hope is that this will stimulate interest in other materials that could also display layertronics.

Supplementary Material
Refer to Web version on PubMed Central for supplementary material.    Unit-cell averaged LDOS at the bottom (red) and the top layer (blue) plotted for both topologically protected states for V = 0.5 eV above and below the Fermi level; δ = 0.1 eV.
Schemes of the two topological states follow the notation given in figure 3. The layer localization is now distinct at voltages accessible in experiments.