Equilibrium currents in a Corbino graphene ring

We address the description of a graphene Corbino disk in the context of a tight binding approach that includes both kinetic and Rashba spin-orbit coupling due to an external out-of-plane electric field. Persistent equilibrium currents are induced by an external magnetic field breaking time reversal symmetry. By direct diagonalization, we compute the spectrum and focus on the dispersion near the $K$ points at the Fermi level. The dispersion keenly reproduces that of a continuum model in spite of the complexity of the boundary conditions. We validate the assumptions of the continuum model in terms of predominant zig-zag boundaries conditions and weak sub-band coupling. The wave functions displaying the lowest transverse modes are obtained, showing the predominance of edge states with charge density at the zig-zag edges. The persistent charge currents, nevertheless, do not follow the traditional argument of current cancellation from levels below the Fermi level, and thus they depart in the tight-binding from those found in the continuum model.


Introduction
Graphene has been a source of recent interest because it is a stable two dimensional material, with linear dispersion around the Fermi level, and chiral edges states exhibiting a variety of exciting and tuneable properties [1]. Another reason for the interest is that in the vicinity of the Dirac cones, many exact results can be obtained and the tight-binding (TB) approach is perfectly suited, differing very little from first principles calculations.
In the long-wavelength limit, expansion of the kinetic energy around the Dirac points leads to a set of two Dirac Hamiltonians where τ z is the Pauli matrix acting on the valley index which labels the Dirac points K = 4π 3a (1, 0) and K = −K and σ i are the Pauli matrices which encode for the two sublattices α = A, B . When the Rashba spinorbit (RSO) interaction, due to an external out-of-plane electric field, is taken into account, the additional term coupled to the electron spin s appears [2], namely H = H 0 ⊗ 1l s + λ R (τ z ⊗ σ x ⊗ s y − 1l τ ⊗ σ y ⊗ s x ). (2) Written in a Corbino ring geometry of average radius R which selects only the lowest transverse modes along the ring, the energy levels can be calculated as [3], E κ,δ m,λ = κ 2 (ħv F /R) 2 (1 + 4m 2 ) + 8λ 2 R − 4δ (mħv F /R) 2 + λ 2 R (ħv F /R) 2 + 4λ 2 R . ( The parameters κ = ±1 and δ = ±1 refer to the particle-hole and the SO branch indices, respectively. In order to generate persistent charge currents, time reversal symmetry should be broken, e.g. by a magnetic field flux Φ piercing the ring. The energy levels are modified according to the substitution m → m +Φ/Φ 0 (with Φ 0 = 2πħ/|e| the flux quantum) which follows from the ordinary U (1) minimal coupling gauge substitution. The equilibrium charge currents then follow from the relation Such a continuum model is particularly simple and the wave functions satisfying the topology of the ring can be easily obtained [3]. Nevertheless, there are a few important assumptions regarding the boundary conditions of the Corbino ring that must be appropriately addressed: i) The internal and external edges of the ring are neither zig-zag nor armchair, but continuously change as a function of the radial angle.
ii) The transverse modes can be coupled to the longitudinal ones even for the purely kinetic term, as is the case for armchair ribbons and iii) the spin-orbit coupling is also a source of coupling to transverse modes, so that a single ground state mode theory could be compromised.
In this work we study the validity of the continuum approximation by treating the Corbino geometry through a tight-binding Hamiltonian approach, thereby correctly incorporating all boundary conditions and verifying the continuum model assumptions. The paper is organised as follows: We first discuss the tight-binding Hamiltonian contemplating both the kinetic term and the Rashba coupling. Due to the geometry of the Corbino ring, it is not straightforward to address the azimuthal momentum within the tight-binding picture. Thus, we make use of a mapping of azimuthal momentum to the corresponding gauge field, and we can derive the dispersion as its function. The closing of the wave functions, a condition that must be satisfied for the delocalised electron case, can then be used to compute the full charge currents.

The tight-binding Hamiltonian
Through the Tight-Binding approach, one can obtain the dispersion in the whole Brillouin zone (not only the vicinity of the Dirac points), and keep track of the hexagonal structure of the lattice. We consider the model Hamiltonian consisting of two terms. First the kinetic energy, here written as a hopping term between nearest neighbours (nn) lattice sites on the lattice 〈i , j 〉, comprises a phase that encodes the magnetic flux. To evaluate this Aharonov-Bohm phase, in the basis where Bravais lattice vectors are given by a 1 = ax and a 2 = (1/2)ax + ( 3/2)aŷ and with the Landau gauge A = −B yx, one computes, for each link on the lattice, the contribution i e ħ na 1 +ma 2 +δ with r(η) = na 1 + ma 2 + ηδ, δ being one of the n, m site's neighbors.
The second term corresponds to the Rashba spin-orbit interaction due to an external gate voltage producing an electric field perpendicular to the graphene sheet and to the strength of the intrinsic spinorbit interaction [4][5][6]. The leading contribution results from nearest neighbour terms.
This Hamiltonian can be obtained as follows. First, in the absence of any spin-orbit interactions, the Hamiltonian only contains kinetic energy, where the sum is over the lattice sites x of the Bravais lattice and, in general, different orbital types α (α = A, B in the case of graphene). The hopping parameter −t α,α x,x = 〈x, α|H 0 |x , α 〉 is assumed homogeneous (= −t 0 ), and creation and annihilation operators can be defined c † x,α c x ,α = |x, α〉〈x , α | in the 33803-2 localized states basis. We usually take into account only nn sites x and x , like in equation (5). The dominant contribution to the electron transfer between sites come from the (ppπ), contributions from the unbonded |p z 〉 orbitals perpendicular to the graphene plane, so the transfer is between orbitals of the same type, and thus α = α in equation (6).
The Rashba Spin-Orbit interaction operating in graphene depends on two microscopic contributions; the coupling of atomic |p z 〉 to |s〉 both through the external electric field and the intrinsic SO coupling.
The coupling is thus first order in external electric field and first order in the internal atomic electric field sampled by transport electrons. Thus, both the external electric field and the atomic SO couple the atomic |s〉 and |p z 〉, and being combined, generate the lowest order coupling [4].
The intrinsic SO coupling generates a spin flipping mechanism that can be understood from the form of the bare Pauli Hamiltonian, where the electric field involved is that coming from the intrinsic atomic field. In a material, the denominator with the bare mass of the electron should be replaced by the gap in the band structure. This term may be encoded in the kinetic energy term of the (non-relativistic) electrons through the minimal coupling to a non-Abelian gauge field W a [7][8][9][10][11] with g = ħ, the coupling constant to W a . To this gauge field corresponds an unitary operator For a space-independent gauge field, expanding the exponentials to first order, we get Keeping only the leading contribution which comes from nearest-neighbouring lattice sites, the quantity (E × s) · (x − x) reduces to E (s ×d nn ) z a/ 3, where a/ 3 = 1.42 Å is the C−C distance andd nn is the unit vector between nearest neighbours x and x . The second term in equation (5) follows, with the Rashba parameter λ R = eE a(−t 0 ) 4 3mc 2 . This value of the SO interaction has been derived using the coupling in the vacuum when in fact it should be that of electrons within the carbon atoms of graphene. Different authors have derived this form from a Slater-Koster approach [4][5][6] and second order perturbation theory, resulting in the form λ R = eE z 0 ξ/3(spσ), where z 0 is proportional to the atomic size, ξ is the atomic SO strength, (spσ) is the Slater-Koster integral for the s − p matrix element along the σ C−C bond. This is the way atomic SO strengths translate into conduction electron SO strengths.

TB results on the Corbino disk
A Corbino disk is built from the graphene plane, keeping lattice sites which sit between two concentric circles of given radii (figure 2). All dangling bonds (sites of coordination one) are eliminated and no account is made for edge deformations [12] due to the reduction in coordination or surface ondulation. Note that there is a dominance of zig-zag edges versus armchair which is general for disks of any size [13,14].
The TB Hamiltonian of equation (5) is now restricted to the lattice sites forming the disk. The Hamiltonian, restricted to these sites is diagonalised and eigenvalues and eigenfunctions are derived using Kwant software [15].

Dispersion relation
When a magnetic flux is applied across the ring, the magnetic flux being varied may serve as a probe to scan the Brillouin zone (BZ), as can be seen from the minimal coupling prescription p → p − eA. The changes in p are symmetric to the changes in A. In the Corbino geometry we have zero radial momentum p r as there is a fixed transverse state, so the momentum associated to currents will be purely azimuthal p ϕ . We can then use the corresponding component of the vector potential A ϕ to sweep the azimuthal wave vectors up to a constant. The wave vectors in the TB approach are measured with respect to the center of the BZ (Γ point), but Φ (A ϕ ) being changed enables one to approach the Dirac points (K and K points) by changing the crystal momentum analog A ϕ i.e. the magnetic field [16,17]. At the K points, the valence band and the conduction band cross the Fermi level with linear dispersion in the absence of SO interactions. The dispersion close to the Fermi energy, obtained by following the vector potential, is depicted in figure 3. We have also depicted two closest upper and lower sub-bands to assess the energy gap in order to excite those sub-bands. One can see that this gap is five times the energy scale associated with the width of the energy structure close to the Fermi level. This energy separation justifies the assumptions made in [3] In order to incorporate the quantized values of the energies due to the closure conditions on the wave function (whose trace we follow with the flux) we recognize that different values for these quantum numbers follow along this same structure simultaneously, as the field is changed. To reproduce this evolution, one needs to superpose dephased traces on the same pattern [3]. This procedure is illustrated in figure 4, where we depict the raw dispersion, on the left, and the dephased repeated traces on the middle panel, identifying the corresponding quantum numbers. From the continuum model, the value for the separation between quantised values is ∆Φ/Φ 0 = 2 2 + 4λ 2 R ), where = ħv F /R.

Charge persistent currents
Once the energy dependence on the flux is established, we have evidence of delocalised states and thus the possibility of persistent currents. By referring to figure 4, we use equation (4) in order to compute the charge current. While the continuum treatment only contemplated the spectrum depth depicted in the figure (almost identical to the continuum description), here we have a detailed access to all occupied levels, built from the 2s and 2p orbitals. Figure 5 shows the persistent currents derived from the spectrum in the right-hand panel of figure 4. The continuum model by definition has infinite bands and the assumption is that the deep level bands produce current cancelling contributions. The finite Corbino disk treated here, on the other hand, has a finite spectrum, and the deep level responds very strongly to the magnetic field. The sum of all contributions then has more a signature of the deep levels than that of the levels close to the Fermi surface. Whether this makes sense for finite systems might be a very interesting possibility, since that would mean persistent currents protected from thermal fluctuations [3].

Charge density
In the continuum treatment [3], no attention was paid to the details of the lateral confining potential, since the universal features of the spectrum and the persistent currents, do not depend on it. Nevetheless, the persistent currents are carried by the edge states, and their configuration is accessible to the TB solution. Figure 6 shows the charge density on the disk for the ground state, when solely the kinetic term is considered at the K point (left-hand panel). For the top of the valence band states, the charge piles up on the edges of the disk, almost exclusively on the outer zigzag edges [18], with no expected spin distinction. The non-bonding character of the charge distribution seen for zig-zag edges in nano-ribbons [19], does not survive in the Corbino geometry, even in the bulk. This distribution of charge describes the nature of the ground state transverse mode of the disk. Once the Rashba SO term is turned on, for sufficiently large values of the coupling, we see a noticeable shift of the charge density on the edges (figure 6 right-hand panel). While an overwhelming amount of charge is piled up on the zig-zag edges, there is a modulation in the charge density as compared to the free electron case. Thus, there appears to be a way to tune the charge distribution on the edges by changing the electric field perpendicular to the graphene plane. This concentration of charge is more an indication that additional terms in the Hamiltonian must be contemplated, such as electron-electron interactions and electron-phonon interactions, which have been associated to edge magnetism [19] and lattices distortions, respectively.

Conclusions
We have assessed a tight-binding version of the Corbino disk for graphene in the presence of Rashba SO interactions and discussed in the light of analytical solutions for the continuum approximation. The assumptions regarding boundary conditions based on the works in references [13,14] and [19], where these are predominantly zig-zag and the fact that a single mode approximation is appropriate, are justified. The full wave functions for the Corbino disk were drawn for the free and SO coupled cases, where the transverse mode structure revealed both internal and external radius edge modes.
Regarding the charge current on the disk, there are, nevertheless, qualitative differences arising from the finite spectrum in the TB approach (while the continuum has an infinite number of states). Current states from a finite set of levels being added up, give rise to contributions from deep in the Fermi sea that change the currents qualitatively and allow for the possibility of enhanced robustness of persistent current against thermal fluctuations.