Spin-orbit interaction and snake states in graphene $p$-$n$ junctions

TL;DR

This paper investigates spin-orbit interactions and snake states in graphene p-n junctions under magnetic fields, revealing zero modes with spin polarization and interference effects useful for spintronic devices.

Contribution

It provides an exact spectrum analysis of spin-resolved Landau levels and identifies wave-vector-dependent zero modes as quantum snake states influenced by Rashba spin-orbit interaction.

Findings

Discovery of linearly dispersing zero modes with spin polarization

Rashba spin-orbit interaction causes observable wave vector shifts

Potential application in spin field-effect transistors

Abstract

We study a model of a $p$-$n$ junction in single-layer graphene in the presence of a perpendicular magnetic field and spin-orbit interactions. By solving the relevant quantum-mechanical problem for a potential step, we determine the exact spectrum of spin-resolved dispersive Landau levels. Close to zero energy, we find a pair of linearly dispersing zero modes, which possess a wave-vector-dependent spin polarization and can be regarded as quantum analogous of spinful snake states. We show that the Rashba spin-orbit interaction, in particular, produces a wave vector shift between the dispersions of these modes with observable interference effects. These effects can in principle provide a way to detect the presence of Rashba spin-orbit interaction and measure its strength. Our results suggest that a graphene $p$-$n$ junction in the presence of strong spin-orbit interaction could be used as a building block in a spin field-effect transistor.

Equations

$${\mathcal{H}}_{\tau }$$

$${\mathcal{H}}_{\text{Z}}=\frac{g_s{\mu }_B}{2}B{s}_z.$$

$${\mathcal{H}}_{\text{RSO}}$$

$${\mathcal{H}}_{\text{ISO}}$$

$${\mathcal{H}}_{-}=(i{\sigma }_y){\mathcal{H}}_{+}(-i{\sigma }_y),$$

$$V(x)=V_0\text{sign}(x),V_0>0.$$

$${\Psi }_k(x,y)=\frac{1}{\sqrt{L}}{e}^{iky}{\psi }_k(x),$$

$${\psi }_k^T(x)=\left({\psi }_{\text{A}\uparrow },{\psi }_{\text{B}\uparrow },{\psi }_{\text{A}\downarrow },{\psi }_{\text{B}\downarrow }\right).$$

$$\left[\mathcal{H}(k)-E\right]\psi (x)=0,$$

$$\hat{a}={\partial }_q+\frac{q}{2},{\hat{a}}^{†}=-{\partial }_q+\frac{q}{2},[a,a^{†}]=1,$$

$$\mathcal{H}(k)-E=\left(\begin{array}{cccc}-{\mu }_{+}& -i\sqrt{2}\hat{a}& 0& 0\\ i\sqrt{2}{\hat{a}}^{†}& -{\mu }_{-}& -i\lambda & 0\\ 0& i\lambda & -{\nu }_{-}& -i\sqrt{2}\hat{a}\\ 0& 0& i\sqrt{2}{\hat{a}}^{†}& -{\nu }_{+}\end{array}\right),$$

$${\mu }_{\alpha }$$

$${\nu }_{\alpha }$$

$${\psi }^{\gtrless }(x)=\left(\begin{array}{c}{d}_1{D}_{p-2}(\pm q)\\ i{d}_2{D}_{p-1}(\pm q)\\ d_3{D}_{p-1}(\pm q)\\ i{d}_4{D}_p(\pm q)\end{array}\right),$$

$$\hat{a}{D}_p(\pm q)$$

$${\hat{a}}^{†}{D}_p(\pm q)$$

$$\left(\begin{array}{cccc}-{\mu }_{+}& \pm \sqrt{2}(p-1)& 0& 0\\ \pm \sqrt{2}& -{\mu }_{-}& -\lambda & 0\\ 0& -\lambda & -{\nu }_{-}& \pm \sqrt{2}p\\ 0& 0& \pm \sqrt{2}& -{\nu }_{+}\end{array}\right)\left(\begin{array}{c}{d}_1\\ d_2\\ d_3\\ d_4\end{array}\right)=\left(\begin{array}{c}0\\ 0\\ 0\\ 0\end{array}\right).$$

$$(2p-2-\mu )(2p-\nu )-{\lambda }^2{\mu }_{+}{\nu }_{+}=0,$$

$$\mu \equiv {\mu }_{+}{\mu }_{-}$$

$$\nu \equiv {\nu }_{+}{\nu }_{-}$$

$$p_{\pm }(E)=\frac{1}{4}\left[\mu +\nu +2\pm \sqrt{(\mu +2-\nu {)}^2+4{\lambda }^2{\mu }_{+}{\nu }_{+}}\right].$$

$${\psi }_p^{\gtrless }(x)$$

$$\psi =c_1{\psi }_{p_{+}}^{>}+c_2{\psi }_{p_{-}}^{>}+c_3{\psi }_{p_{+}}^{<}+c_4{\psi }_{p_{-}}^{<}.$$

$$E_{0,\pm }$$

$$E_{n,\pm }$$

$$E_{n,\pm }=\text{sign}(n){\left[2|n|-1+\frac{{\lambda }^2}{2}\pm \sqrt{{\left(1-\frac{{\lambda }^2}{2}\right)}^2+2|n|{\lambda }^2}\right]}^{\frac{1}{2}},$$

$$\psi (x)=\{\begin{array}{ll}{c}_1{\psi }_{p_{+}}^{<}(x)+c_2{\psi }_{p_{-}}^{<}(x)& x<0,\\ c_3{\psi }_{p_{+}}^{>}(x)+c_4{\psi }_{p_{-}}^{>}(x)& x>0,\end{array}$$

$$c_1{\psi }_{p_{+}}^{<}(0)+c_2{\psi }_{p_{-}}^{<}(0)=c_3{\psi }_{p_{+}}^{>}(0)+c_4{\psi }_{p_{-}}^{>}(0),$$

$$\mathbf{W}\mathbf{c}=0,$$

$$\mathbf{W}=\left[\begin{array}{cccc}{\psi }_{p_{+}}^{<}(0)& {\psi }_{p_{-}}^{<}(0)& -{\psi }_{p_{+}}^{>}(0)& -{\psi }_{p_{-}}^{>}(0)\end{array}\right],$$

$${\mathbf{c}}^T=\left(\begin{array}{cccc}{c}_1& c_2& c_3& c_4\end{array}\right).$$

$$\det \left[\mathbf{W}(k,E)\right]=0.$$

$$\Delta k\equiv k_{+}-k_{-}=\frac{\lambda }{\hslash v_{\text{F}}},$$

Full text

I Introduction

Fifteen years after the isolation of graphene Novoselov et al. (2004, 2005), fundamental research on this exceptional two-dimensional material is still very active. Spin-orbit interaction (SOI) effects, in particular, continue to attract much attention. Graphene is, in fact, a promising material for spintronic applications Žutić et al. (2004); Han et al. (2014); Bercioux and Lucignano (2015); Offidani et al. (2017), where control of the SOI plays a crucial role. Moreover, SOI is at the basis of fundamental phenomena as the quantum spin Hall effect, whose theoretical discovery in graphene Kane and Mele (2005) started the field of topological matter Hasan and Kane (2010); Chiu et al. (2016); Bercioux et al. (2018).

Due to the small atomic number of carbon, the SOI in pristine graphene is very weak (on the order of $10\,\mu$eV)  Gmitra et al. (2009); Konschuh et al. (2010). However, in the last few years there has been an increasing body of experimental and theoretical evidence that SOI can be artificially enhanced either by chemical functionalization or by placing graphene on strong-SOI substrates Castro Neto and Guinea (2009); Weeks et al. (2011); Hu et al. (2012); Marchenko et al. (2012); Balakrishnan et al. (2013, 2014); Otrokov et al. (2018). The latter approach is particularly promising, as it preserves the electronic properties of graphene. While the first experimental reports presented signatures of enhanced SOI only in the spectral properties, a breakthrough came with the advent of two-dimensional van der Waals materials Geim and Grigorieva (2013), in particular the transition metal dichalcogenides (TMDC) WSe2, WS2, and MoS2 Roldán et al. (2014). Graphene/TMDC heterostructures allow for high quality interfaces, which preserve the large mobility of graphene and provide ideal platforms for transport measurements, showing uniform enhancement of SOI of up to three orders of magnitude Avsar et al. (2014); Gmitra and Fabian (2015); Gmitra et al. (2016); Wang et al. (2015, 2016); Yang et al. (2017); Benítez et al. (2017); Cummings et al. (2017); Wakamura et al. (2018); Zihlmann et al. (2018); Wakamura et al. (2019); Garcia et al. (2018). For example, the Rashba SOI for graphene on WS2 is estimated to have a magnitude of up to $5$ meV Wang et al. (2015, 2016), and even larger for graphene on a metallic substrate, where it can reach a value of the order of $50$ meV Otrokov et al. (2018); Marchenko et al. (2012).

Motivated by these developments, we study here the effects of strong SOI on the electronic properties of graphene p-n junctions in a perpendicular magnetic field. From the early days of graphene research, these systems have attracted much attention, in part due to the peculiar transport channels that form at the interface between the $p$- and $n$-doped regions Williams et al. (2007); Abanin and Levitov (2007); Özyilmaz et al. (2007). These chirally propagating electronic states can be visualized as snake states, whose quasiclassical trajectories describe a snaking motion along the interface. They play an important role in the magnetotransport properties of these devices Williams and Marcus (2011); Barbier et al. (2012); Liu et al. (2015); Taychatanapat et al. (2015); Rickhaus et al. (2015); Cohnitz et al. (2016); Fräßdorf et al. (2016); Tóvári et al. (2016); Kolasiński et al. (2017); Makk et al. (2018). Similar snake states are expected also for inhomogeneous magnetic fields containing a line separating $B>0$ and $B<0$ regions, without a $p$-$n$ junction Park and Sim (2008); Ghosh et al. (2008); Oroszlány et al. (2008). These purely magnetic snake states might give rise to effects analogous to those presented in this paper, but we do not further consider them here. We note in passing that the effects of SOI on the transport properties of p-n junctions in the absence of magnetic field have been studied in Yamakage et al. (2009).

In this article, we investigate how the presence of finite SOI changes the interface states in a graphene p-n junction. We focus on the low-energy physics described by the Dirac-Weyl Hamiltonian, and include in our model two SOI terms 111Recently, a third type of SOI has been argued to exist in graphene/TMDC heterostructures, the so-called valley-Zeeman SOI Gmitra and Fabian (2015); Wang et al. (2015, 2016); Gmitra et al. (2016); Cummings et al. (2017); Zihlmann et al. (2018). This term has the form of a Zeeman coupling to an effective perpendicular magnetic field having opposite sign at the two valleys. Thus, within the single-valley description we use in this paper, one could in principle incorporate it in the Zeeman coupling.: the Rashba SOI Kane and Mele (2005); Rashba (2009), due to inversion symmetry breaking by the substrate or, possibly, by a perpendicular electric field; and the Kane and Mele intrinsic SOI Kane and Mele (2005). We present the exact solution of the corresponding quantum mechanical problem in the presence of a potential step, and thereby provide a full description of the energy spectrum and the eigenstates of the system. We find infinitely many branches of dispersive Landau levels (LLs), with flat portions separated by regions with approximately linear dispersion. The shape of the levels can be understood as follows. The electron states are parameterized by the wave vector along the junction, $k$, which defines the centre around which the wave function is localized within a magnetic length. At large $k$, the states are localized deep in the bulk, away from the junction, and their energy is essentially independent of $k$. These states correspond to graphene’s LLs modified by the SOI. At small $k$, the states are localized close to the interface and feel the potential step. Their energy acquires a $k$-dependence and, as a consequence, the corresponding states have a finite group velocity along the junction. Due to the presence of the magnetic field, the system is chiral and the group velocity has the same sign for all branches. We will focus in particular on the two chiral modes evolving from the pristine zero-energy LLs. These modes are responsible for the low-energy transport properties when the Fermi level is close to the charge neutrality point. We refer to them as the snake modes. We will show that their group velocity depends on the size of the potential step, it is proportional to the classical drift velocity in crossed electric and magnetic fields, and it is only weakly affected by the SOI.

In the absence of SOI and Zeeman coupling, the snake modes are degenerate. We show that the Rashba SOI removes this degeneracy and leads to a relative wave vector shift between their dispersions. This resembles the relative horizontal shift of the energy dispersions of spin-up and spin-down electrons which occurs in a single-channel quantum wire in the presence of Rashba SOI, and suggests that graphene p-n junctions could be exploited to realize the Datta-Das effect Datta and Das (1990); Bercioux and Lucignano (2015). Indeed, electrons injected at a given energy at some point along the junction with a definite spin polarization will populate the two modes and will acquire different phase shifts while propagating. Monitoring the current downstream, one would observe a modulation due to quantum interference effects between the two modes.

Before closing this section, let us briefly describe how this paper is organized. In Sec. II we formulate the model for a graphene p-n junction in a perpendicular magnetic field in the presence of finite SOIs, and review its exact solution for a uniform potential. In Sec. III we determine the exact spectrum of eigenstates of the quantum mechanical problem with the potential step. In Sec. IV we discuss the corresponding observables, namely group velocity and spin polarization. Finally, in Sec. V we present some conclusions and discuss possible extensions of our work. Some technical details are relegated to the appendix.

II Model and general solution

Our starting point is the Dirac-Weyl Hamiltonian for graphene in a perpendicular magnetic field ${\bf B}=B\hat{z}$ with SOI  Castro Neto et al. (2009); Rashba (2009); Bercioux and De Martino (2010); Liu et al. (2012):

[TABLE]

where $\bm{\Pi}=\textbf{p}+\frac{e}{c}\textbf{A}$, A and $V$ are respectively vector and scalar potential, and ${\mathcal{H}}_{\text{Z}}$ is the Zeeman term

[TABLE]

We include both the Rashba SOI and the Kane-Mele intrinsic SOI:

[TABLE]

The symbol ${\bm{\sigma}}=(\sigma_{x},\sigma_{y},\sigma_{z})$ (resp. ${\bf s}=(s_{x},s_{y},s_{z})$) denotes as usual the Pauli matrices in sublattice (resp. spin) space, while $\tau=\pm$ refers to the valleys $K/K^{\prime}$. In Eq. (1) we have assumed that the two valleys are decoupled, which holds provided all potentials are smooth on the scale of graphene’s lattice constant 222Our analysis does not apply to situations in which the inter-valley coupling might play an important role as, for instance, in the presence of short-range impurities or defects.. The Hamiltonians at $K$ and $K^{\prime}$ are related by a unitary transformation,

[TABLE]

which amounts to a change of sublattice basis at the valley $K^{\prime}$. Equation (4) implies that the spectra of eigenvalues of ${\mathcal{H}}_{+}$ and ${\mathcal{H}}_{-}$ are the same and that the eigenstates at $K^{\prime}$ are obtained from the eigenstates at $K$ by multiplication by $i\sigma_{y}$. This mapping does not affect any of the observables discussed below. We can then restrict ourselves to a single valley, say $\tau=+1$, and suppress the valley index from now on.

Throughout this paper we shall use dimensionless variables, where lengths are expressed in units of magnetic length $\ell_{B}=\sqrt{\hbar c/eB}$, wave vectors in units of $\ell_{B}^{-1}$, and energies in units of the relativistic cyclotron energy $\hbar\omega_{c}=\hbar v_{\text{F}}/\ell_{B}$. For a magnetic field of $2$ T, one finds $\ell_{B}\simeq 18$ nm and $\hbar\omega_{c}\simeq 34$ meV, while the Zeeman coupling strength is approximately $100\,\mu$eV.

We shall study a linear geometry, see the sketch in Fig. 1, where the potential profile is a sharp antisymmetric step,

[TABLE]

The potential creates an interface at $x=0$ separating a $n$-doped region for $x<0$ from a $p$-doped region for $x>0$. The local in-plane electric field due to the potential step might in principle induce an additional local SOI, but this is negligible due to the assumed smoothness of the electrostatic potential on the scale of the lattice constant 333To estimate the effect, we follow Ref. Kane and Mele (2005) and replace the Coulomb potential with the potential of the p-n junction. Assuming that the potential difference at the junction is of the order of the cyclotron energy and that the potential varies on the scale of the magnetic length, we estimate the induced SOI to be of the order of $10^{-6}$ meV..

In the Landau gauge $\textbf{A}=(0,Bx,0)$ the Hamiltonian (1) is translationally invariant in the $y$-direction, and its eigenstates take the form

[TABLE]

where $\psi_{k}(x)$ is a four component wave function describing the electron amplitudes on A/B sublattices for spin up/down:

[TABLE]

(Here, $T$ denotes transposition and, from now on, we shall omit the index $k$ on the spinors.) After projection onto the plane waves in the $y$-direction, we need to solve the equation

[TABLE]

where ${\cal H}(k)=e^{-iky}{\mathcal{H}}e^{iky}$. It is convenient to introduce the ladder operators

[TABLE]

where $q=\sqrt{2}(x+k)$ is the shifted and rescaled spatial coordinate and $\partial_{q}=\frac{1}{\sqrt{2}}\partial_{x}$ the corresponding differential operator. Then we have

[TABLE]

where we use the notation ($\alpha=\pm$)

[TABLE]

and $b$ denotes the dimensionless Zeeman coupling.

The general solution to Eq. (7) for piecewise constant external fields was obtained in De Martino et al. (2011) (see also Rashba (2009)). To make this paper self-contained, we concisely review the main steps of the solution. The structure of Eq. (9) suggests the following ansatz:

[TABLE]

where $D_{p}(\pm q)$ are the parabolic cylinder functions Olver et al. (2010) and $p$ is a real parameter. The parabolic cylinder functions are the eigenfunctions of the operator $\hat{a}^{\dagger}\hat{a}$, $\hat{a}^{\dagger}\hat{a}D_{p}(\pm q)=pD_{p}(\pm q)$, and satisfy the following recurrence relations:

[TABLE]

The superscript $>$ (resp. $<$) in Eq. (10) indicates that the wave function vanishes for $x\rightarrow+\infty$ (resp. $x\rightarrow-\infty$). Inserting this ansatz into Eq. (7), we get a homogeneous linear system of equations for the coefficients $\{d_{i}\}$:

[TABLE]

The existence of a non-trivial solution requires the vanishing of the determinant of the coefficient matrix, which results in the equation

[TABLE]

where we have introduced the notation

[TABLE]

The solutions are

[TABLE]

By solving Eq. (11) we then obtain

[TABLE]

where $p$ must be replaced by either $p_{+}$ or $p_{-}$. The set $\left\{\psi^{>}_{p_{+}},\psi^{>}_{p_{-}},\psi^{<}_{p_{+}},\psi^{<}_{p_{+}}\right\}$ is a basis in the space of eigenfunctions of ${\mathcal{H}}(k)$ with eigenvalue $E$, and thus the general eigenfunction can be written as

[TABLE]

For a system homogeneous in the $x$-direction (with $V(x)=V_{0}$) the wave functions $\psi^{\gtrless}_{p_{\pm}}(x)$ must be normalizable. This happens only if $p_{\pm}$ are non-negative integers Olver et al. (2010). In this case, $\psi^{<}_{p_{\pm}}$ is proportional to $\psi^{>}_{p_{\pm}}$ and can be omitted. Then, for each given $p_{+}=0,1,2,\dots$ (resp. $p_{-}=0,1,2,\dots$ ) the solution of Eq. (13) for $E$ gives a full spectrum of LLs, in one-to-one correspondence with spinless graphene’s LLs. Therefore, overall we obtain two sets of levels, which reduce to the usual twofold spin-degenerate LLs when $\lambda=\Delta=b=0$. This spectrum is discussed in detail in De Martino et al. (2011). Here, we report the results in two simple limits, where explicit analytical expressions can be obtained. These results will be useful for comparison with the exact spectrum of the p-n junction discussed in the next Section. In the case $\lambda=0$, omitting the obvious shift by $V_{0}$, one finds

[TABLE]

where the index $\pm$ can be identified with the spin up/down index. We note that for $b=0$ the intrinsic SOI splits the spin degeneracy of the zero-energy states, while all other states remain doubly degenerate. In the case $b=\Delta=0$, one finds the non-zero eigenvalues given by ($n=\pm 1,\pm 2,\dots$)

[TABLE]

in addition to two zero-energy levels 444For $n=\pm 1$ one only need to keep the solutions $E_{\pm 1,+}$.. We observe that the Rashba SOI splits the degeneracy of all levels except the zero energy ones. We note, in passing, that when $b=\Delta=0$ the Hamiltonian in Eq. (9) is unitarily equivalent to the Hamiltonian of bilayer graphene in Bernal stacking for a single valley McCann and Fal’ko (2006), with the spin index identified with the layer index and the Rashba SOI identified with the interlayer hopping amplitude. In the limit $|E|\ll\lambda$, from Eq. (12) we obtain $E\approx\pm\frac{2}{\lambda}\sqrt{p(p-1)}$ for the two smallest solutions, which indeed coincide with the LLs of bilayer graphene in the effective two-band approximation when $p$ is a non-negative integer McCann and Fal’ko (2006).

III Exact solution of the p-n junction

In this section we use the general solution in Eq. (15) to construct the exact eigenstates for the p-n junction problem. In the presence of a non-uniform potential the LLs become dispersive, and we shall determine the dispersion relations $E=E(k)$. For the potential step given in Eq. (5) we need to impose the continuity of the wave function at the interface $x=0$ between the $n$-doped region ($x<0$, $V=-V_{0}$) and the $p$-doped region ($x>0$, $V=V_{0}$). Taking into account the normalizability properties of the parabolic cylinder functions, the general form of the wave function is given by

[TABLE]

and the matching condition at $x=0$ reads

[TABLE]

where $p_{\pm}$ are related to the energy $E$ as given in Eq. (13). Therefore we get a homogeneous linear system

[TABLE]

with the matrix of coefficients given by

[TABLE]

where the spinors $\psi^{\gtrless}_{p_{\pm}}(0)$ define the columns of the matrix ${\bf W}$, and

[TABLE]

This linear system admits a non-trivial solution only if the determinant vanishes:

[TABLE]

This condition implicitly defines the exact spectral branches of the p-n junction problem. We study this equation numerically. For given $k$, we find infinitely many solutions $E=E_{n,\pm}(k)$ that can be labelled by a Landau level index $n\in{\mathbb{Z}}$ and an index $\pm$ which, in the limit of vanishing $\lambda$, reduces to a spin projection quantum number 555 We note that, due to the spin degeneracy of all LLs at $\lambda=\Delta=b=0$, the limits $\lambda\rightarrow 0$ and $(b,\Delta)\rightarrow 0$ do not commute: the limit $\lambda\rightarrow 0$ at finite $b$ or $\Delta$ produces eigenstates of $s_{z}$, whereas the limit $\lambda\rightarrow 0$ at $b=\Delta=0$ gives eigenstates of $s_{x}$.. Then, for any given solution $\left\{k,E_{n,\pm}(k)\right\}$ of Eq. (22), we solve the linear system (20) to obtain the corresponding eigenfunction, which can then be normalized. Let us now discuss in detail our findings.

The spectrum consists of infinitely many branches of dispersive LLs. This is illustrated in Fig. 2a for a generic case, with $\lambda$, $\Delta$, and $b$ all finite. Two special cases of particular interest are illustrated in Fig. 2b, where $b=\Delta=0$, and $\lambda$ is finite; and in Fig. 2c, where $\lambda=0$ and $\Delta$ and $b$ are finite and equal.

The flat portions of each level (for $|k|\gg 1$) in Fig. 2 correspond to bulk states. The approximately linearly dispersing portions visible throughout the spectrum correspond to states mainly localized close to the interface. Due to the Rashba SOI, which mixes spin-up and spin-down, these levels have $k$-dependent spin polarization, which we discuss in Sec. IV below. Adjacent levels generally exhibit anticrossings, see Figs. 2a and 2b, unless $\lambda=0$, in which case the projection of spin along the $z$-direction is a good quantum number, and one finds level crossings instead, as in Fig. 2c.

For a not too large amplitude of the potential step, there exists an energy window around zero energy in which there is a single pair of interface modes. This is the situation illustrated in the three panels of Fig. 2. If the Fermi level lies in this energy window, the observable properties of the system will be determined by these linearly dispersing chiral zero modes. In the rest of the paper we will mainly focus on them. As they evolve from the two zero-energy modes of the uniform system, we will refer to them as the zero modes of the p-n junction. Close to zero energy their dispersion is linear to a very good approximation, and they can be regarded as the quantum states corresponding to spin-resolved quasiclassical snake states.

Let us consider now the case $b=\Delta=0$, illustrated in Fig. 2b. Remarkably, the dispersions of the two zero modes exhibit a relative shift in wave vector reminiscent of the horizontal shift between the energy dispersions in a single-channel quantum wire with Rashba SOI Datta and Das (1990); Bercioux and Lucignano (2015). The modes become again degenerate (with exponential accuracy) for $|k|\rightarrow\infty$. We can rationalize this behaviour as follows. In the gauge we are working in, the electronic states can be visualized as stripes along the $y$-direction, localized in the $x$-direction around $x=-k$ within a distance of the order of the magnetic length. For $|k|$ very large, their centre is far from the interface and they are not affected by the potential step. The solution for the uniform case then shows that there are two degenerate zero modes, whose energies do not depend on the Rashba coupling. This occurs because they both have the same vanishing (or more precisely, exponentially small) sublattice component. Since the Rashba SOI couples the amplitudes of electrons with opposite spin on different sublattices, its matrix element will be at most exponentially small, and it will not lift the spin degeneracy. For $|k|\approx 0$, instead, the electronic states are close to the interface and the zero modes have a significant amplitude on both sublattices. Then the Rashba SOI couples them, lifting the degeneracy.

The wave vector splitting between the two zero modes can be calculated analytically in perturbation theory in $\lambda$ (see App. A). We find the simple result

[TABLE]

where $k_{\pm}$ are defined as the solutions of $E_{0,\pm}(k)=E$, and we have reintroduced the units. Interestingly, the splitting is independent of energy, of magnetic field and even of the potential difference across the junction. Using our exact dispersions, we have also computed the exact splitting as a function of $\lambda$ numerically. The exact result is compared with the perturbative one in Fig. 3, which shows that the two are in good agreement up to $\lambda\approx 0.2$.

The p-n junction then hosts at given energy two degenerate modes at different wave vector with finite overlap, and we expect interference effects in the observables induced by the Rashba SOI, which we discuss in the next section. A relative shift of the zero-mode dispersions can be observed also for vanishing $\lambda$ and finite $b$ and $\Delta$. In this case, however, the degeneracy is restored at large $k$ only in the non generic situation that $b=\Delta$ illustrated in Fig. 2c. More importantly, for $\lambda=0$ the two modes are eigenstates of $s_{z}$, thus their overlap vanishes and interference effects are not possible.

IV Observables

From now on, we focus solely on the effects of the Rashba SOI ($\Delta=b=0$) and we restrict to an energy window in which there exist only the two dispersive zero modes. At given energy $E$, the electronic state is then a superposition of snake states at two different wave vectors

[TABLE]

All observables, bilinear in the wave function, will contain a term uniform in $y$ and an interference term oscillating along $y$, with period given by $2\pi\Delta k^{-1}$. These oscillating terms manifest themselves in correlation functions, as we discuss below, provided the overlap between the states at $k_{\pm}$ does not vanish. This happens only if the Rashba SOI is different from zero. If $\lambda=0$, the two modes are spin eigenstates with opposite spin projection and hence are orthogonal. As a consequence, the presence of such oscillating terms in correlation functions is a hallmark of Rashba SOI and provides a way of measuring the strength of the coupling.

The probability density associated to the states at $k_{\pm}$, $P_{\alpha}(x)=\psi^{\dagger}_{k_{\alpha}}(x)\psi_{k_{\alpha}}(x)$ ($\alpha=\pm$) is illustrated in Fig. 4a for the states at $E=0$. One observes that the profiles are localized at $x=0$ within a magnetic length, are skewed, with maxima on opposite sides of the interface. Although the wave functions seem to have a substantial spatial overlap, the overlap integral

[TABLE]

(which turns out to be real) is rather small, as shown in Fig. 4b. It reaches a minimum at $E=0$, but remains finite in the whole energy window. This can be rationalised with the observation that a finite Rashba SOI produces a spin texture in the eigenfunctions, which would otherwise have opposite spin orientations and would be orthogonal. As discussed, above, a finite overlap is important as it is the weight of the oscillating term in the density-density correlator.

Next, we consider the mode velocities. The velocity operator is ${\bf v}=s_{0}\otimes{\bm{\sigma}}$, and the group velocities of the two zero modes are given by

[TABLE]

We find $v_{\pm,x}=0$, while the $y$-components $v_{\pm,y}$, which coincide with the slopes of the zero-mode dispersions $E^{\prime}_{0,\pm}(k_{\pm})$, are illustrated in Fig. 5 as function of energy. In the whole energy range the velocities practically coincide. They are finite and reach their maximum (in absolute value) at $E=0$. Close to $E=0$ they are almost constant, which implies that the dispersions are almost perfectly linear. We have checked numerically that they depend only weakly on $\lambda$. This holds also in the limit $V_{0}\ll\lambda$. In fact, a straightforward perturbative calculation in $V_{0}$ (but exact in $\lambda$) shows that for $k\approx 0$ the velocities can be expressed as

[TABLE]

where $C(x)=\frac{2+x^{2}/2}{\sqrt{\pi}(1+x^{2}/2)}$ and we have restored units. This value does not depend on graphene’s Fermi velocity. Since $V_{0}/e\ell_{B}$ can be interpreted as the electric field $E_{x}$ induced by the potential step in the $x$-direction and $\ell^{2}_{B}\propto 1/B$, $V_{0}\ell_{B}/\hbar$ corresponds to the classical drift velocity of a charged particle in crossed magnetic and electric fields, $cE_{x}/B$ Cohnitz et al. (2016). The effect of the Rashba SOI is a renormalization factor which decreases for increasing $\lambda$ from $2/\sqrt{\pi}$ (for $\lambda=0$) to $\frac{1}{\sqrt{\pi}}$ (for $\lambda\gg\hbar\omega_{c}$) and thus tends to slightly reduce the velocities.

Next, we consider the spin polarization. The spin operator is given by ${\bf S}={\bf s}\otimes\sigma_{0}$ (for simplicity we omit the factor $1/2$), so the expectation value of the spin for the two zero modes is given by

[TABLE]

We find that $S_{y,\pm}=0$ and that ${\bf S}_{\pm}$ points predominantly in the $x$-direction. This can be rationalized as due to the fact that the Rashba SOI induces a spin precession in a direction perpendicular to the direction of motion and, in the geometry considered here, it selects the eigenstates of $s_{x}$ (see discussion at the end of App. A). The $z$-component switches sign across $E=0$ and we observe the symmetry properties $S_{z,-}(E)=S_{z,+}(-E)$ and $S_{x,-}(E)=-S_{x,+}(-E)$.

The matrix element of the spin operator between the two zero modes is given by

[TABLE]

and is illustrated in Fig. 6. It plays an important role in the spin-spin correlation function, because it is the weight of the interference term. We find that $\Sigma_{x,\pm}$ and $\Sigma_{z,\pm}$ are real, with $\Sigma_{x(z),+}=\Sigma_{x(z),-}$, while $\Sigma_{y,\pm}$ are imaginary, with $\Sigma_{y,+}=-\Sigma_{y,-}$.

We turn now to the density-density and spin-spin correlations. While we have access in principle to the full exact correlators, for the purpose of a study of the low-energy properties around $E=0$, it is sufficient to focus on the linear part of the spectrum and describe our system as a pair of one-dimensional non-interacting chiral electronic modes, with (second-quantized) Hamiltonian given by

[TABLE]

and the fermion field operator given by

[TABLE]

where $\psi_{\alpha}(x)$ are the exact eigenstates of the zero modes at $E=0$ and $L_{\alpha}(y)$ fermionic operators with usual anti-commutation relations:

[TABLE]

Using the expression for the field operator (31), the effective one-dimensional density operator, obtained by integrating the two-dimensional density operator over $x$, is given by

[TABLE]

with

[TABLE]

and the overlap integral ${\cal O}$, defined in Eq. (25), is evaluated at $E=0$. Therefore, the density-density correlator will contain two contributions:

[TABLE]

The uniform term is the usual correlator of one-dimensional chiral fermions:

[TABLE]

The interference term is given by

[TABLE]

As anticipated, the density-density correlator exhibits spatial oscillations with a period proportional to the inverse of $\Delta k$. The structure factor, which is essentially its Fourier transform, will then exhibit a pronounced observable feature at $q=\Delta k$.

The spin-spin correlator has a similar structure. The effective one-dimensional spin-density operator is given by

[TABLE]

where the matrix elements ${\bf S}_{\alpha}$ and ${\bf\Sigma}_{\alpha}$, defined in Eqs. (28) and (29), are evaluated at $E=0$. The spin-spin correlator then will be

[TABLE]

Therefore the spin structure factor will present a feature at momentum $q=\Delta k$, indicating the existence of an interference term.

V Conclusions and Outlook

In this paper we have studied the effects of finite spin-orbit interaction on the low-energy electronic properties of p-n junctions in graphene subject to a perpendicular magnetic field. We have found the exact solution of the appropriate Dirac-Weyl equation and thereby determined the exact spectrum of dispersive LLs. The spectrum contains a pair of propagating modes which cross zero energy with approximately linear dispersion, and can be interpreted as quantum limit of snake states. When the Fermi level lies close to $E=0$ between two bulk LLs, the transport properties are dominated by these snake states, similar to those forming in usual graphene p-n junctions in magnetic field. However, the Rashba SOI induces a mismatch in the Fermi momenta of these states, which has observable consequences. In fact, we predict that the density-density and the spin-spin correlation functions have a spatially oscillating behaviour along the junction, which could be tested experimentally.

The relative horizontal shift of the dispersions of the interface zero modes and the $k$-dependence of their spin polarization are strongly reminiscent of the analogous effects in a single-channel quantum wire with Rashba SOI. This analogy suggests that the p-n junction works as a Datta-Das channel. If one injects electrons with definite spin polarization at some point along the junction, and measure correlations downstream, in the presence of a substantial Rashba SOI one should observe oscillations due to interference effects between the two modes. The junction could then provide a realization of a spin-FET. With respect to the case of a quantum wire, this realization would have the advantage that the system is chiral, and therefore it is practically immune to backscattering by impurities. Moreover, although one would expect, in principle, large effects of electron-electron interactions on these effectively one-dimensional electronic states, the chirality of the system suppresses the backscattering processes and prevents the opening of a gap at the Fermi level. The e-e interactions might then only renormalize the velocities and modify the exponents in the power law decay of correlations, but could not spoil the effects discussed here. We leave the detailed study of the e-e interaction effects for future works.

We hope that our work will stimulate further experimental and theoretical activity on the properties of snake states in graphene’s p-n junction, and that the predictions put forward here can be soon tested experimentally.

Acknowledgements.

We thank Reinhold Egger and Tineke van den Berg for a careful reading of the manuscript and useful suggestions. The work of DB is supported by Spanish Ministerio de Ciencia, Innovation y Universidades (MICINN) under the project FIS2017-82804-P, and by the Transnational Common Laboratory Quantum-ChemPhys.

Appendix A Perturbation theory

In this appendix we briefly discuss the calculation of the spectrum of the Hamiltonian in Eq. (1) by perturbation theory in $\lambda$. We use this result to find an explicit formula for the dependence of the splitting $\Delta k$ on the Rashba SOI.

For $\lambda=0$ the Hamiltonian (1) is block-diagonal in spin. Each block coincides with the Hamiltonian studied in Cohnitz et al. (2016) (up to the intrinsic SOI term and the Zeeman term, which can be easily incorporated) and we can borrow their results. The wave functions can be expressed as

[TABLE]

where $|s\rangle$ are the eigenstates of $s_{z}$ and

[TABLE]

Here $\epsilon(x)=\text{sign}(x)$, ${\cal N}$ is a normalization constant and

[TABLE]

The wave function (37) is continuous at $x=0$ provided the matching condition

[TABLE]

is satisfied. This condition determines the eigenvalues $E$ and plugging the eigenvalues back into Eq. (37) one obtains the corresponding eigenstates $\phi_{k,\uparrow}$ . The eigenstates $\phi_{k,\downarrow}$ and the corresponding eigenvalues are obtained by replacing $\Delta$ and $b$ with $-\Delta$ and $-b$.

Henceforth, we set $\Delta=b=0$. Then $\phi_{k,\uparrow}=\phi_{k,\downarrow}$ and we can omit the spin index. All eigenstates (at given $k$) are doubly degenerate. We focus on the two zero modes, with unperurbed energy $E_{0}(k)$. Using standard degenerate perturbation theory, we find that the first order perturbative correction to the energy due to the Rashba SOI is given by

[TABLE]

where $h(x)$ is the matrix element of the perturbing operator

[TABLE]

Up to the factor $-\lambda/2$, this expression is nothing else than the $y$-component of the group velocity $v^{(0)}_{y}(k)$ of the unperturbed states. Therefore, $h(k)=-\frac{\lambda}{2}v^{(0)}_{y}(k)$ and the mode dispersions read

[TABLE]

The difference in wave vector between the two modes at energy $E$ is then simply given by

[TABLE]

where $k_{\pm}$ are defined by the equation $E_{0,\pm}(k_{\pm})=E$. In other words, $\lambda$ produces a relative horizontal shift of the dispersions. Remarkably, at this order in $\lambda$, the shift is independent of the energy, the potential strength, and the magnetic field. The unperturbed eigenstates that diagonalize ${\mathcal{H}}_{\text{RSO}}$ are the eigenstates of $s_{x}$:

[TABLE]

This is consistent with the fact that our exact eigenstates have spin polarization predominantly along the $x$-direction, see Fig. 6, and explains the observation that at $b=\Delta=0$ they reduce to eigenstates of $s_{x}$ in the limit $\lambda\rightarrow 0$.