Helical currents in metallic Rashba strips

We study the texture of helical currents in metallic planar strips in the presence of Rashba spin-orbit coupling (RSOC) on the lattice at zero temperature. In the noninteracting case, and in the absence of external electromagnetic sources, we determine by exact numerical diagonalization of the single-particle Hamiltonian, the distribution across the strip section of these Rashba helical currents (RHC) as well as their sign oscillation, as a function of the RSOC strength, strip width, electron filling, and strip boundary conditions. Then, we study the effects of charge currents introduced into the system by an Aharonov-Bohm flux for the case of rings or by a voltage bias in the case of open strips. The former setup is studied by variational Monte Carlo, and the later, by the time-dependent density-matrix-renormalization group technique. Particularly for strips formed by two, three and four coupled chains, we show how these RHC vary in the presence of such induced charge current, and how their differences between spin-up and spin-down electron currents on each chain, help to explain the distribution across the strip of charge currents, both of the spin conserving and the spin flipping types. We also predict the appearance of polarized charge currents on each chain. Finally, we show that these Rashba helical currents and their derived features remain in the presence of an on-site Hubbard repulsion as long as the system remains metallic, at quarter filling, and even at half-filling where a Mott-Hubbard metal-insulator transition occurs for large Hubbard repulsion.


I. INTRODUCTION
There is an ongoing effort to take advantage of electron spin, in addition to its charge, for computing and information technologies, that has evolved into the new field of spintronics. [1][2][3] An important ingredient for spintronics devices, is the spin-orbit (SO) interaction, which provides the essential mechanisms for spin polarization and spin currents without the need of magnetic materials. [4][5][6] In particular, at interfaces or surfaces of materials, the breaking of structure inversion symmetry, combined with the atomic spin-orbit interaction gives rise to the Rashba (or Bychkov-Rashba) spin-orbit coupling. 7,8 It is well-known that the Rashba SO coupling (Rashba SOC or RSOC) in two-dimensional (2D) systems leads to the presence of pure spin currents, which are transversal to the direction of an applied charge current. This spin current is considered the intrinsic mechanism for the spin-Hall effect, which in many systems coexists and competes with extrinsic mechanisms. 6,[9][10][11][12][13] For finite width 2D systems such as wires or strips, this spin current is the origin of the phenomenon of spin accumulation. 14,15 The Rashba SOC in noninteracting 2D systems has been widely studied in first quantization, 5,16 where the finite width of strips is enforced by an harmonic potential or by an infinite potential well. 17 Much less has been done using second quantization on lattices, which provides a road to study the crossover from one-dimensional (1D) to 2D behavior, [18][19][20] and to include realistic hoppings and SO couplings, as well as various types of electron correlations. These electron correlations are present in a number of surfaces and interfaces involving transition metal oxides, such as SrTiO 3 21-25 , particularly LaAlO 3 /SrTiO 3 , [26][27][28][29][30] and BaTiO 3 , for example in the heterostructure BaTiO 3 /BaOsO 3 . 31 In addition, second quantization and discrete lattices provide a scale of detail at a microscopic level, not accessible by continuum coarse grained techniques. The main physical property we want to address within this approach is the distribution of spin-up and spin-down electron currents that appear on strips due to the presence of edges, even in the absence of external electromagnetic potentials or fields. Helical edge currents where electrons with up and down spins move in opposite directions, such as those illustrated in Fig. 1, are well-known to be present in topological insulators of the quantum spin-Hall (QSH) type. [32][33][34][35][36][37] In the present work, we examine helical currents which appear in metallic planar strips in the presence of Rashba SOC. These helical currents are slowly decaying from the edges, and hence for a narrow strip, depending on the RSOC strength, they are present over all the section of the strip, not just on its edges as in the QSH case. This feature should be ultimately related to the fundamental difference between the Rashba SOC and the SO term that consists of two copies of the Haldane term and is the essence of the Kane-Mele model 36 This Haldane-like term of the Kane-Mele model, contrary to the Rashba SO term, conserves the z-component of the total spin, while breaking the timereversal symmetry. We emphasize that what we call in the following Rashba helical currents (RHC) (defined in section II), are charge currents for each spin projection, and should not be confused with the already mentioned transversal spin currents, which are the ones related to the spin-Hall conductivity 11,38 . In the absence of appplied magnetic or electric fields, at any point of the system, the spin-up Rashba helical current is exactly opposite to the spin-down one, and hence of course there is no net charge current.
One of the main features we will study is the texture of these Rashba induced helical currents, that is, the amplitude and sign oscillations of spin-up and spindown electron currents on different chains as a function of their distance to the edge. These sign oscillations of helical currents were not so far detected experimentally, and hence did not receive too much theoretical attention, since in semiconductor devices the electron densities are low and RSOC is weak. However, some traces of this texture have been reported theoretically on studies on semiconductor wires. 14,39 New materials which are being considered for spintronic devices can perform at higher electron densities and present much larger RSOC, and hence the results of the present study may be relevant for such systems. 31,40 Hence, the first purpose of the present article is to numerically determine the helical currents that appear due to the RSOC in two different systems, as shown in Fig. 1. We will first consider a ring geometry, that is, a strip with open boundary conditions (OBC) in the transversal direction and periodic boundary conditions (PBC) in the longitudinal one ( Fig. 1(a)). Then, we will analyze a strip with full OBC (Fig. 1(b)), which could be regarded as an idealized version of the wires employed in spintronic devices.
In Section III we start our study with noninteracting systems in the absence of external electromagnetic sources. By using a combination of momentum and realspace techniques on both kinds of strips, the Rashba helical current distribution across the transversal section will be determined.
Next, we examine how the Rashba helical currents determine the net charge currents that appear in the presence of external electromagnetic sources. This part of the work, performed in Section IV, although certainly is the most relevant to physical systems or devices, could only be properly understood by resorting to the behavior of the Rashba helical currents studied in the previous section. In the case of rings, shown in Fig. 1(a), persistent currents can be induced by applying an Aharonov-Bohm flux through it. For strips such as those of Fig. 1(b), these currents are induced by connecting the trips to two leads with a voltage difference. As shown in Fig. 1, both kinds of strips belong in principle to different topologies and transport physics, and the purpose of the present work is to examine how the Rashba helical currents give rise to net longitudinal currents in both geometries.
In this second part of our effort, we will also include the effects of electron correlations, introduced by a Hubbard term in the Hamiltonian, again relevant for a number of recently introduced materials and devices, as mentioned before. So far, few studies have been accomplished in models including the Rashba SOC and Hubbard repulsion, and most of them were done on quasionedimensional lattices. 19,20 We will show that the analysis of RHC gives a fresh approach on this kind of correlated models, particularly on its transport properties.

II. MODEL AND METHODS
The Rashba SO Hamiltonian on the square lattice is given by: 41,42 corresponds to the Hubbard model. The notation is standard. We will consider the isotropic case t x = t y = t, and V SO,x = V SO,y = V SO , where x (y) is the longitudinal (transversal) directions of the strips shown in Fig. 1, except otherwise stated. In the following, we have adopted the normalization, t 2 + V 2 SO = 1, which we adopt as the unit of energy.
It is well-known that with the Rashba term, the zprojection of the total spin is no longer conserved, as well as the rotational symmetry in the spin space, while the time-reversal symmetry is still conserved.
In this work, we will study planar strips, with L sites on the longitudinal direction, with periodic ( Fig. 1(a)) or open ( Fig. 1(b)), boundary conditions, and W sites in the transversal direction, with OBC. In general, we will take L ≫ W , so that the strips could be considered as W coupled chains. The total number of sites is N = L × W , and the total number of electrons, N e , such that the filling is n = N e /N . We will not consider the presence of impurities in our system.
The main quantities we are interested in are the currents associated to the hopping term in Eq. (2), termed the spin-conserving current, J σ,x,y,μ , on each link between (x, y) and (x, y) +μ, whereμ =x,ŷ, which is the expectation value of the operator: and the contribution from the RSOC term, which leads to the so-called spin-flipping current, J SO,x,y,μ which is the expectation value of a similar current operator that can be also similarly derived from Eq. (1) as the first order perturbation induced by a magnetic potential A added to Eq. (1) by the Peierls factors, γμ → γμ exp (−ieAμ) (in units of ), where γμ are the appropriate coupling constant in theμ-direction. In the following we will also adopt units where e = 1. Let us write down only the x-component of the SO current, which we will need for future discussions: The expectation value of each term in parenthesis will be denoted as J ′ SO,x,y,x and J ′′ SO,x,y,x , respectively, and similarly for the y-direction. The total current on each link is then J tot,x,y,μ = J ↑,x,y,μ + J ↓,x,y,μ + J SO,x,y,μ .
In the following we will compute and discuss the RHC as the "spin-up currents" ("spin-down currents"), defined as the hopping currents J ↑,x,y,x (J ↓,x,y,x ) described above.
Spin accumulation is a well-known consequence of the spin-Hall effect in the presence of boundaries in the direction transversal to the direction of the charge current. 14,15,19 That is, a charge current in the longitudinal direction that appears due to an external electromagnetic field, induces a transversal spin current, which then leads to spin polarization with opposite sign on opposite edges of the conducting strip. In principle, the spin accumulation is defined as the relative out-of-plane spin polarization between the two sides of the strip that is: where S z (j) is the sum of the z-component of the electron spins on the chain j. For wide strips it is expected that the spin accumulation is concentrated on the strip edges as it was observed experimentally 43 , although there is no consensus about the interpretation of these experimental results. In the narrow strips we will consider in Section IV, the spin polarization could take appreciable values across the whole section of the strip, so we will compute the polarization or spin accumulation on each chain, ∆S z (j) = S z (j) − S z (W − j + 1), j = 1, . . . , W/2. Hence the total spin acumulation can be expressed as ∆S z = 1≤j≤W/2 ∆S z (j). It has been noticed that spin polarization could also appear on strips due to the inverse spin galvanic effect (ISGE). 1 However, the current-induced spin polarization due to ISGE refers to in-plane (x, y) components of the spin, and has been studied in two-dimensional electron systems for varying x, y-directions of the currents. 44 In contrast, we would like to emphasize that the spin polarization we compute and discuss in Section IV, is the out-of-plane (z) component related to the Spin-Hall effect on strips, and hence not to the ISGE.
For the noninteracting case, we solve numerically the Hamiltonian for a single electron. For the case of rings, Fig. 1(a), translation symmetry is preserved along the longitudinal direction (x-axis), while it is broken in the transversal direction (y-axis) due to the open BC. In this case, for each k x one has to solve a 2W × 2W matrix, and the single particle eigenvalues and eigenvectors will be a function of (k x , y, σ). In the absence of external electromagnetic fields, static currents on each chain both of the hopping and SO types, are computed by straightforward ground state averaging. These zero temperature results are presented in Section III B.
For the geometry of Fig. 1(b), translation symmetry is broken also in the x direction, so the single particle Hamiltonian matrix is formulated in real space, and it has dimensions 2N × 2N . In this case, one has in principle currents of spin-σ electrons varying on each bond of the lattice, and in both directions, that is J σ,x,y,μ . In Section III C, only results for the currents at the center of the strip, J σ,L/2,y,μ will be discussed.
To deal with both external electromagnetic sources and nonzero electron-electron interactions, that is a finite value of U in the Hubbard term, we will resort to much more involved many-body techniques to study the zero temperature behavior. For the case of periodic strips, we will study persistent currents along the ring, which are induced by piercing the ring by an Aharonov-Bohm flux. 19,45,46 The study of these currents is accomplished by a simple variational Monte Carlo (VMC) where configurations generated by a trial wave function which consists of a product of Slater determinants, are weighted by a Gutzwiller factor which depends on the number of double-occupied sites and the value of U . 47,48 This technique allows us to study rings up to W = 4, and electronelectron Hubbard repulsion up to U = 8. Although this simple variational wave function underestimates the effect of U , the obtained results are qualitatively reliable as compared with those obtained by other techniques. Typically, for each set of parameters, and for each value of the Aharonov-Bohm flux, Φ, we take Monte Carlo averages over 10 6 sweeps.
In the case of an Aharonov-Bohm flux Φ piercing a ring, the total Hamiltonian H R + H H has to be treated in the presence of Peierls factors in both the hopping and the RSOC terms along the x-direction, and Ax should be replaced by Φ/L, where in the adopted units, 0 ≤ Φ ≤ 2π. In this case, the ground state eigenvalue and eigenvector will be E 0 = E 0 (Φ), Ψ 0 = Ψ 0 (Φ), and hence the ground state expectation value of any physical quantity will also depend on Φ.
To study currents induced by an external electric field, as well as electron correlations, in strips where OBC has been imposed also in the longitudinal direction, we resort to density matrix-renormalization group (DMRG) technique, 49 which is in principle exact for static properties. In this case, currents appear as a time dependent feature after applying a voltage bias between the ends of the strip at a given time, after the iterative process has converged to the ground state. This first standard DMRG was carried out keeping up 800 states, and performing at least 14 sweeps, enough to achieve truncation errors of the order O(10 −6 ). The voltage bias is applied between the two halves of the strip in the longitudinal direction. Alternatively, we have applied this bias to the two ends of the strip without affecting the two central columns of the strip. The results obtained with both recipes are virtually identical. The voltage bias is small, ∆V = 0.01, so as to avoid strong nonequilibrium processes. This calculation requires an extension of the original DMRG algorithm, and the results are considered "quasi-exact". Among different DMRG approaches to deal with time evolution problems 49-51 , we choose the algorithm introduced by Manmanna et al 52 , which is a variant based on a Krylov-space representation of the time evolution operator 53 . Thereby, we are able to evolve systems with interactions involving many coupled chains, and obtain very accurate results for long time periods. We observed that using at each time slice dt = 0.1, an Krylov base of ∼ 30 states, the errors become smaller than the symbol sizes. For each strip, and for a given set of parameters, a single value for each type of current was taken as the amplitude of the first maximum reached in its time evolution, considered as an approximation of the expected plateau.
The application of computational techniques, particularly DMRG, is severely limited by the non-conservation of the total S z , which implies a much larger Hilbert space for the system. The need to work with complex numbers, also implied by the Rashba SO coupling, further increases CPU time and memory requirements.

A. Rashba helical currents
Let us start with the noninteracting case, U = 0. As we will extensively show in the next sections, the RHC consist of spin-up or spin-down hopping currents which are exactly opposite on each chain across the strip section in the absence of external electromagnetic sources (fields or potentials). In addition, the currents for a given spin orientation, have opposite direction on chains symmetrically located with respect to the middle of the section, and their sign and strength vary over the section. In principle, counter propagating currents for each spin projection are guessed from general arguments similar to those leading to the prediction of the spin-Hall effect 11,37 and which ultimately stem from requirements of time-reversal symmetry.
An alternative and instructive way of understanding the presence of these currents with opposite circulation of spin up and down electrons in the absence of external sources, can be done in real space in systems with Rashba SOC. At an effective level, these Rashba helical currents appear due a term in the Hamiltonian, whereμ =x orŷ. The origin of this term is schematically shown in Fig. 2. A spin up electron at site "1" could reach site "4" by moving to site "2" involving a spin-orbit process, and then to site "3" with another SO process. Finally, it could reach site "4" by an ordinary hopping. The full process would involve σ y × σ x = −iσ z , with a factor λ proportional to tV 2 SO . An analogous process along sites "2' " and "3' " would have an opposite sign but this contribution is absent if that line is on a system edge. It is easy to see that the reverse process from "4" to "1" would lead to σ x × σ y = iσ z , thus providing the sign for the current in Eq. 6. In addition, processes involving more sites further away from the edge would contribute to that effective current on the edge, and also on legs located in the bulk due to the lack of compensation. If sites "1" and "4" belong to a leg located at y, processes through sites "2" and "3", located at "y+a" and process through sites "2' " and "3' ", located at "y-a" would not cancel each other because legs "y+a" and "y-a" are in general not symmetric due to the proximity to the edge. For example, an effective current involving five sites would lead to a current between sites at the line next to the edge, and so forth. Of course, each hopping would on average have a probability of (1 − n/2), where n is the electron filling, and assuming an equal number of electrons with up and down spins. Notice that for the so-called persistent helix state, 54 where the Rashba and Dresselhaus, 55 couplings are present with equal strength, the SO coupling is mediated by σ x −σ y in both directions, and hence the effective terms illustrated in Fig. 2 turns out to be of the kinetic type, not of the current type.
Of course, for very large strip widths, processes with longer paths become unlikely, and far from the edges, processes with shorter path would cancel because legs become approximately equivalent. It is expected then that for large strip widths, these currents become essentially nonzero only at the strip edges. It will be tempting to assign the changes in sign and strength as a function of the distance of the chains to the edge, to a sort of Friedel oscillations, as it will be explored below.
It should also be noticed that the Rashba Hamiltonian Eq. (1) can be written in terms of currents, and the same kind of currents here studied could be the one seen in Ref. 56 in the presence of impurities. That is, the impurities provide the breaking of translational symmetry that appear in our case due to the open boundaries of strip geometry.
Let us remark that according to the previous arguments, the RHC are of the hopping or spin-conserving type, given by Eq. (3), not of the SO or spin-flipping type, corresponding to Eq. (4).
In Fig. 3 the evolution of energy bands as a function of k x for various values of V SO /t is shown for rings. Results were obtained for the 64×2000 ring. The top panel, corresponding to V SO /t ≫ 1, resembles the typical situation found in topological Dirac semimetals 57-59 , with its characteristic Dirac nodes in the Fermi surface. 60,61 However, the purpose of the present work is to study the evolution of Rashba helical currents when V SO /t < 1, when the electron density is reduced below half-filling, and hence the physics of the topological Dirac semimetals is not relevant. As V SO /t is reduced, a more conventional band metallic behavior starts to appear. At V SO /t = 0.4, a single conduction band appears. It is interesting to note that the behavior shown in Fig. 3 is also similar to the one found in the ferromagnetic Kondo lattice model with Rashba SOC in 2D, with classical localized spins, for large values of the Hund coupling at quarter filling. 62 Of course, this band structure is different to the one characteristic of QSH systems, reflecting the insulating character of the bulk, and the metallic character of the edge. 36,63 The density of states for the same 64 × 2000 ring in the noninteracting case is shown in Fig. 4, for various values of V SO /t. It can be clearly seen a suppression of energy states around zero energy, corresponding to half-filling, as V SO /t is increased from 0. In this Subsection, we study planar rings such as the one depicted in Fig. 1(a). Since as said earlier, Rashba helical currents of spin-down electrons are exactly opposite to the spin-up ones, in the following we will only consider the spin-up electron currents on each chain. Moreover, since currents are spatially antisymmetric on chains with respect to the central longitudinal axis, we will only show the currents on the chains in half the section of the strip. Finally, due to translation invariance along the longitudinal direction, the currents will not depend on x. Hence, in all figures below, we will depict J ↑,l , where l = 1, . . . , W/2 or J ↑,ν , where ν is the normalized distance of the chain to the edge, or chain depth, defined in such a way that ν = 0 at the edge and ν = 1 at the center, independently of the width W .
As discussed in Subsection III A, the RHC, computed here as J ↑,l , is a spin-conserving current, and spinflipping currents are zero everywhere, which is remark-  able since the SO coupling is needed to be nonzero in both directions for the RHC to exist. Moreover, notice that the conservation of the average charge on each site in the stationary regime (Kirchhoff's law) is trivially satisfied for up and down spin electrons separately. This is compatible not only with the net SO currents being zero but also with the two contributions J ′ SO,x,y,μ and J ′′ SO,x,y,μ , possibly being zero in both directions,μ =x,ŷ.
A direct evidence of the existence of the RHC and their predicted dependence with the square of V SO /t can be observed for strips of W = 2, 3, 4 and 8 at n = 0.5 at small V SO /t, as shown in Fig. 5. This dependence becomes increasingly complex as the strip width is in- creased. For all W , there is a clear change of behavior around V SO /t ≈ 0.5, and for W > 2, this behavior presents some cusps, which are typical of level crossings. In this, and in the following calculations in this Subsection, results for L = 800 are virtuable indiscernible from those for L = 8000, hence we are confident that the results are representative of the continuum limit in the longitudinal direction. As suggested by the arguments of Subsection III A, the RHC should concentrate near the strip edges, and only due to its slow decay from the edges, for small widths, the RHC appears to pervade the whole system. To examine this behavior we will consider a broad range of strip widths. As it can be seen in Fig. 6, the module of J ↑ as a function of the chain's depth ν indicates that only for W ≥ 16 and V SO /t ≥ 0.4 it is apparent that these currents become essentially different from zero near to the strip edges. In general, wider strips have to be considered to show an edge-like behavior for smaller SO couplings. Two more features become apparent. First, consistently with the predictions of Subsection III A, and with the behavior already shown in Fig. 5, the amplitude of these currents increases almost quadratically with V SO /t, and second, and more interesting, J ↑,ν presents an oscillatory behavior as a function of the chain depth, reminiscent of a Friedel oscillation behavior.
To quantify the decay of RHC as the chains are located away from the edge, we have fitted the curves for the widest ring, with W = 64, with an exponential law, Results from these fits for V SO /t = 0.1 are shown in Fig. 6(a), and a systematic study of this fitting will be presented below in Fig. 9. In Fig. 7, it can be noticed that the amplitude oscillation observed in the previous Figure, actually comes from a sign oscillation. In fact, close to the edge, and more clearly for larger widths, it seems that the oscillations would have a wavenumber π, for these results at quarter filling (n = 0.5).
Let us now consider the case of half-filling (n = 1). The Rashba helical currents on each chain are shown in Fig. 8 for the same values of W and V SO /t as before, now including also V SO /t = 0.8. Similarly to the behavior observed for n = 0.5, the RHC tend to concentrate near the strip edges as W and V SO /t increase. As a difference with the case of n = 0.5, at n = 1 the RHC is maximal at the edges, smoothly decaying towards the strip center. The absence of sign and amplitude oscillations for this density n = 1, is, as in the n = 0.5 case, compatible with a wavenumber 2nπ (n = 1 at half-filling), which would correspond to 4k F in a quasi-onedimensional system.
In Fig. 9 the coefficients B of the above mentioned exponential and power law fits of the RHC for (a) n = 0.5 and (b) n = 1, are shown. The results in (a) were obtained from the values of J ↑,ν shown in Fig. 6, and the results in (b) were obtained form the corresponding values in Fig. 8. In both cases, the largest ring width, W = 64 was considered. Although these results have somewhat large error bars, estimated by adopting different subsets of points on each curve, it is clear in both cases that B grows as a function of V SO,x /t, indicating a faster decay from the edges. For n = 0.5, and the values V SO,x /t considered, the quality of the fittings using both laws, as measured by the correlation coefficient, is similar, although in principle a power law behavior would seem more natural for a gapless system. On the other hand, for n = 1, the exponential fit is definitely poor, and only the power law fit leads to minimally acceptable correlation coefficient.  Let us emphasize that the behavior discussed in this Section occurs in the metallic phase. In Section IV B, the effect of the Hubbard term and a possible Mott-Hubbard type of metal-insulator transition (MIT) will be considered. However, before leaving this Subsection, we can explore the effect of a MIT due to the opening of a band gap which may occur for sufficiently large ratio t y /t x , and at specific electron fillings. To this end, we should consider an spatially anisotropic extension of the Hamiltonian, where we keep the ratio V SO,x /t x = V SO,y /t y = 1, and we do not use the normalization defined in Section II, but we used t x = 1 instead. In Fig. 10, a change of behavior in the Rashba helical currents is observed when a band gap opens in strips of widths W = 2 and 4, at n = 0.5. It is apparent that in the gaped region the RHC decrease approximately following a linear dependence, but clearly a more systematic research with more realistic band structures should be performed.

C. Results on open strips
We now turn to the geometry depicted in Fig. 1(b), with open boundary conditions in both directions. In this case, as shown in that Figure, the RHC would follow closed circuits. However, the currents on vertical segments would also present oscillatory behaviors similar to those shown previously for the currents on rings in the longitudinal direction. The result of the interference between horizontal and vertical oscillations leads to complex patterns, as the one illustrated in Fig. 11, that are little resembling the closed loops illustrated in Fig. 1(b), although some segments of closed loops are still visible near the edges.
One important difference between the full open strips considered in this Subsection, and rings, studied in the preceding one, is that spin-up and spin-down electron densities at each site, are no longer conserved separately, that is, the sum over the links connected to a given site of the hopping currents for each spin is no longer zero. This can be checked for example with the values of the currents in the zoom of Fig. 11. However, the Kirchhoff law for the total electron density at each site is still satisfied since on each link, J ↑,x,y,μ = −J ↓,x,y,μ , and the net SO current on each link is zero. It is interesting to notice, that in order to conserve the density of up electrons on each site, the currents J ′ SO,x,y,x defined after Eq. (4) must be nonzero in such a way that its sum over the links connected to a site must be equal to the corresponding sum of spin-up hopping currents. A similar conclusion can be reached for the sum of the spin-down currents and the other term of the SO currents, J ′′ SO,x,y,x . From the above discussion, to study the dependence of RHC on chains as a function of their distance to the edge, we have considered for simplicity the currents on the central link of these chains, that is J ↑,L/2,ν , where we have always adopted L even. To reduce finite size effects, for each width we have taken the largest value of L accessible to our computational facilities. to that found previously for rings for the same sets of parameters, that is, RHC is concentrated close to the edge for large V SO /t and W , while for small RSOC strength and widths, the overlapping of the tails of the decay from teh edges leads to RHC finite over the whole strip section. Results for W = 64 are probably affected by finite size effects on the length L.
Similarly also to the case of rings, the oscillations actually stem from sign oscillations, as it can observed in Figs. 12(c),(d), for the same parameters as in Figs. 12(a),(b), respectively. Again, it seems that this oscillation has a wavenumber of π at this density. Overall, the results for J ↑,L/2,ν seem more erratic than the results for J ↑,ν considered before for rings: although the later is a global quantity, the former is a local one.
In addition, a similar dependence of the results with the electron filling can be inferred by observing the results depicted in Fig. 13, corresponding to half-filling. Effects of V SO /t and W are shown in Fig. 13(a)-(c), with an overall behavior similar to that shown before for rings, Fig. 8, that is, an apparent absence of oscillations except for the 128 × 64, probably due to finite size effects. A somewhat more systematic study of the oscillations dependence on the electron filling is presented in Fig. 13(d), where results for V SO /t = 0.2, on the 256 × 32 strip, for the fillings n = 0.5, 0.75 and 1. For n = 0.75, it seems that the oscillation has period four, consistent with a 4k F value of π/2 in the first Brillouin zone.
To close this Subsection, let us emphasize that although it is not possible for open strips to reach the continuum limit as it was chieved in rings, the fact that one obtains essentially the same behavior indicates that finite size effects, and in lattice systems, the consequent discretized energy spectrum, are neligible.

IV. INTERACTING CASE: CURRENTS DUE TO EXTERNAL ELECTROMAGNETIC SOURCES
In this section, we will not only address the evolution of net carrier currents out of helical currents in Rashba metallic strips when external electromagnetic fields are applied, but also we will examine how electron correlations in the form of a Hubbard repulsion could affect those Rashba helical currents.

A. Persistent currents
We will consider first how the RHC in the presence of persistent currents in rings will lead to net charge currents on each chain, both of the spin-conserving and the spin-flipping types, and to polarized currents, as well as spin accumulation. Persistent currents induced by an Aharonov-Bohm (AB) flux, are studied in this section by VMC, up to U = 8. The value of each physical quantity obtained for a given V SO /t and U , density, and lattice size, is adopted as the average over the AB flux between 0 and π. Of course other criterions could be adopted, for example to take the maximum of the absolute value of each quantity on that interval of the AB flux. The applied AB flux installs a net persistent hopping current, J hop /2 = l J ↑,l,x = l J ↓,l,x = 0, where the sum extends over all the chains of the strip. Since there is a net charge current along the longitudinal direction, it is well-known that a torque will appear on the spins making them rotate in the (x, z) plane 37 , that is, there will be a net SO or spin-flipping current along the xdirection, J SO . Then, the total current in the x-direction is J tot = J hop + J SO Let us start by discussing the case of three-chain rings. In Fig. 14, the evolution of RHC currents on each chain, J ↑,y , with the AB flux, Φ is shown for V SO /t = 0.6, on the 16 × 3 ring. In this, and in similar plots, we display properties as a function of the AB flux in the range 0 ≤ Φ ≤ π, taking advantage of the symmetry where O is any physical quantity. The relation involving the outer chains, J ↑,1 = −J ↑,3 , is broken as soon as Φ is different from zero, as illustrated in Fig. 14(a), for n = 0.5 and U = 0. However, the conditions J ↓,1 = J ↑,3 and J ↑,1 = J ↓,3 , still hold for any Φ. On the inner chain, it is observed, for n = 0.5, (Fig. 14(b)) that J ↓,2 = J ↑,2 for all Φ. The same relations are present at half-filling, n = 1, as shown in Fig. 14(c) for the outer chains, and in Fig. 14(d) for the inner chain. The dependence of J ↑,y with the Hubbard repulsion U is also included in this plot. One should remember that in general VMC underestimates the effect of U . Nonetheless, it is quite apparent that the effects of U at n = 1 are stronger than at n = 0.5. In fact, as we will see in the next Subsection, where the DMRG technique exactly treats the electron correlations, the model considered experiences a Mott-Hubbard type of metal-insulating transition, given by the vanishing of the total currents. The results of a systematic study as a function of both V SO /t and U , are depicted in Fig. 15. In Fig. 15(a) it can be observed that the flux-averaged total currents decrease as a function of both V SO /t and U . For a fixed V SO /t, this decreasing of J tot as a function of U is much stronger for n = 1 than for n = 0.5, as expected on general grounds for the Hubbard repulsion. For both densities, and for any fixed value of U , the decrease of J tot with V SO /t is in general more important than the one due to U . In particular at U = 0, this result for W = 3 seems in contradiction with some suggestion in the literature that the longitudinal conductivity is independent of V SO /t and density. 65 Further numerical support for the observed variation of J tot with V SO /t and density will be given in the next Subsection, and further discussion will be provided in the final Section. Results for the flux-averaged spin accumulation, is shown in Fig. 15(b). There is a general trend that ∆S z reaches its maximum between V SO /t = 0.4 and 0.6 (the result for V SO /t = 1.0 is a finite size effect, since it is expected that ∆S z vanishes as V SO /t → ∞), and more importantly, there is a general trend towards an enhancement of the spin accumulation with the Hubbard coupling U . Both features were already observed for the W = 2 strip 19 , and together with the results provided below for the W = 4 strip, suggest that they could be valid for general strip widths W .
In Fig. 16, we show the dependence of various quantities as a function of the AB flux Φ on the 16 × 4 ring, n = 0.5, and V SO /t = 0.6. In Fig. 16(a) the current of spin-up electrons is shown for each chain. Following the behavior above observed for the W = 3 strip, the currents no longer satisfy the conditions J ↑,y = −J ↑,4−y+1 as soon as the AB flux is nonzero, but the relations J ↓,y = J ↑,4−y+1 still hold. The SO currents on each chain are shown in Fig. 16(b), where it can be verified the symmetry condition, J so,y = J so,4−y+1 . The total S z (y) in each chain, which provides a measure of partial spin accumulation on each chain, is shown in Fig. 16(c). It should be noticed that the largest difference in S z (y) occurs between the symmetrical inner chains y = 2, 3. This largest contribution to the total spin accumulation from the inner chains occurs for all V SO /t (except the pathological case of V SO /t = 0.8) as it will be shown in Fig. 17(c). Finally, in Fig. 16(d) we show the total spin accumulation, defined as ∆S z = S z (1) + S z (2) − (S z (3) + S z (4)), confirming that its maximum value around Φ = 0.4 is mainly originated by the difference in S z (y) in the inner legs. In Fig. 17(a), the flux-averaged spin-up electron currents on the external and internal legs of the four-chain ring are shown as a function of V SO /t and for the values of U indicated on the plot. To be more precise, we plot the difference J ↑,y − J ↑,4−y+1 , or equivalently, J ↑,y − J ↓,y , y = 1, 2. This average is equal to J ↑,y at Φ = 0. The behavior of the flux averaged differences have the same behavior of J ↑,y (Φ = 0), as it can be assessed by comparing it with Fig. 5(c) (differences between the results in both Figures are just finite size effects). It is interesting to notice that the current on the inner leg has the opposite sign of the current on the outer legs, and that both currents reach their maximum absolute value at V SO /t ∼ 0.5. After this point the inner current starts to increase becoming positive at V SO /t ∼ 0.7 while the outer current decreases. At V SO /t ∼ 0.8, the current in the inner leg becomes larger than the one of the outer leg.
In Fig. 17(b) SO currents in inner and outer legs are shown. In spite of an irregular behavior due to the finite size involved in the calculations, it can be seen observed that SO currents roughly follow the behavior of the up currents in the inner leg but their behavior is opposite for the outer legs. More relevant for spintronics applications is the behavior observed for the absolute value of ∆S z n = S z n − S z 4−n+1 , shown in Fig. 17(c). It can be observed that ∆S z 2 , corresponding to the inner legs, definitely dominates over ∆S z 2 , corresponding to the outer legs, and determines the overall behavior of the resulting spin accumulation, ∆S z . The absolute value of this quantity is plotted in Fig. 17(d), and it can be observed that its maximum is obtained at intermediate values of V SO /t, and also in agreement with previous calculations on the 2-leg ladder 19 , it is enhanced by the Hubbard repulsion U.

B. Voltage bias
In this subsection we present results obtained with DMRG on strips (or generalized "ladders") with two and three coupled chains (or "legs") with OBC in both directions. The purpose is to determine the time evolution of the Rashba helical currents following the sudden application of a bias voltage between both halves of the strip, as explained in Section II.
Similarly to the previous case, the applied voltage bias installs a net hopping current in the longitudinal x-direction, J hop , together with a spin-flipping current, J SO , and a total current, J tot = J hop + J SO . This net SO current implies in turn that the two contributions J ′ SO,x,y,x and J ′′ SO,x,y,x will no longer cancel each other, and as a consequence, the relation on each link, J ↑,x,y,μ = −J ↓,x,y,μ no longer holds. The antisymmetry of the system still implies that y = 1, . . . , W . Of course, no net SO current, which is a charge current, in the transversal direction will appear. A more subtle analysis would show the existence of transversal spin currents 11 , leading in turn to the spin accumulation that we will describe below. As in Subsection III C, the currents on each chain are measured at the center of each chain, x = L/2, and taking advantage of the condition Eq. (7), only the currents J ↑,l,x will be shown in the figures below.
Let us start to examine the case of a strip with W = 2, or "two-leg ladder". Although some results for the total hopping or spin-conserving, and SO or spin-flipping currents, were already reported, 19 the main goal here is to explain their behavior in terms of the RHC studied in the previous Section, which are only slightly modified after applying the external voltage bias. Fig. 18 shows the time-evolution of the up-spin electron currents on both chains of a 24 × 2 strip, following the application of a voltage bias at time τ = 0. Figs. 18(a) and (b) correspond to electron fillings n = 0, 5 and n = 1 respectively. For the sake of clarity, we have inverted one of them, J ↑,2 , since J ↑,2 = −J ↑,1 = J ↓,1 = −J ↓,2 in the static case, that is, at the start of the time evolution. As it can be observed, as τ increases, both currents evolve following different oscillatory behaviors around their τ = 0 values. This difference leads to a nonzero total hopping current along each chain, J ↑,1 + J ↓,1 , and a net total current on the strip, as it will be further discussed below. During the whole time evolution the spin-down electron currents on each leg remain exactly equal to their spinup counterpart on the opposite leg, that is J ↑,2 = J ↓,1 , J ↑,1 = J ↓,2 . Hence, on each chain there will appear a nonzero polarized current, J pol,l = J ↑,l − J ↓,l , l = 1, 2 that cancel each other leading to an overall zero polarized hopping current along the strip. For both densities considered, the effect of the Hubbard correlation is to suppress both the RHC currents as well as their difference. For n = 0.5 ( Fig. 18(a)) the system remains metallic for all U considered, as it can be inferred from (Fig. 19(c)), and the RHC are smoothly decreasing while keeping their sign. On the other hand, at half-filling, (Fig. 18(b)) the RHC, as well as their differences, decrease more rapidly, and at U ≈ 3 the RHC change their sign and they virtually cancel each other, leading to a vanishing overall current (Fig. 19(d)).
The result of a systematic analysis of the behavior illustrated in Fig. 18 for the 24×2 strip is presented in Fig. 19. At n = 0.5, U = 0, the current of spin-up electrons in leg "1" at τ = 0 ( Fig. 19(a)), which as discussed in the previous paragraph is approximately equal to J pol,1 /2, follows the general behavior illustrated in Fig. 5(a) for rings, with an approximate initial quadratic behavior as a function of V SO /t, and becoming increasingly suppressed for V SO /t > 0.6. The effect of the Hubbard repulsion is to suppress these RHC, and this suppression is more pronounced near V SO /t = 0.5. The total currents for this filling, shown in Fig. 19(c) for various values of V SO /t, are suppressed following the relation predicted for Luttinger liquids, 66 J tot (U )/J tot (0) = K ρ where K ρ is the one computed for Hubbard two-leg ladders 67 although the Rashba term could lead to corrections as observed in onedimensional systems. [68][69][70] It is also worth to emphasize the striking similarity between the behavior of the RHC (Fig. 19(a)), as a function of V SO /t, with the one of the spin accumulation as it can be seen in Fig. 4(c) of Ref. 19. Of course, for U = 0, this similarity could be eventually traced back to their common origin in the spin-Hall effect. 11 The fact that their behavior with in- creasing U is the opposite is typical of Hubbard models, where magnetic orderings are in general opposed to transport.
Results at half-filling are strikingly different. The Rashba helical currents at τ = 0, shown in Fig. 19(b), are not only suppressed by U but they also change their sign for large Hubbard repulsion. This behavior is somewhat erratic, particularly for V SO /t = 0.4. Notice that for large U , RHC become quite large in absolute value and their behavior is similar to the one for n = 0, 5. This behavior of the helical current is loosely correlated with the vanishing of J tot with U ( Fig. 19(d)), which indicates the Mott-Hubbard type of metal-insulator transition that is present in Hubbard two-leg ladders close to half-filling. 71 We try now to see if the previously discussed behaviors could be also found in wider strips by studying the strip with W = 3. As it is well-known, numerical difficulties of applying DMRG to wider strips increase exponentially, hence most of the results shown below corresponds to a small strip, 8 × 3, and finite size effects would be estimated by comparing with results for the 16 × 3 strip, although in most cases we could only make qualitative statements.
Let start with the analysis of the time evolution of the RHC, computed on the 16 × 3 strip, shown in Fig. 20(a). In this case, at τ = 0 the RHC in the external chains, y = 1 and 3, have the same antisymmetric properties as the two chains in the W = 2 strip, and these currents vanish on the central chain. As in the W = 2 case, the voltage bias leads to oscillations of these currents around their τ = 0 values. Relatively small differences appearing between antisymmetric currents lead now to the appearance of a net charge hopping current as well as of a net charge SO current. The currents in the central chain, J ↑,2 and J ↓,2 evolve with τ around zero but having the same sign, that is, with net charge current and vanishing polarized current. The time evolution of RHC in the external chains, are similar to the ones observed for the W = 2 chain, where not only a much stronger suppression with U is observed, but also reversal of direction.
As discussed above, differences appearing between J σ,y currents, σ =↑, ↓, y = 1, . . . , 3, lead to net total SO currents, shown in Fig. 20(b), net total hopping currents, shown in Fig. 20(c) which summed give the total net current. Notice that J hop is, as before, about an order of magnitude smaller than the RHC, but still it is possible to reliably assign a plateau value to each of them although it is clearly seen that they are suppressed by the electron Hubbard U term. The total SO currents are in turn an order of magnitude smaller than J hop , and clearly, from Fig. 20(b), much more precision and consequent computer power would be required to assign to them an overall value. However, it seems that it can be observed a trend with U of not only suppressing these SO currents but also to reversing its direction, which starts to occur just at the end of the time interval considered. This behavior is more clear for the smaller 8 × 3 strip, where larger times can be achieved with our currrent computer facilites, and this direction reversal is consistent, although not as definite, with the behavior observed for the W = 2 strip. 19 In the same fashion, it is not reliable to assign a single value to the evolution of the spin accumulation, which is defined as ∆S z = S z 1 − S z 3 , where S z y is the sum of all the average z-projection of the spin on each site on chain y. However, as it can be seen in Fig. 20(d), results point to a trend of enhancement of ∆ with U , which again is consistent with the behavior observed for the W = 2 strip. However, for V SO /t = 0.4 and 0.6, the behavior is not monotonic with U , and this is another indication that the W = 3 strip may be not in the same class of the W = 2 strip, which is not unexpected. Clearly, far more computational effort on the W = 3 strip would be needed to set this question. Another interesting feature that can be observed in Fig. 20 is that ∆S z reaches its maximum with a delay with respect to the time at which J hop and J SO reach their maximum. This delay is an indication of the presence of transversal spin currents appearing after the charge currents have been installed and leading to the spin accumulation, as expected in the SHE mechanism with open boundaries.
Finally, a more systematic dependence of the RHC as a function of V SO /t and U is presented in Figs. 21(a) and (b), for n = 0.5 and n = 1, respectively. Although J ↑,1 have nonmonotonic variations with V SO /t, their variation with U is noticeable weaker and smoother for n = 0.5 than for n = 1. Actually, for n = 1, for most values of V SO /t, J ↑,1 changes its sign as U is increased, similarly to the previous analyzed case of W = 2 strips. Again, this behavior occurs while the total current, depicted in Fig. 21(d), goes to zero, indicating a Mott-Hubbard metal-insulation transition. On the other hand, the be- havior of J tot , for a given V SO /t decreases with U following the general metallic behavior of quasi-onedimensional systems. We are not aware of Luttinger Liquid studies on the Hubbard model on W = 3 strips, and consequent estimations of K ρ .

V. CONCLUSIONS
In this manuscript, we have described in the first place the spatial texture of Rashba helical currents in planar strips, as a function of various parameters, such as RSOC strength, strip width, electron filling, in the absence of external electromagnetic sources, and for the noninteracting case. By numerically solving the single particle Hamiltonian Eq. (1), we found that these RHC are mostly concentrated at the strip edges, but as the RSOC strength and/or the width are small enough, they appear as spread over the whole system. This is a different behavior to the one theoretically predicted for the helical edge currents, characteristic of the QSH, and topological insulators, as resulting for example from the Kane-Mele model. Certainly, this difference could be traced to the difference between the Haldane-like term for counter propagating spin up and down electrons, present in that model, which leads to the opening a gap in the energy spectrum, and hence to an insulator state in the bulk, and the Rashba SOC which favours the metallic state.
The results also suggest a periodicity of the sign oscillations of the RHC that could originate as Friedel oscillations induced by the edges, thus explaining its filling dependence. As mentioned in the Introduction, the natural source of these oscillations are local impurities, and in fact, the current inhomogeneities reported in Ref. 56 could be of the same nature as the present RHC. It would be interesting to perform similar calculations on strips of the Kane-Mele model, for the sake of comparison.
In the second part of the present work, we consider the effects of charge currents induced by external electromagnetic sources, which could be relevant for actual spintronic devices. Our main results are that the RHC remain in the presence of such induced net currents, and their differences on each chain with respect to their stationary equilibrium values, help to explain how charge currents, both of the spin conserving or the spin flipping types, are distributed on each chain of the strip. These RHC also allow us to predict that on each chain there could be polarized charge currents, although the total current on the strip would be unpolarized. 64 Since, as predicted by our handwaving model presented in Section III A in systems where Rashba and Dresselhaus SOC have equal strength, the RHC would not exist, and taking into account the relation between RHC and the net currents, it would be tempting to propose an explanation for the vanishing of optical conductivity at the persistent helix point. 72 In this study, we also include the presence of on-site electron correlations given by the Hubbard repulsion. In the context of topological insulators, it has been emphasized that the robustness of helical edge states with respect to interactions is a result of time-reversal symmetry, and in particular, the effect of a Hubbard term added to the Kane-Mele Hamiltonian has been thoroughly analyzed. 63,73 In the model here considered, the RHC remain for finite values of U , although in the quasi-onedimensional metallic state, their are suppressed roughly following the prediction of Luttinger liquid theory. In fact, even if the system undergoes a Mott-Hubbard type of metal-insulator transition, as detected by the vanishing of the total current, the RHC are still present.
In principle, the Rashba helical currents would be as difficult to detect as the helical edge currents of the QSH. For the case of QSH, it has been shown that experimental evidence for helical edge currents can be achieved by multiterminal probes. [74][75][76] However, this scheme depends heavily on the strict edge character of QSH helical currents, and/or on the quantized nature of these currents. Clearly, a lot of ingenuity will be needed to adapt those experimental setups to measure the Rashba helical currents we observe in the present effort, and particularly the presence of polarized currents varying across the strip section. In addition, the study of a dissipationless character of the RHC is out of the scope of the present work. Notice that our model is dissipationless in the sense that it is not coupled to external degrees of freedom such as phonons.
Finally, we would like to emphasize apparent differences with some previous theoretical studies on this model. The band we consider is a cosine one, and we considered fillings n ≥ 0.5, hence it is not simply reducible to a parabolic band, and, more importantly, the finite width of our system, and eventually also its finite length, sets a departure from previous studies on the infinite twodimensional system. Certainly, finite size systems are per se relevant for nanoscale applications.