# Magnetoelectric torque and edge currents in spin-orbit coupled graphene
nanoribbons

Matheus S. M. de Sousa  Department of Physics, PUC-Rio, 22451-900 Rio de
Janeiro, Brazil  Manfred Sigrist  Institute for Theoretical Physics, ETH
Zurich, 8093 Zurich, Switzerland  Wei Chen  Department of Physics, PUC-Rio,
22451-900 Rio de Janeiro, Brazil

###### Abstract

For graphene nanoribbons with Rashba spin-orbit coupling, the peculiar
magnetic response due to the presence of a magnetization and geometric
confinement are analyzed within a tight-binding model. We observe a sizable
transverse susceptibility that can be considered as a gate voltage-induced
magnetoelectric torque without the need of a bias voltage, with different
directions for zigzag and armchair ribbons. The local torque generates non-
collinear spin polarization between the two edges and/or along the ribbon, and
the net torque averages to zero if the magnetization is homogeneous.
Nevertheless, a nonzero net torque can appear in partially magnetized
nanoribbons or in nanoflakes of irregular shapes. The equilibrium spin current
produced by the spin-orbit coupling also appears in nanoribbons, but the
component flowing in the direction of confinement is strongly suppressed. Even
without the magnetization, an out-of-plane polarized chiral edge spin current
is produced, resembling that in the quantum spin Hall effect. Moreover, a
magnetization pointing perpendicular to the edge produces a laminar flow of
edge charge currents, whose flow direction is symmetric (non chiral) or
antisymmetric (chiral) between the two edges depends on whether the
magnetization points in-plane or out-of-plane.

##  I Introduction

The celebrated Rashba spin-orbit coupling (RSOC) has a strong impact on the
physical properties of two-dimensional (2D) metals with parabolic bandsBychkov
and Rashba (1984a, b), especially their magnetic response and charge to spin
interconversion. A number of these features originate from the RSOC induced
spin-momentum locking, such as the Edelstein effect that causes a bias
voltage-induced nonequilibrium in-plane spin polarizationEdelstein (1990). In
magnetized 2D systems, this effect can be exploited to induce magnetization
dynamics known as the spin-orbit torqueManchon and Zhang (2008, 2009); Haney
and Stiles (2010); Gambardella and Miron (2011); Pesin and MacDonald (2012a);
Haney  _et al._ (2013) , whose feasibility has been demonstrated extensively
in experimentsMiron  _et al._ (2010, 2011); Manchon  _et al._ (2019). In
addition to these properties, RSOC also modifies the equilibrium properties of
2D metals, most notably causing an in-plane polarized spin current flowing
throughout the systemRashba (2003); Sonin (2007); Tokatly (2008).

Besides these effects occurring in infinitely extended 2D metals or electron
gases (2DEG), RSOC also causes peculiar effects at the boundary of
geometrically confined mesoscopic 2D systems. A notable example is the
generation of out-of-plane polarized equilibrium edge spin currentsReynoso
_et al._ (2004); Grigoryan  _et al._ (2009). Furthermore, in-plane spin
polarization perpendicular to the edge can induce a persistent charge current
decaying and oscillating in sign away from the edgeUsaj and Balseiro (2005).
This behavior is very similar to that found in other spin-momentum locking
systems in proximity to a magnet, such as for topological
insulator/ferromagnetic metal junctions (TI/FMMs)Zegarra  _et al._ (2020);
Chen (2020a). The tunability of RSOC by a gate voltage may provide a means to
control such phenomena and engineer dedicated devices.

In this paper, we use a tight-binding model to explore the spin torque and
equilibrium currents caused by RSOC and geometric confinement in graphene
nanoribbons. Numerous peculiarities make graphene nanoribbons a particularly
intriguing system to study: (1) the linear dispersion near the two Dirac
points that causes a unique spin-momentum locking profileNovoselov  _et al._
(2005); Zhou  _et al._ (2006); Geim and Novoselov (2007); Castro Neto  _et
al._ (2009); Geim (2009); Das Sarma  _et al._ (2011); (2) the emergence of
zero energy edge states in zigzag terminated ribbons Fujita  _et al._ (1996);
Nakada  _et al._ (1996); Brey and Fertig (2006); Ezawa (2006); Peres  _et al._
(2006); Akhmerov and Beenakker (2008); Tao  _et al._ (2011); Ruffieux  _et
al._ (2011): (3) the opening of energy gap depending on the ribbon width for
armchair ribbons Son  _et al._ (2006a, b); Barone  _et al._ (2006); Yang  _et
al._ (2007); Han  _et al._ (2007); Ritter and Lyding (2009), and (4) the
pronounced magnetic response at some sample edges Wakabayashi  _et al._
(1999); Lee  _et al._ (2005); Pisani  _et al._ (2007). Moreover, it is
important to note that an enhanced and gate-tunable spin-orbit coupling (SOC)
has been demonstrated in graphene/transition-metal dichalcogenide (TMD)
heterostructuresYang  _et al._ (2016); Wang  _et al._ (2016); Yang  _et al._
(2017); Safeer  _et al._ (2019); Ghiasi  _et al._ (2019); Benítez  _et al._
(2020). In graphene/ yttrium iron garnet (YIG)Mendes  _et al._ (2015); Wang
_et al._ (2015); Dushenko  _et al._ (2016); Leutenantsmeyer  _et al._ (2016)
and graphene/CoRybkin  _et al._ (2018) heterostructures, both SOC and
ferromagnetism are induced, which might be used as the stage for the phenomena
we are going to discuss. Finally, the RSOC and magnetization are also known to
cause quantum spin Hall (QSHE) and quantum anomalous Hall effects (QAHE) in
graphene when combined with the intrinsic SOCDiniz  _et al._ (2013); Guassi
_et al._ (2015).

We first elaborate on the feature that, in contrast to the usual longitudinal
in-plane susceptibility in an infinitely large graphene, a transverse in-plane
susceptibility develops in nanoribbons due to the RSOC and geometric
confinement. The controllability of the RSOC suggests that this transverse
susceptibility serves as a gate-voltage-induced magnetoelectric torque on the
magnetization without the need of a bias voltage, contrary to the current-
induced spin-orbit torqueDyrdał and Barnaś (2015); Li and Manchon (2016);
Rodriguez-Vega  _et al._ (2017); Rybkin  _et al._ (2018); Zollner  _et al._
(2020); Rybkina  _et al._ (2020). While this torque averages to zero for a
homogeneously spin polarized system, well designed local variations in the
magnetization can yield a net torque signal available to practical purpose. We
further investigate the pattern of equilibrium spin currents in nanoribbons,
especially concerning the symmetry of flow directions between the two edges,
e.g. distinguished chiral and non-chiral structure. Charge currents induced by
spin polarization exist in nanoribbons too, whose variability upon the
direction of the magnetization will be a point of attention.

We start by introducing the lattice model for graphene with RSOC and spin
magnetization, and then calculate the persistent spin current to show how the
spin-momentum locking is modified by the presence of the magnetization in Sec.
II. In Sec. III, we use a two-site toy model to demonstrate analytically the
existence of a transverse susceptibility due to geometric confinement, which
is then transferred to the discussion of both zigzag and armchair ribbons. In
addition, we characterize the patterns and chiralities of persistent charge
and spin currents, and also use an specific example to demonstrate that these
phenomena also survive in graphene nanoflakes. Section IV summarizes our
results.

##  II Extended tight-binding model of graphene

In order to emphasize the impact of geometric confinement, we start by
addressing the spintronics related effects due to RSOC for an infinite
graphene sheet, before discussing the same physics for graphene ribbons and
nanoflakes. As indicated in Fig. 1 (a), three lattice vectors (in units of
bond length $a=1$) characterize the honeycomb lattice,

$\displaystyle{\boldsymbol{\delta}}_{1}=\left(\frac{1}{2},\frac{\sqrt{3}}{2}\right),\;\;\;{\boldsymbol{\delta}}_{2}=\left(\frac{1}{2},-\frac{\sqrt{3}}{2}\right),\;\;\;{\boldsymbol{\delta}}_{3}=\left(-1,0\right).\;\;\;$
(1)

connecting neighbouring lattice sites belonging to the two different
sublattices, $A$ and $B$. We now formulate the following tight-binding model
that incorporates both RSOC and spin polarizationKane and Mele (2005a, b); Min
_et al._ (2006); Dyrdał and Barnaś (2017); Peralta  _et al._ (2019)

$\displaystyle H$ $\displaystyle=$ $\displaystyle-t\sum_{\langle
ij\rangle,\sigma}c_{i\sigma}^{{\dagger}}c_{j\sigma}+J_{ex}\sum_{i,\alpha,\beta}{\bf
S}\cdot c_{i\alpha}^{{\dagger}}{\boldsymbol{\sigma}}_{\alpha\beta}c_{i\beta}$
(2) $\displaystyle+$ $\displaystyle i\lambda_{R}\sum_{\langle
ij\rangle,\alpha,\beta}c_{i\alpha}^{{\dagger}}\left({\boldsymbol{\sigma}}_{\alpha\beta}\times{\bf
d}_{ij}\right)^{z}c_{j\beta}-\mu\sum_{i\sigma}c_{i,\sigma}^{{\dagger}}c_{i\sigma}.$

Here $c^{{\dagger}}_{i\sigma}(c_{i\sigma})$ creates (annihilates) an electron
of spin $\sigma$ on the lattice site $i$, where sum with $\langle ij\rangle$
run over nearest neighbor lattice sites. The hopping matrix element is $-t$,
$\lambda_{R}$ is the RSOC coupling constant,
${\boldsymbol{\sigma}}=(\sigma^{x},\sigma^{y},\sigma^{z})$ are the spin Pauli
matrices, ${\bf d}_{ij}$ is the vector connecting the site $i$ to $j$,
$J_{ex}$ is the exchange coupling between the magnetization ${\bf
S}=S(\sin\theta\cos\phi,\sin\theta\sin\phi,\cos\theta)$ and the spin, and
$\mu$ is the chemical potential. The bipartite structure of our system allows
us to use a sublattice formulation with $I=\left\\{A,B\right\\}$, where we
define the basis
$\psi=\left(A\uparrow,B\uparrow,A\downarrow,B\downarrow\right)$, and the
electron operators in Eq. (2) can be split into
$c_{i\sigma}\rightarrow\left\\{c_{Ai\sigma},c_{Bi\sigma}\right\\}$, where $i$
now denotes the position of each unit cell containing an $A$ and $B$ lattice
site. After Fourier transformation $c_{Ii\sigma}=\sum_{\bf k}e^{i{\bf k\cdot
r}_{i}}c_{I{\bf k}\sigma}$, the Hamiltonian $H=\sum_{\bf
kIJ\alpha\beta}c_{I{\bf k}\alpha}^{{\dagger}}H_{I\alpha J\beta}({\bf
k})c_{J{\bf k}\beta}$ defines the matrix

$\displaystyle H_{I\alpha J\beta}({\bf k})=$
$\displaystyle\left(\begin{array}[]{cccc}J_{ex}S^{z}&tZ^{\ast}&J_{ex}S_{\perp}e^{-i\phi}&\lambda_{R}Y^{\ast}\\\
tZ&J_{ex}S^{z}&\lambda_{R}X^{\ast}&J_{ex}S_{\perp}e^{-i\phi}\\\
J_{ex}S_{\perp}e^{i\phi}&\lambda_{R}X&-J_{ex}S^{z}&tZ^{\ast}\\\
\lambda_{R}Y&J_{ex}S_{\perp}e^{i\phi}&tZ&-J_{ex}S^{z}\end{array}\right),\;\;\;\;$
(7)

with

$\displaystyle Z\equiv e_{1}+e_{2}+e_{3},$ $\displaystyle
X\equiv\frac{-1-i\sqrt{3}}{2}e_{1}^{\ast}+\frac{-1+i\sqrt{3}}{2}e_{2}^{\ast}+e_{3}^{\ast},$
$\displaystyle
Y\equiv\frac{1+i\sqrt{3}}{2}e_{1}+\frac{1-i\sqrt{3}}{2}e_{2}-e_{3},$

where we have defined $e_{a}=e^{i{\bf k}\cdot{\boldsymbol{\delta}}_{a}}$ and
$S_{\perp}=\sqrt{(S^{x})^{2}+(S^{y})^{2}}=S\sin\theta$ ($\theta\in[0,\pi]$).
In the following we choose the parameters

$\displaystyle
J_{ex}=0.2,\;\;\lambda_{R}=0.2,\;\;\mu=0.5,\;\;k_{B}T=0.03,\;\;\;\;$ (9)

with $t=1$ as energy unit, such that the magnetic response is more pronounced,
but we emphasize that their patterns are fairly robust against variation of
these parameters.

Figure 1: (a) The definition of coordinate in our honeycomb lattice. (b) The
resulting band structure in the presence of RSOC without magnetization ${\bf
S}={\bf 0}$. (c) The spin-momentum locking profile of the lowest eigenstate
$|u_{{\bf k}1}\rangle$ in the absence ${\bf S}={\bf 0}$ and presence ${\bf
S}\parallel\left\\{{\hat{\bf x}},{\hat{\bf y}},{\hat{\bf z}}\right\\}$ of the
magnetization.

Diagonalizing the Hamiltonian yields two filled bands $E_{{\bf k}1}<E_{{\bf
k}2}<0$ and two empty bands $0<E_{{\bf k}3}<E_{{\bf k}4}$, with the
corresponding four eigenstates $\left\\{|u_{{\bf k}1}\rangle,|u_{{\bf
k}2}\rangle,|u_{{\bf k}3}\rangle,|u_{{\bf k}4}\rangle\right\\}$. These four
bands are shown in Fig. 1 (b) for zero magnetization. In this case the spin
expectation value $\langle{\boldsymbol{\sigma}}\rangle_{{\bf k}\eta}=\langle
u_{{\bf k}\eta}|{\boldsymbol{\sigma}}|u_{{\bf k}\eta}\rangle$ for each
eigenstate has only inplane components $\left\\{\langle\sigma^{x}\rangle_{{\bf
k}\eta},\langle\sigma^{y}\rangle_{{\bf k}\eta}\right\\}$, as displayed in the
vector plots in Fig. 1 (c) for $S=0$, where we also see
$\langle{\boldsymbol{\sigma}}\rangle_{{\bf
k}1}=-\langle{\boldsymbol{\sigma}}\rangle_{{\bf k}2}$. We see that the spin
texture forms a vortex around each Dirac point, ${\bf
K}=(2\pi/3,2\pi/3\sqrt{3})$ and ${\bf K^{\prime}}=(2\pi/3,-2\pi/3\sqrt{3})$
with the same vorticity, but opposite vorticity appears for the vortex around
the origin ${\bf k}=(0,0)$. This spin pattern represents the spin-momentum
locking due to the RSOC.

As common in 2D Rashba metals, the presence of a magnetization ${\bf S}$ is
expected to distort the profile of spin pattern Manchon and Zhang (2009);
Gambardella and Miron (2011), as demonstrated for magnetizations along the
three principle directions ${\bf S}\parallel\left\\{{\hat{\bf x}},{\hat{\bf
y}},{\hat{\bf z}}\right\\}$ in Fig. 1 (c), which shows the planar components
of $\langle{\boldsymbol{\sigma}}\rangle_{{\bf k}1}$ only. For ${\bf
S}\parallel{\hat{\bf x}}$, the $\langle{\boldsymbol{\sigma}}\rangle_{{\bf
k}1}$ becomes asymmetric between $+k_{y}$ and $-k_{y}$, which also causes an
asymmetry in the dispersion in this direction. Analogously, ${\bf
S}\parallel{\hat{\bf y}}$, yields an asymmetry between $+k_{x}$ and $-k_{x}$,
whereas for ${\bf S}\parallel{\hat{\bf z}}$, the planar components basically
remain the same as that in the absence of magnetization. By calculating the
spin polarization in our lattice model with periodic boundary condition (PBC),
we also verify that the spin polarization on either sublattice is always
longitudinal, i.e., along the direction of magnetization.

The charge and spin current density operators are constructed from the
equations of motion of the charge and spin density operator, respectively,
which correspond to continuity equationsZegarra  _et al._ (2020); Chen (2020a)

$\displaystyle\dot{n}_{i}$ $\displaystyle=$
$\displaystyle\frac{i}{\hbar}\left[H,n_{i}\right]=\frac{i}{\hbar}\left[H_{t}+H_{R},n_{i}\right]$
$\displaystyle=$ $\displaystyle-{\boldsymbol{\nabla}}\cdot{\bf
J}_{i}^{0}=-\frac{1}{a}\sum_{\eta}J_{i,i+\eta}^{0},$
$\displaystyle\dot{m}_{i}^{a}$ $\displaystyle=$
$\displaystyle\frac{i}{\hbar}\left[H,m_{i}^{a}\right]=\frac{i}{\hbar}\left[H_{t}+H_{R},m_{i}^{a}\right]+\frac{i}{\hbar}\left[H_{J},m_{i}^{a}\right]$
(10) $\displaystyle=$ $\displaystyle-{\boldsymbol{\nabla}}\cdot{\bf
J}_{i}^{a}+\tau_{i}^{a}$ $\displaystyle=$
$\displaystyle-\frac{1}{a}\sum_{\eta}J_{i,i+\eta}^{a}+\frac{2J_{ex}}{\hbar}\left({\bf
S}\times c_{i\alpha}{\boldsymbol{\sigma}}_{\alpha\beta}c_{i\beta}\right)^{a},$

where
${\boldsymbol{\eta}}=\left(-{\boldsymbol{\delta}_{1}},-{\boldsymbol{\delta}_{2}},-{\boldsymbol{\delta}_{3}}\right)$
if $i\in A$, and
${\boldsymbol{\eta}}=\left({\boldsymbol{\delta}_{1}},{\boldsymbol{\delta}_{2}},{\boldsymbol{\delta}_{3}}\right)$
if $i\in B$. Moreover, $\tau_{i}^{a}=\frac{2J_{ex}}{\hbar}\left({\bf S}\times
c_{i\alpha}{\boldsymbol{\sigma}}_{\alpha\beta}c_{i\beta}\right)^{a}$ denotes
the spin torque term as originating from the usual Landau-Lifshitz dynamics.
The calculation of the commutators in Eq. (10) is detailed in Appendix A.

##  III Graphene nanoribbons with magnetization and Rashba spin-orbit coupling

###  III.1 Two-site toy model

For a deeper insight into the impact of geometrical confinement on properties
of graphene with RSOC and magnetization, we first present an exactly solvable
two-site toy model to demonstrate the feature of the transverse spin
susceptibility and the persistent spin current. Similar models have been
proposed to explain the microscopic mechanism of noncollinear magnetic order
and equilibrium spin currentsKatsura  _et al._ (2005); Bruno and Dugaev
(2005); Chen  _et al._ (2014), while we will put emphasis on the spin torque
here. The two-site model Hamiltonian reads,

$\displaystyle H$ $\displaystyle=$
$\displaystyle\sum_{\sigma}t\left(c_{A\sigma}^{{\dagger}}c_{B\sigma}+c_{B\sigma}^{{\dagger}}c_{A\sigma}\right)+J_{ex}\sum_{i=A,B}{\bf
S}_{i}\cdot
c_{i\alpha}^{{\dagger}}{\boldsymbol{\sigma}}_{\alpha\beta}c_{i\beta}$ (11)
$\displaystyle+$ $\displaystyle
i\lambda_{R}\sum_{ij}c_{i\alpha}^{{\dagger}}\left({\boldsymbol{\sigma}}_{\alpha\beta}\times{\bf
d}_{ij}\right)^{z}c_{j\beta},$

where the two sites $i=\left\\{A,B\right\\}$ are assumed to be connected along
${\bf d}_{AB}=-{\bf d}_{BA}\parallel{\hat{\bf x}}$, and we consider the same
magnetization on the two sites ${\bf S}_{A}={\bf S}_{B}$. The two sites
feature the $A$ and $B$ sublattice sites in the graphene unit cell connected
along ${\hat{\bf x}}$ direction, as implied in Fig. 1. The $4\times 4$
Hamiltonian in the basis of
$(c_{A\uparrow},c_{B\uparrow},c_{A\downarrow},c_{B\downarrow})$ is

$\displaystyle
H=\left(\begin{array}[]{cccc}J_{ex}S^{z}&t&J_{ex}S_{\perp}e^{-i\phi}&-\lambda_{R}\\\
t&J_{ex}S^{z}&\lambda_{R}&J_{ex}S_{\perp}e^{-i\phi}\\\
J_{ex}S_{\perp}e^{i\phi}&\lambda_{R}&-J_{ex}S^{z}&t\\\
-\lambda_{R}&J_{ex}S_{\perp}e^{i\phi}&t&-J_{ex}S^{z}\end{array}\right),$ (16)

analogous to Eq.(7). Assuming $\left\\{J_{ex},\lambda_{R},t\right\\}>0$ and
$\left\\{J_{ex},\lambda_{R}\right\\}\ll t$, for any magnetization direction
there are two negative ($E_{1},E_{2}$) and two positive ($E_{3},E_{4}$)
eigenenergies. Suppose the two negative eigenstates $|u_{1}\rangle$ and
$|u_{2}\rangle$ are occupied. Then the $\alpha$ component spin expectation
value on site $i$ is given by $\langle\sigma_{i}^{\alpha}\rangle=\langle
u_{1}|\sigma_{i}^{\alpha}|u_{1}\rangle+\langle
u_{2}|\sigma_{i}^{\alpha}|u_{2}\rangle$, which can be expanded in powers of
$J_{ex}S$ to obtain the susceptibility along principle directions
$\chi_{i}^{\alpha\beta}=\partial\langle\sigma_{i}^{\alpha}\rangle/\partial(J_{ex}S^{\beta})$.
We find

$\displaystyle\chi_{A}^{\alpha\beta}=\left(\begin{array}[]{ccc}\chi_{A}^{xx}&\chi_{A}^{xy}&\chi_{A}^{xz}\\\
\chi_{A}^{yx}&\chi_{A}^{yy}&\chi_{A}^{yz}\\\
\chi_{A}^{zx}&\chi_{A}^{zy}&\chi_{A}^{zz}\end{array}\right)$ (21)
$\displaystyle=\frac{1}{(t^{2}+\lambda_{R}^{2})^{3/2}}\left(\begin{array}[]{ccc}-\lambda_{R}^{2}&0&t\lambda_{R}\\\
0&0&0\\\
-t\lambda_{R}&0&-\lambda_{R}^{2}\end{array}\right)=\left(\chi_{B}^{\alpha\beta}\right)^{T}.$
(25)

Thus there is a transverse response $\left\\{\chi^{xz},\chi^{zx}\right\\}$
between $x$ and $z$ directions, and consequently a spin torque on both sites
due to Landau-Lifshitz dynamics,

$\displaystyle\frac{d{\bf
S}_{i}}{dt}=\frac{J_{ex}}{\hbar}\langle{\boldsymbol{\sigma}}_{i}\rangle\times{\bf
S}_{i}.$ (26)

However, because the transverse susceptibility is opposite on the two sites
$\chi_{A}^{xz}=-\chi_{B}^{xz}$, the net torque $\sum_{i=A,B}d{\bf S}_{i}/dt$
vanishes. If we assume a coupling of the two magnetic moments by some exchange
interaction $J_{AB}{\bf S}_{A}\cdot{\bf S}_{B}$, then the transverse spin
polarization causes a canting angle $\theta_{AB}\approx
2\cos^{-1}\left(J_{ex}^{2}\lambda_{R}/2J_{AB}t^{2}\right)$ between them in the
$xy$-plane, realizing a Dzyaloshinskii–-Moriya interaction (DMI)Dzyaloshinsky
(1958); Moriya (1960). Note that $\chi_{A}^{xx}$ and $\chi_{A}^{zz}$ in Eq.
(25) represent the corrections to the longitudinal spin polarization due to
RSOC, which must be negative, since the spins would be fully polarized
$\langle{\boldsymbol{\sigma}}_{i}\rangle=1$ along ${\bf S}$ if RSOC were
absent. The $\chi_{A}^{yy}=0$ means the magnetic response in the ${\hat{\bf
y}}$ direction is unaffected.

The local charge and spin current operators in this toy model can be
constructed from the same formalism used in Sec. II. Denoting the current
operator $J_{AB}^{a}$ as the one flowing from site $A$ to $B$ and $J_{BA}^{a}$
for the opposite direction, we find that for either ${\bf S}\parallel{\hat{\bf
x}}$ or ${\bf S}\parallel{\hat{\bf z}}$, the only nonzero expectation value is
the spin current $\langle J^{y}\rangle$, and it satisfies

$\displaystyle\langle J_{AB}^{y}\rangle=-\langle J_{BA}^{y}\rangle$
$\displaystyle=\frac{J_{ex}S\lambda_{R}}{\sqrt{(t-J_{ex}S)^{2}+\lambda_{R}^{2}}}-\frac{J_{ex}S\lambda_{R}}{\sqrt{(t+J_{ex}S)^{2}+\lambda_{R}^{2}}}.$
(27)

There is no spin current for ${\bf S}\parallel{\hat{\bf y}}$. For all cases
the charge current is always absent, but all three components of the spin
current can be non-zereo, in general. Moreover, the $x$ and $z$ components of
the spin currents flowing in the two directions are not negative of each other
$\langle J_{AB}^{x,z}\rangle\neq-\langle J_{BA}^{x,z}\rangle$, but the
continuity equation
$\langle\dot{m}_{A}^{a}\rangle=\langle\dot{m}_{B}^{a}\rangle=0$ is explicitly
satisfied, if the spin torque is taken into account as in Eq. (10).

Figure 2: The band structure $E(n_{x},k_{y})$, persistent charge and spin
currents $J^{\alpha}$, and transverse spin polarization $\sigma^{\alpha}$ (red
positive and blue negative) in the mirror symmetric zigzag ribbon with RSOC at
different magnetization directions ${\bf S}$. The pattern of these quantities
repeats along the ribbon direction ${\hat{\bf y}}$, and is either symmetric
(labeled SY) or antisymmetric (labeled AS) under $P_{\rm zig}^{\rm mir}$ in
Eq. (28), so we only present the left half of a unit cell, except for ${\bf
S}={\bf 0}$ we show the entire unit cell. For the absence of magnetization
case ${\bf S}={\bf 0}$, the colors in the band structure indicate the spin up
and down components at the left half space (LE) for the eigenstates
$|n_{x},k_{y}\rangle$, and for ${\bf S}\neq{\bf 0}$ cases they indicate the
probability of the eigenstates at the left and right (RE) half space. The
largest arrow in the figure has magnitude $0.015$, and the largest disk
$0.011$.

###  III.2 Persistent currents and spin torques in Zigzag ribbons

The results in the two-site model suggests the existence of spin torque and
equilibrium spin current in geometrically confined Rashba systems, which
motivates us to investigate our graphene model in Sec. II in the nanoribbon
geometry. We first consider the zigzag ribbon by taking open boundary
condition (OBC) in the ${\hat{\bf x}}$ direction and PBC along ${\hat{\bf y}}$
direction for our lattice model in Eq. (2). There are two kinds of zigzag
ribbons as far as the symmetry between the two edges is concerned, namely the
mirror symmetric and the glide plane symmetric ones, each invariant under the
corresponding symmetry operation

$\displaystyle P_{\rm zig}^{\rm mir}=\left\\{P_{x}|\,{\bf
0}\right\\},\;\;\;P_{\rm zig}^{\rm gli}=\left\\{P_{x}|\,\sqrt{3}/2{\hat{\bf
y}}\right\\},$ (28)

using the common notation of space group operation, where $P_{x}$ denotes the
mirror reflection in ${\hat{\bf x}}$ directions with respect to the central
axis of the ribbon. We find that the symmetries of the patterns of the
currents $J^{\alpha}$ and spin polarization $\sigma^{\alpha}$ in the mirror
symmetric ribbon defined with respect to $P_{\rm zig}^{\rm mir}$ are the same
as those in the glide plane symmetric ribbon defined with respect to $P_{\rm
zig}^{\rm gli}$, so we only present the former ones for simplicity. In
particular, we use the zigzag ribbon of 24-site width as an example, as shown
in Fig. 2. The patterns are translationally invariant along the ribbon
direction ${\hat{\bf y}}$, and are either symmetric (labeled SY) or
antisymmetric (labeled AS) under $P_{\rm zig}^{\rm mir}$, such that we only
show the left half of a unit cell, with size of the arrows and disks
indicating the magnitude of the currents and spin polarization, respectively.
In addition, to clarify the origin of the edge charge current, for each
eigenstate $|n_{x},k_{y}\rangle$ we calculate the weight of the wave function
closer to the left (LE) and right (RE) edges

$\displaystyle n_{n_{x},k_{y}}^{LE}=\sum_{1\leq x\leq
N_{x}/2}|\psi_{n_{x},k_{y}}(x)|^{2},$ $\displaystyle
n_{n_{x},k_{y}}^{RE}=\sum_{N_{x}/2+1\leq x\leq
N_{x}}|\psi_{n_{x},k_{y}}(x)|^{2}.$ (29)

Likewisely, to understand the spin current $J^{z}$ at the left edge, we
calculate the spin polarization $\sigma^{z}$ at the left edge for each
eigenstate

$\displaystyle m_{n_{x},k_{y}}^{z,LE}=\sum_{1\leq x\leq
N_{x}/2}\langle\sigma^{z}_{n_{x},k_{y}}(x)\rangle.$ (30)

These quantities are represented by colors in the band structure. The results
are summarized below according to different magnetization directions.

(i) ${\bf S}={\bf 0}$ : In the absence of magnetization, the
$\left\\{J^{x},J^{y}\right\\}$ flowing along $\left\\{{\hat{\bf y}},{\hat{\bf
x}}\right\\}$ induced by RSOC are present in the zigzag ribbon, whereby
$J^{y}$ is largely suppressed due to OBC in ${\hat{\bf x}}$Sun and Xie (2007).
A helical (spin chiral) edge spin current $J^{z}$ is also produced, as hinted
by the results in 2DEGReynoso  _et al._ (2004); Grigoryan  _et al._ (2009);
Nakhmedov and Alekperov (2012), although in zigzag ribbons it exists even
without a magnetic field and demonstrates helicity. In Fig.2 the colored band
structure clarifies the origin of $J^{z}$: for every left edge (LE) spin up
$k_{y}$ state (blue) there exists a corresponding LE spin down $-k_{y}$ state
(green), yielding counter propagating spins at the LE. This is true for all
eigenstates, so the finite $J^{z}$ is not only the result of low energy
states, similar to that has been discussed recently for the QSHE in
topological insulatorsChen (2020b). These features for the spin currents
remain true even in the presence of a finite magnetization ${\bf S}\neq{\bf
0}$.

(ii) ${\bf S}\parallel{\hat{\bf x}}$ : It is known that in 2DEG with RSOC, an
in-plane magnetization pointing perpendicular to the edge produces an edge
charge current, whose flow direction depends on the distance away from the
edgeUsaj and Balseiro (2005). ${\bf S}\parallel{\hat{\bf x}}$ corresponds to
this situation, in which we indeed see a charge current $J^{0}$ that is
symmetric (nonchiral) between the two edges. The eigenstates are not
particularly localized at either edge (the band structure is not particularly
red or green), but the band structure becomes asymmetric
$E(n_{x},k_{y})=-E(n_{x},-k_{y})$, which can cause $J^{0}$ along the ribbon, a
mechanism that has been pointed out for a superconductor/noncollinear
magnetChen and Schnyder (2015) and TI/FMM heterostructuresZegarra  _et al._
(2020); Chen (2020a). A large transverse spin polarization $\sigma^{z}$ is
induced near the edge without any bias voltage (in contrast to that produced
by a bias voltageNomura  _et al._ (2005); Bokes and Horváth (2010); Erlingsson
_et al._ ; Zhang (2014); Khaetskii (2014); Khaetskii and Egues (2017); Cysne
_et al._ ), signaling the existence of a local torque. The torque is
antisymmetric between the two edges and hence averages to zero, and is
expected to cause a non-collinear order between the two edges as in the two-
spin model in Sec. III.1. For the parameters in Eq. (9), the largest edge
currents are of the order of $\sim 0.01$, in units of $et/\hbar\sim 10^{-4}$A
for $J_{y}^{0}$ and $\mu_{B}t/\hbar\sim 10^{15}\mu_{B}/$s for
$\left\\{J^{x},J^{y},J^{z}\right\\}$, where $\mu_{B}$ is the Bohr magneton,
and the largest edge spin polarizations are of the order of $\sim
0.01\mu_{B}$.

(iii) ${\bf S}\parallel{\hat{\bf y}}$ : For the case of a magnetization along
the ribbon, we find no charge current $J^{0}$ and no transverse spin
polarization $\sigma^{x}=\sigma^{z}=0$, and hence there is no local torque.
The energy spectrum is half-metallic.

(iv) ${\bf S}\parallel{\hat{\bf z}}$ : Interestingly, we find that an out-of-
plane magnetization also produces a charge current $J^{0}$, but it is
antisymmetric (chiral) between the two edges. Comparing with that in the ${\bf
S}\parallel{\hat{\bf x}}$ case, this suggests the chirality of $J^{0}$ can be
controlled by the orientation of the magnetization. This current only occurs
when the chemical potential is finite $\mu\neq 0$, similar to that in the QAHE
in Chern insulatorsChen (2020b), although our spectrum remains gapped. The
band structure is symmetric between $+k_{y}$ and $-k_{y}$ at ${\bf
S}\parallel{\hat{\bf z}}$, but the wave function distribution is not: If
$+k_{y}$ state is mostly localized at the left edge (red), then the $-k_{y}$
state is more localized at the right edge (green), suggesting counter
propagating charge currents at the two edges. Note that some of the low energy
states are inherited from the flat band edge states of the pristine zigzag
ribbon, which become dispersive and chiral under the influence of
magnetization and RSOC. Finally, a transverse spin polarization $\sigma^{x}$
occurs near the edge, which are opposite at the two edges and hence average to
zero.

We conclude that for magnetization along principle directions ${\bf
S}\parallel\left\\{{\hat{\bf x}},{\hat{\bf y}},{\hat{\bf z}}\right\\}$, only
$\left\\{\chi^{xz},\chi^{zx}\right\\}$ of the transverse susceptibility are
nonzero, which then yield a local spin torque that requires no bias voltage
according to Eq. (26). However, the transverse spin polarization is alway
antisymmetric between the two edges and hence integrates to zero
$\sum_{i}d{\bf S}_{i}/dt=0$, indicating no net torque on a macroscopic scale.
These features are similar to those in the two-site model in Sec. III.1.

Figure 3: Same quantities as in Fig. 2 but for glide plane symmetric armchair
ribbon. The pattern is symmetric (labeled SY) or antisymmetric (labeled AS) is
defined with respect to $P_{\rm arm}^{\rm gli}$ in Eq. (31). The largest arrow
and disk in the figure are of the magnitude $0.024$ and $0.011$, respectively.

###  III.3 Persistent currents and spin torques in armchair ribbons

The armchair ribbons are simulated by imposing PBC in ${\hat{\bf x}}$ and OBC
in ${\hat{\bf y}}$ to Eq. (2). The symmetry between the two edges
distinguishes mirror symmetric and the glide plane symmetric armchair ribbons,
with the corresponding symmetry operators

$\displaystyle P_{\rm arm}^{\rm mir}=\left\\{P_{y}|\,{\bf
0}\right\\},\;\;\;P_{\rm arm}^{\rm gli}=\left\\{P_{y}|\,3/2{\hat{\bf
x}}\right\\},$ (31)

where $P_{y}$ denotes the mirror operation along ${\hat{\bf y}}$ with respect
to the central axis of the ribbon. The symmetry properties of the currents and
spin polarization for the mirror and glide plane symmetric armchair ribbons,
defined under $P_{\rm arm}^{\rm mir}$ and $P_{\rm arm}^{\rm gli}$,
respectively, are identical, so we only present the glide plane symmetric case
in Fig. 3. The edge spin and charge operators are analogous to those
introduced in Eq. (29) and (30), with an exchange of coordinates
$x\leftrightarrow y$. The results for different magnetization directions are
summarized below.

(i) ${\bf S}={\bf 0}$ : The case $S=0$ has both the $J^{x}$ and $J^{y}$ caused
by RSOC. However, $J^{x}$ is strongly suppressed due to OBC in ${\hat{\bf y}}$
direction. The out-of-plane polarized $J^{z}$ at the edge is clearly visible,
and also originates from the counter propagating spins as can be deduced from
the color codes of the band structure. These features of spin currents remain
valid in the ${\bf S}\neq{\bf 0}$ cases below.

(ii) ${\bf S}\parallel{\hat{\bf x}}$ : For the magnetization pointing along
the ribbon, there is no charge current $J^{0}$, and a spin polarization in
both transverse directions is induced. The $\sigma^{y}$ component is
antisymmetric and, hence, constitutes a non-collinear magnetic order between
the two edges. The $\sigma^{z}$ component is symmetric for the two edges, but
has an alternating signs along the ribbon, causing non-collinear magnetic
order along the ribbon.

(iii) ${\bf S}\parallel{\hat{\bf y}}$ : For magnetization in-plane but
perpendicular to the ribbon direction, an asymmetric band structure causes a
nonchiral $J^{0}$, and the spin polarization in both transverse directions
$\sigma^{x}$ and $\sigma^{z}$ are induced and are asymmetric between the two
edges.

(iv) ${\bf S}\parallel{\hat{\bf z}}$ : The out-of-plane magnetization again
causes a symmetric band structure. For $+k_{x}$ state that is more localized
at the right edge, we find a corresponding $-k_{x}$ state at the left edge,
which suggests the existence of a chiral $J^{0}$. Interestingly, this feature
is not obvious at low energy, but more prominent for the higher energy states.
Comparing this with the ${\bf S}\parallel{\hat{\bf y}}$ result, we again see
that the chirality of $J^{0}$ can be controlled by the orientation of the
magnetization. The $\sigma^{z}$ component is antisymmetric between the two
edges, while $\sigma^{y}$ is symmetric but with sign-alternating along the
ribbon. Finally, comparing the spin polarization $\sigma^{\alpha}$ in all
three situations ${\bf S}\parallel\left\\{{\hat{\bf x}},{\hat{\bf
y}},{\hat{\bf z}}\right\\}$, we conclude that every component of the local
susceptibility tensor $\chi^{\alpha\beta}$ is nonzero, different from that of
the zigzag ribbon in Sec. III.2. Nevertheless, the total torque vanishes in
all cases.

###  III.4 Partially magnetized nanoribbons and irregular nanoflakes

To make the magnetoelectric torque observable, one must overcome the overall
cancelling of the torque as seen in our discussion of the nanoribbons. We
suggest two situations that the net torque can be nonzero. The first is to
make the magnetization spatially varying, e.g. only occupy the region closer
to one edge, such that the other edge is idle. Using zigzag ribbon with ${\bf
S}\parallel{\hat{\bf x}}$ as an example, in Fig. 4 we indeed see a nonzero net
transverse spin polarization in this situation, whose average value per site
is of the order of $\overline{\sigma^{z}}\sim 10^{-4}$ when the magnetization
only covers a region near the left edge, and it is dramatically enhanced to
$\overline{\sigma^{z}}\sim 10^{-2}$ if the magnetization only exists on the
edge sites. Following Eq. (26), the $\overline{\sigma^{z}}\sim 10^{-4}$ and
$J_{ex}\sim 0.1$eV would yield a very large spin torque $d{\bf S}/dt\sim
10$GHz. Although we likely overestimated this torque due to the large $J_{ex}$
and $\lambda_{R}$ in Eq. (9), even if the torque is reduced by two orders of
magnitude to $\sim 0.1$GHz, it is still significant.

The second proposal is to use graphene nanoflakes with an irregular shapes,
which may help to generate nonuniform distribution of spin polarization that
does not sum to zeroMa and Sheng (2011); Weymann  _et al._ (2012); Güçlü  _et
al._ (2013); Luo and Sheng (2014); Szałowski (2015), given that the two edges
are not equivalent or it may even be ambiguous to identify two opposite edges.
In the example shown in Fig. 5, we see that the L-shape nanoflake indeed
yields a net $\overline{\sigma^{y}}$ and $\overline{\sigma^{z}}$ per site of
the order of $10^{-3}$ tp $10^{-4}$ when magnetization points at ${\bf
S}\parallel{\hat{\bf x}}$, indicating a net torque of significant strength. In
addition, despite the open boundary in all directions, still equilibrium
currents $\left\\{J^{0},J^{x},J^{y},J^{z}\right\\}$ exist. These currents turn
into networks of local currents, suggesting that they survive even in
realistic experimental situations of small and open boundary nanoflakes.

Figure 4: The unit cell of partially magnetized zigzag ribbon at magnetization
${\bf S}\parallel{\hat{\bf x}}$, where the transverse spin polarization
$\sigma^{z}$ is indicated by the size of disks, with red positive and blue
negative. The region to the left of the dotted line is magnetized, and to the
right is unmagnetized, and hence the panels from top to bottom represent the
magnetization covering $1/2$, $1/3$, $1/6$, and the edge sites of the ribbon.
The average spin polarization per site in the magnetized region is indicated
by $\overline{\sigma^{z}}$.  Figure 5: An L-shape nanoflake with magnetization
${\bf S}\parallel{\hat{\bf x}}$, where we show the transverse spin
polarizations $\sigma^{y}$ and $\sigma^{z}$, charge current $J^{0}$, and out-
of-plane polarized spin current $J^{z}$. The $\overline{\sigma^{y,z}}$
represents the average spin polarization per site.

##  IV Conclusions

In summary, we demonstrate that graphene nanoribbons due to RSOC and geometric
confinement display a peculiar magnetic response, including a bias voltage-
free spin torque, chiral and nonchiral edge currents, and equilibrium spin
currents. Using a two-site toy model, we could analytically show features such
as a transverse susceptibility caused by geometric confinement and linear in
the RSOC strength. Numerical calculation suggests that the same features occur
in nanoribbons. Given the controllability of RSOC by gate voltageYang  _et
al._ (2016); Wang  _et al._ (2016); Yang  _et al._ (2017); Safeer  _et al._
(2019); Ghiasi  _et al._ (2019); Benítez  _et al._ (2020), the transverse
susceptibility can serve as a gate-voltage-induced magnetoelectric torque
distinct from the usual current-induced spin-orbit torque. In a zigzag ribbon,
this torque tends to create noncollinear spin polarization between the two
edges, and in armchair ribbon the noncollinear order is not only between the
edges but can also be along the ribbon. Although the net torque sums to zero
in uniformly magnetized ribbons, it can be nonzero if the magnetization only
covers parts of the ribbon, say, closer to one edge, or in nanoflakes of
irregular shapes. We estimate a net torque that can reach sub-GHz magnitude,
pointing to the possibility of practical applications.

We further confirm that the equilibrium edge charge current $J^{0}$ predicted
for 2DEG, when the magnetization points perpendicular to the edge, also exists
in Rashba nanoribbons. Moreover, we discovered that $J^{0}$ can be nonchiral
due to the asymmetric band structure, or chiral due to counter propagating
eigenstates localized at opposite edges, depending on the magnetization
pointing in-plane or out-of-plane. In addition, the in-plane polarized
persistent spin currents $J^{x}$ and $J^{y}$ that exist ubiquitously in 2D
Rashba systems also manifest in nanoribbons, but the component flowing in the
confined direction is strongly suppressed. Besides, the RSOC causes an out-of-
plane polarized chiral edge spin current $J^{z}$ even in the absence of
magnetization. This discovery poses a challenge to distinguish $J^{z}$ from
that of purely topological origin in graphene-based topological insulators
containing RSOC, such as the Kane-Mele modelKane and Mele (2005a, b).
Generally, the energy gap does not affect the existence of these currents and
spin torques, since they are contributed from all the eigenstates in the Fermi
sea, not only the low energy flat band edge states that in some cases become
chiral. In addition, all these currents survive even in nanoflakes that have
open boundary in every direction. We anticipate that the controllability of
RSOC by gate voltage and the magnetization by magnetic field offers a
practical way to engineer these effects, which may help to realize them for
practical purposes, such as building graphene-based spintronic devicesPesin
and MacDonald (2012b); Han  _et al._ (2014); Roche  _et al._ (2015).

##  V Acknowledgement

We thank exclusively J. C. Egues and A. Zegarra for fruitful discussions. This
work is financially supported by the productivity in research fellowship from
CNPq.

##  Appendix A Explicit form of the charge and spin current operators

To calculate the commutators in Eq. (10), we observe that the hopping and the
Rashba Hamiltonians take the following general form

$\displaystyle H_{t}+H_{R}=\sum_{i\in
A}\sum_{\eta}\sum_{\alpha\beta}\left\\{T_{\alpha\beta}^{\eta}c_{i\alpha}^{{\dagger}}c_{i+\eta\beta}+T_{\alpha\beta}^{\eta\ast}c_{i+\eta\beta}^{{\dagger}}c_{i\alpha}\right\\},$

where $T_{\alpha\beta}^{\eta}$ is the hopping amplitude that the spin $\beta$
at site $i+\eta$ becomes spin $\alpha$ at site $i$, and
$T_{\alpha\beta}^{\eta\ast}$ is the complex conjugate of this amplitude. The
results of the commutators give

$\displaystyle J_{i,i+\eta}^{0}$ $\displaystyle=$
$\displaystyle\frac{ia}{\hbar}\sum_{\alpha\beta}\left\\{T_{\alpha\beta}^{\eta}c_{i\alpha}^{{\dagger}}c_{i+\eta\beta}-T_{\alpha\beta}^{\eta\ast}c_{i+\eta\beta}^{{\dagger}}c_{i\alpha}\right\\},$
$\displaystyle J_{i,i+\eta}^{a}$ $\displaystyle=$
$\displaystyle\frac{ia}{\hbar}\sum_{\alpha\beta}\sum_{\mu\nu}\left\\{T_{\alpha\beta}^{\eta}\delta_{\alpha\nu}\sigma_{\mu\nu}^{a}c_{i\mu}^{{\dagger}}c_{i+\eta\beta}\right.$
(33)
$\displaystyle\left.-T_{\alpha\beta}^{\eta\ast}\delta_{\mu\alpha}\sigma_{\mu\nu}^{a}c_{i+\eta\beta}^{{\dagger}}c_{i\nu}\right\\}.$

For
${\boldsymbol{\eta}}=\left(-{\boldsymbol{\delta}_{1}},-{\boldsymbol{\delta}_{2}},-{\boldsymbol{\delta}_{3}}\right)$,
the complete list of nonzero hopping amplitudes are

$\displaystyle
T_{\uparrow\uparrow}^{-\delta_{i}}=T_{\downarrow\downarrow}^{-\delta_{i}}=t\;\;\;{\rm
for}\;i=\left\\{1,2,3\right\\},$ $\displaystyle
T_{\uparrow\downarrow}^{-\delta_{1}}=i\lambda_{R}\left(-\frac{\sqrt{3}}{2}-\frac{i}{2}\right),\;\;\;T_{\downarrow\uparrow}^{-\delta_{1}}=i\lambda_{R}\left(-\frac{\sqrt{3}}{2}+\frac{i}{2}\right),$
$\displaystyle
T_{\uparrow\downarrow}^{-\delta_{2}}=i\lambda_{R}\left(\frac{\sqrt{3}}{2}-\frac{i}{2}\right),\;\;\;T_{\downarrow\uparrow}^{-\delta_{2}}=i\lambda_{R}\left(\frac{\sqrt{3}}{2}+\frac{i}{2}\right),$
$\displaystyle
T_{\uparrow\downarrow}^{-\delta_{3}}=-\lambda_{R},\;\;\;T_{\downarrow\uparrow}^{-\delta_{3}}=\lambda_{R},$
(34)

and $T_{\alpha\beta}^{\delta_{i}}=-T_{\alpha\beta}^{-\delta_{i}}$ for
${\boldsymbol{\eta}}=\left({\boldsymbol{\delta}_{1}},{\boldsymbol{\delta}_{2}},{\boldsymbol{\delta}_{3}}\right)$.
The explicit forms of the current operators are then given by

$\displaystyle
J_{i,i+\eta}^{0}=\frac{ia}{\hbar}\left\\{t\,c_{i\uparrow}^{{\dagger}}c_{i+\eta\uparrow}-t\,c_{i+\eta\uparrow}^{{\dagger}}c_{i\uparrow}+t\,c_{i\downarrow}^{{\dagger}}c_{i+\eta\downarrow}-t\,c_{i+\eta\downarrow}^{{\dagger}}c_{i\downarrow}\right.$
$\displaystyle\left.+T_{\uparrow\downarrow}^{\eta}\,c_{i\uparrow}^{{\dagger}}c_{i+\eta\downarrow}-T_{\uparrow\downarrow}^{\eta\ast}\,c_{i+\eta\downarrow}^{{\dagger}}c_{i\uparrow}+T_{\downarrow\uparrow}^{\eta}\,c_{i\downarrow}^{{\dagger}}c_{i+\eta\uparrow}-T_{\downarrow\uparrow}^{\eta\ast}\,c_{i+\eta\uparrow}^{{\dagger}}c_{i\downarrow}\right\\},$
$\displaystyle
J_{i,i+\eta}^{x}=\frac{ia}{\hbar}\left\\{t\,c_{i\downarrow}^{{\dagger}}c_{i+\eta\uparrow}-t\,c_{i+\eta\uparrow}^{{\dagger}}c_{i\downarrow}+t\,c_{i\uparrow}^{{\dagger}}c_{i+\eta\downarrow}-t\,c_{i+\eta\downarrow}^{{\dagger}}c_{i\uparrow}\right.$
$\displaystyle\left.+T_{\uparrow\downarrow}^{\eta}\,c_{i\downarrow}^{{\dagger}}c_{i+\eta\downarrow}-T_{\uparrow\downarrow}^{\eta\ast}\,c_{i+\eta\downarrow}^{{\dagger}}c_{i\downarrow}+T_{\downarrow\uparrow}^{\eta}\,c_{i\uparrow}^{{\dagger}}c_{i+\eta\uparrow}-T_{\downarrow\uparrow}^{\eta\ast}\,c_{i+\eta\uparrow}^{{\dagger}}c_{i\uparrow}\right\\},$
$\displaystyle
J_{i,i+\eta}^{y}=\frac{a}{\hbar}\left\\{-t\,c_{i\downarrow}^{{\dagger}}c_{i+\eta\uparrow}-t\,c_{i+\eta\uparrow}^{{\dagger}}c_{i\downarrow}+t\,c_{i\uparrow}^{{\dagger}}c_{i+\eta\downarrow}+t\,c_{i+\eta\downarrow}^{{\dagger}}c_{i\uparrow}\right.$
$\displaystyle\left.-T_{\uparrow\downarrow}^{\eta}\,c_{i\downarrow}^{{\dagger}}c_{i+\eta\downarrow}-T_{\uparrow\downarrow}^{\eta\ast}\,c_{i+\eta\downarrow}^{{\dagger}}c_{i\downarrow}+T_{\downarrow\uparrow}^{\eta}\,c_{i\uparrow}^{{\dagger}}c_{i+\eta\uparrow}+T_{\downarrow\uparrow}^{\eta\ast}\,c_{i+\eta\uparrow}^{{\dagger}}c_{i\uparrow}\right\\},$
$\displaystyle
J_{i,i+\eta}^{z}=\frac{ia}{\hbar}\left\\{t\,c_{i\downarrow}^{{\dagger}}c_{i+\eta\uparrow}-t\,c_{i+\eta\uparrow}^{{\dagger}}c_{i\downarrow}-t\,c_{i\uparrow}^{{\dagger}}c_{i+\eta\downarrow}+t\,c_{i+\eta\downarrow}^{{\dagger}}c_{i\uparrow}\right.$
$\displaystyle\left.+T_{\uparrow\downarrow}^{\eta}\,c_{i\uparrow}^{{\dagger}}c_{i+\eta\downarrow}-T_{\uparrow\downarrow}^{\eta\ast}\,c_{i+\eta\downarrow}^{{\dagger}}c_{i\uparrow}-T_{\downarrow\uparrow}^{\eta}\,c_{i\downarrow}^{{\dagger}}c_{i+\eta\uparrow}+T_{\downarrow\uparrow}^{\eta\ast}\,c_{i+\eta\uparrow}^{{\dagger}}c_{i\downarrow}\right\\}.$

We then evaluate the equilibrium charge and spin current numerically by their
expectation values

$\displaystyle\langle J_{i,i+\eta}^{a}\rangle=\sum_{n}\langle
n|J_{i,i+\eta}^{a}|n\rangle\,f(E_{n}),$ (36)

where $|n\rangle$ is the eigenstate of the lattice Hamiltonian, and
$f(E_{n})=\left(e^{E_{n}/k_{B}T}+1\right)^{-1}$ is the Fermi distribution
function. Note that this strategy includes all the eigenstates in the Fermi
sea, in contrast to previous theories that consider only the low energy
statesZhang  _et al._ (2014); Li and Yao (2017). Often times we ignore that
bracket $\langle{\hat{\cal O}}\rangle\equiv{\hat{\cal O}}$ of the expectation
value for convenience.

## References

  * Bychkov and Rashba (1984a) Y. A. Bychkov and E. I. Rashba, JETP Lett. 39, 78 (1984a).
  * Bychkov and Rashba (1984b) Y. A. Bychkov and E. I. Rashba, J. Phys. C: Solid State Phys. 17, 6039 (1984b).
  * Edelstein (1990) V. M. Edelstein, Solid State Commun. 73, 233 (1990).
  * Manchon and Zhang (2008) A. Manchon and S. Zhang, Phys. Rev. B 78, 212405 (2008).
  * Manchon and Zhang (2009) A. Manchon and S. Zhang, Phys. Rev. B 79, 094422 (2009).
  * Haney and Stiles (2010) P. M. Haney and M. D. Stiles, Phys. Rev. Lett. 105, 126602 (2010).
  * Gambardella and Miron (2011) P. Gambardella and I. M. Miron, Phil. Trans. R. Soc. A 369, 3175 (2011).
  * Pesin and MacDonald (2012a) D. A. Pesin and A. H. MacDonald, Phys. Rev. B 86, 014416 (2012a).
  * Haney  _et al._ (2013) P. M. Haney, H.-W. Lee, K.-J. Lee, A. Manchon, and M. D. Stiles, Phys. Rev. B 88, 214417 (2013).
  * Miron  _et al._ (2010) I. M. Miron, G. Gaudin, S. Auffret, B. Rodmacq, A. Schuhl, S. Pizzini, J. Vogel, and P. Gambardella, Nat. Mater. 9, 230 (2010).
  * Miron  _et al._ (2011) I. M. Miron, K. Garello, G. Gaudin, P.-J. Zermatten, M. V. Costache, S. Auffret, S. Bandiera, B. Rodmacq, A. Schuhl, and P. Gambardella, Nature 476, 189 (2011).
  * Manchon  _et al._ (2019) A. Manchon, J. Železný, I. M. Miron, T. Jungwirth, J. Sinova, A. Thiaville, K. Garello, and P. Gambardella, Rev. Mod. Phys. 91, 035004 (2019).
  * Rashba (2003) E. I. Rashba, Phys. Rev. B 68, 241315 (2003).
  * Sonin (2007) E. B. Sonin, Phys. Rev. B 76, 033306 (2007).
  * Tokatly (2008) I. V. Tokatly, Phys. Rev. Lett. 101, 106601 (2008).
  * Reynoso  _et al._ (2004) A. Reynoso, G. Usaj, M. J. Sánchez, and C. A. Balseiro, Phys. Rev. B 70, 235344 (2004).
  * Grigoryan  _et al._ (2009) V. L. Grigoryan, A. Matos Abiague, and S. M. Badalyan, Phys. Rev. B 80, 165320 (2009).
  * Usaj and Balseiro (2005) G. Usaj and C. A. Balseiro, EPL 72, 631 (2005).
  * Zegarra  _et al._ (2020) A. Zegarra, J. C. Egues, and W. Chen, Phys. Rev. B 101, 224438 (2020).
  * Chen (2020a) W. Chen, Phys. Rev. B 102, 144442 (2020a).
  * Novoselov  _et al._ (2005) K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, M. I. Katsnelson, I. V. Grigorieva, S. V. Dubonos, and A. A. Firsov, Nature 438, 197 (2005).
  * Zhou  _et al._ (2006) S. Y. Zhou, G.-H. Gweon, J. Graf, A. V. Fedorov, C. D. Spataru, R. D. Diehl, Y. Kopelevich, D.-H. Lee, S. G. Louie, and A. Lanzara, Nat. Phys. 2, 595 (2006).
  * Geim and Novoselov (2007) A. K. Geim and K. S. Novoselov, Nat. Mater. 6, 183 (2007).
  * Castro Neto  _et al._ (2009) A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. Novoselov, and A. K. Geim, Rev. Mod. Phys. 81, 109 (2009).
  * Geim (2009) A. K. Geim, Science 324, 1530 (2009).
  * Das Sarma  _et al._ (2011) S. Das Sarma, S. Adam, E. H. Hwang, and E. Rossi, Rev. Mod. Phys. 83, 407 (2011).
  * Fujita  _et al._ (1996) M. Fujita, K. Wakabayashi, K. Nakada, and K. Kusakabe, Journal of the Physical Society of Japan 65, 1920 (1996).
  * Nakada  _et al._ (1996) K. Nakada, M. Fujita, G. Dresselhaus, and M. S. Dresselhaus, Phys. Rev. B 54, 17954 (1996).
  * Brey and Fertig (2006) L. Brey and H. A. Fertig, Phys. Rev. B 73, 235411 (2006).
  * Ezawa (2006) M. Ezawa, Phys. Rev. B 73, 045432 (2006).
  * Peres  _et al._ (2006) N. M. R. Peres, F. Guinea, and A. H. Castro Neto, Phys. Rev. B 73, 125411 (2006).
  * Akhmerov and Beenakker (2008) A. R. Akhmerov and C. W. J. Beenakker, Phys. Rev. B 77, 085423 (2008).
  * Tao  _et al._ (2011) C. Tao, L. Jiao, O. V. Yazyev, Y.-C. Chen, J. Feng, X. Zhang, R. B. Capaz, J. M. Tour, A. Zettl, S. G. Louie, H. Dai, and M. F. Crommie, Nat. Phys. 7, 616 (2011).
  * Ruffieux  _et al._ (2011) P. Ruffieux, S. Wang, B. Yang, C. Sánchez-Sánchez, J. Liu, T. Dienel, L. Talirz, P. Shinde, C. A. Pignedoli, D. Passerone, T. Dumslaff, X. Feng, K. Müllen, and R. Fasel, Nature 531, 489 (2011).
  * Son  _et al._ (2006a) Y.-W. Son, M. L. Cohen, and S. G. Louie, Phys. Rev. Lett. 97, 216803 (2006a).
  * Son  _et al._ (2006b) Y.-W. Son, M. L. Cohen, and S. G. Louie, Nature 444, 347 (2006b).
  * Barone  _et al._ (2006) V. Barone, O. Hod, and G. E. Scuseria, Nano Lett. 6, 2748 (2006).
  * Yang  _et al._ (2007) L. Yang, C.-H. Park, Y.-W. Son, M. L. Cohen, and S. G. Louie, Phys. Rev. Lett. 99, 186801 (2007).
  * Han  _et al._ (2007) M. Y. Han, B. Özyilmaz, Y. Zhang, and P. Kim, Phys. Rev. Lett. 98, 206805 (2007).
  * Ritter and Lyding (2009) K. A. Ritter and J. W. Lyding, Nat. Mater. 8, 235 (2009).
  * Wakabayashi  _et al._ (1999) K. Wakabayashi, M. Fujita, H. Ajiki, and M. Sigrist, Phys. Rev. B 59, 8271 (1999).
  * Lee  _et al._ (2005) H. Lee, Y.-W. Son, N. Park, S. Han, and J. Yu, Phys. Rev. B 72, 174431 (2005).
  * Pisani  _et al._ (2007) L. Pisani, J. A. Chan, B. Montanari, and N. M. Harrison, Phys. Rev. B 75, 064418 (2007).
  * Yang  _et al._ (2016) B. Yang, M.-F. Tu, J. Kim, Y. Wu, H. Wang, J. Alicea, R. Wu, M. Bockrath, and J. Shi, 2D Mater. 3, 031012 (2016).
  * Wang  _et al._ (2016) Z. Wang, D.-K. Ki, J. Y. Khoo, D. Mauro, H. Berger, L. S. Levitov, and A. F. Morpurgo, Phys. Rev. X 6, 041020 (2016).
  * Yang  _et al._ (2017) B. Yang, M. Lohmann, D. Barroso, I. Liao, Z. Lin, Y. Liu, L. Bartels, K. Watanabe, T. Taniguchi, and J. Shi, Phys. Rev. B 96, 041409 (2017).
  * Safeer  _et al._ (2019) C. K. Safeer, J. Ingla-Aynés, F. Herling, J. H. Garcia, M. Vila, N. Ontoso, M. R. Calvo, S. Roche, and F. Hueso, Luis E.and Casanova, Nano Lett. 19, 1074 (2019).
  * Ghiasi  _et al._ (2019) T. S. Ghiasi, A. A. Kaverzin, P. J. Blah, and B. J. van Wees, Nano Lett. 19, 5959 (2019).
  * Benítez  _et al._ (2020) L. A. Benítez, W. Savero Torres, J. F. Sierra, M. Timmermans, J. H. Garcia, S. Roche, M. V. Costache, and S. O. Valenzuela, Nat. Mater. 19, 170 (2020).
  * Mendes  _et al._ (2015) J. B. S. Mendes, O. Alves Santos, L. M. Meireles, R. G. Lacerda, L. H. Vilela-Leão, F. L. A. Machado, R. L. Rodríguez-Suárez, A. Azevedo, and S. M. Rezende, Phys. Rev. Lett. 115, 226601 (2015).
  * Wang  _et al._ (2015) Z. Wang, C. Tang, R. Sachs, Y. Barlas, and J. Shi, Phys. Rev. Lett. 114, 016603 (2015).
  * Dushenko  _et al._ (2016) S. Dushenko, H. Ago, K. Kawahara, T. Tsuda, S. Kuwabata, T. Takenobu, T. Shinjo, Y. Ando, and M. Shiraishi, Phys. Rev. Lett. 116, 166102 (2016).
  * Leutenantsmeyer  _et al._ (2016) J. C. Leutenantsmeyer, A. A. Kaverzin, M. Wojtaszek, and B. J. van Wees, 2D Mater. 4, 014001 (2016).
  * Rybkin  _et al._ (2018) A. G. Rybkin, A. A. Rybkina, M. M. Otrokov, O. Y. Vilkov, I. I. Klimovskikh, A. E. Petukhov, M. V. Filianina, V. Y. Voroshnin, I. P. Rusinov, A. Ernst, A. Arnau, E. V. Chulkov, and A. M. Shikin, Nano Lett. 18, 1564 (2018).
  * Diniz  _et al._ (2013) G. S. Diniz, M. R. Guassi, and F. Qu, J. Appl. Phys. 114, 243701 (2013).
  * Guassi  _et al._ (2015) M. R. Guassi, G. S. Diniz, N. Sandler, and F. Qu, Phys. Rev. B 92, 075426 (2015).
  * Dyrdał and Barnaś (2015) A. Dyrdał and J. Barnaś, Phys. Rev. B 92, 165404 (2015).
  * Li and Manchon (2016) H. Li and A. Manchon, Phys. Rev. B 93, 235317 (2016).
  * Rodriguez-Vega  _et al._ (2017) M. Rodriguez-Vega, G. Schwiete, J. Sinova, and E. Rossi, Phys. Rev. B 96, 235419 (2017).
  * Zollner  _et al._ (2020) K. Zollner, M. D. Petrović, K. Dolui, P. Plecháč, B. K. Nikolić, and J. Fabian, Phys. Rev. Research 2, 043057 (2020).
  * Rybkina  _et al._ (2020) A. A. Rybkina, A. G. Rybkin, I. I. Klimovskikh, P. N. Skirdkov, K. A. Zvezdin, A. K. Zvezdin, and A. M. Shikin, Nanotechnology 31, 165201 (2020).
  * Kane and Mele (2005a) C. L. Kane and E. J. Mele, Phys. Rev. Lett. 95, 226801 (2005a).
  * Kane and Mele (2005b) C. L. Kane and E. J. Mele, Phys. Rev. Lett. 95, 146802 (2005b).
  * Min  _et al._ (2006) H. Min, J. E. Hill, N. A. Sinitsyn, B. R. Sahu, L. Kleinman, and A. H. MacDonald, Phys. Rev. B 74, 165310 (2006).
  * Dyrdał and Barnaś (2017) A. Dyrdał and J. Barnaś, 2D Mater. 4, 034003 (2017).
  * Peralta  _et al._ (2019) M. Peralta, E. Medina, and F. Mireles, Phys. Rev. B 99, 195452 (2019).
  * Katsura  _et al._ (2005) H. Katsura, N. Nagaosa, and A. V. Balatsky, Phys. Rev. Lett. 95, 057205 (2005).
  * Bruno and Dugaev (2005) P. Bruno and V. K. Dugaev, Phys. Rev. B 72, 241302 (2005).
  * Chen  _et al._ (2014) W. Chen, P. Horsch, and D. Manske, Phys. Rev. B 89, 064427 (2014).
  * Dzyaloshinsky (1958) I. Dzyaloshinsky, J. Phys. Chem. Solids 4, 241 (1958).
  * Moriya (1960) T. Moriya, Phys. Rev. 120, 91 (1960).
  * Sun and Xie (2007) Q.-F. Sun and X. C. Xie, Int. J. Mod. Phys. B 21, 3687 (2007).
  * Nakhmedov and Alekperov (2012) E. Nakhmedov and O. Alekperov, Phys. Rev. B 85, 153302 (2012).
  * Chen (2020b) W. Chen, Phys. Rev. B 101, 195120 (2020b).
  * Chen and Schnyder (2015) W. Chen and A. P. Schnyder, Phys. Rev. B 92, 214502 (2015).
  * Nomura  _et al._ (2005) K. Nomura, J. Wunderlich, J. Sinova, B. Kaestner, A. H. MacDonald, and T. Jungwirth, Phys. Rev. B 72, 245330 (2005).
  * Bokes and Horváth (2010) P. Bokes and F. Horváth, Phys. Rev. B 81, 125302 (2010).
  * (78) S. I. Erlingsson, J. Carlos Egues, and D. Loss, Phys. Stat. Sol. C 3, 4317.
  * Zhang (2014) W. Zhang, Sci. Rep. 4, 6320 (2014).
  * Khaetskii (2014) A. Khaetskii, Phys. Rev. B 89, 195408 (2014).
  * Khaetskii and Egues (2017) A. Khaetskii and J. C. Egues, EPL 118, 57006 (2017).
  * (82) T. P. Cysne, F. S. M. Guimarães, L. M. Canonico, T. G. Rappoport, and R. B. Muniz, arXiv:2101.05049 .
  * Ma and Sheng (2011) Z. Ma and W. Sheng, Appl. Phys. Lett. 99, 083101 (2011).
  * Weymann  _et al._ (2012) I. Weymann, J. Barnaś, and S. Krompiewski, Phys. Rev. B 85, 205306 (2012).
  * Güçlü  _et al._ (2013) A. D. Güçlü, P. Potasz, and P. Hawrylak, “Graphene-based integrated electronic, photonic and spintronic circuit,” in  _Future Trends in Microelectronics_ (John Wiley and Sons, Ltd, 2013) pp. 308–318.
  * Luo and Sheng (2014) K. Luo and W. Sheng, J. Appl. Phys. 115, 053705 (2014).
  * Szałowski (2015) K. Szałowski, J. Magn. Magn. Mater. 382, 318 (2015).
  * Pesin and MacDonald (2012b) D. Pesin and A. H. MacDonald, Nat. Mater. 11, 409 (2012b).
  * Han  _et al._ (2014) W. Han, R. K. Kawakami, M. Gmitra, and J. Fabian, Nat. Nanotechnol. 9, 794 (2014).
  * Roche  _et al._ (2015) S. Roche, J. Åkerman, B. Beschoten, J.-C. Charlier, M. Chshiev, S. P. Dash, B. Dlubak, J. Fabian, A. Fert, M. Guimarães, F. Guinea, I. Grigorieva, C. SchÃ¶nenberger, P. Seneor, C. Stampfer, S. O. Valenzuela, X. Waintal, and B. van Wees, 2D Mater. 2, 030202 (2015).
  * Zhang  _et al._ (2014) H. Zhang, Z. Ma, and J.-F. Liu, Sci. Rep. 4, 6464 (2014).
  * Li and Yao (2017) Z.-X. Li and H. Yao, Phys. Rev. B 96, 241101 (2017).