Robustness of Helical Edge States Under Edge Reconstruction

The helical edge states of time-reversal invariant two-dimensional topological insulators are protected against backscattering in idealized models. In more realistic scenarios with a shallow confining potential at the sample boundary, additional strongly interacting edge states may arise, that could interfere with the topological protection of edge conduction. We find that interaction effects within the reconstructed edges are well described by the Luttinger liquid model. While interactions between this Luttinger liquid and the helical edge states can in principle give rise to dynamical spin polarization and the breaking of time-reversal symmetry, we demonstrate that random spin-orbit coupling strongly suppresses such dynamical spin polarization, resulting in the persistence of near quantized edge conduction.

Despite the amazing observation of quantized ballistic edge conductance, such behavior has yet to be observed in a macroscopic regime over distances beyond tens of microns. There have been many attempts at explaining how TRS-breaking perturbations can limit the mean free path for edge transport. This includes dynamical polarization of nuclear spins interacting with the helical states [36][37][38][39][40], the presence of Rashba spinorbit coupling [28,[41][42][43] in combination with interaction mediated mechnanisms [44,45], scattering off magnetic [46][47][48][49][50][51][52][53], and non-magnetic impurities or electron puddles [6,[54][55][56]. The physical conditions for the occurrence of TRS breaking were considered from a new point of view through a scenario with reconstructed edge states [9], including spontaneous ferromagnetism at the edge. For integer quantum Hall systems, it is well known that a smooth confinement potential (as opposed to a model of sharp confinement due to the termination of the sample) can lead to the emergence of additional edge states, an effect dubbed edge reconstruction [57,58]. This effect can give rise to a non-monotonic variation of the filling fraction in combination with ferromagnetism. Consequences of such edge reconstruction in QSH systems include en-hanced back-scattering due to the local breaking of time reversal, as well as the possibility of engineering spin filters in mesoscopic devices [9]. The formation of edge states at the boundary between Mott insulating stripes and a topological QSH bulk was considered in [59,60].
In this letter, we show that fluctuation effects are likely to protect the QSH edge conductance against TRS breaking due to edge reconstruction. We argue that the reconstructed states exist in a well defined, experimentally accessible regime, and potentially lead to a dynamical spin polarization at the edge. By examining the effects of weak random Rashba spin-orbit coupling, we show that edge transport is generically protected in this scenario. Only in the presence of strong disorder, which strongly localizes the reconstructed edge states, can fluctuating local moments be dynamically polarized via spin flipflop scattering, and thereby mediate backscattering in the helical liquid similar to the mechanism described in Refs. [47][48][49]. Taken together, this analysis demonstrates that generic edge reconstruction in QSH insulators may not be sufficient to induce backscattering in the presence of a moderate amount of Rashba disorder, further strengthening the case for the possibility of observing protected edge transport on macroscopic scales. Edge reconstruction: We consider the Bernevig-Hughes-Zhang (BHZ) Hamiltonian [2] on a strip with N y lattice sites in y-direction, and periodic boundary conditions in the x-direction such that the momentum k x is a good quantum number. We employ open boundary conditions defined by the termination of hopping matrix elements outside the strip. Electron interactions are described by with c nkσ annihilating an electron with band index n = 1, . . . , 2N y , momentum k, and spin σ. V nn ,σσ k1k2k3k4 denotes the interaction matrix element. Since the sharp arXiv:2105.14763v2 [cond-mat.mes-hall] 30 May 2022 shows the selfconsistent band structure in an effective model which includes both the helical and the reconstructed states and shows the splitting between spin up bands (full lines), and spin down bands (dashed lines). Panel b) is a zoom-in into the spectrum close to the Fermi level and for momenta close to the boundary of the BZ. Panel c) shows a zoom-in of the dispersion of the helical modes close to the center of the BZ. Panel d) shows the effective electrostatic potential VES and the confining potential VC close to the left edge. A state may have weight on a lattice point if the total potential less than the bandwidth. Panel e) illustrates that the wave functions of the reconstructed modes do not depend on momentum k such that φ k (y) ≡ φ(y).
edge termination does not serve as a realistic model of the sample edges, we include the effects of the positively charged ions with density n ion which gives rise to an electrostatic potential V ion . Further, we introduce a confining potential V C , which grows linearly with slope m near the edges of the system and vanishes in the bulk. The interplay of both terms affects the presence of additional edge modes at the boundary of the Brillouin zone (BZ). We choose the ion density to drop abruptly from its bulk value to zero at a distance of w lattice sites measured from the sharp edges of the sample. In addition, we assume the confining potential to vanish for lattice sites j ≥ w. We begin by treating interactions within the Hartree-Fock (HF) approximation and apply a self-consistent electrostatic modeling of the sample edges. The Hamiltonian is given by H = H BHZ + V ES + V C (j), with V ES = V H + V ion (j) denoting the total electrostatic potential containing the Hartree potential V H and the ion potential [61]. By varying the slope m of the confining potential, we control the density of the additional edge modes, and we tune the system such that these additional modes do (as in Fig. 1) or do not exist. Similarly to Wang et al. [9], we find that for sufficiently small values of m, the system evolves into a symmetry-broken ground state [9]. This occurs for a lattice constant a 6.5 A appropriate for HgTe/CdTe quantum wells, complementing the finding of Ref. [9] for a 3.14 nm. Moreover, we find that the HF solution is equivalent to that in a reduced model in which we consider the helical and the reconstructed states only. Inspection of the wave function of the reconstructed states φ k (j) shows that its transverse wave function is independent of k, and strongly localized on a single lattice point as shown in Fig. 1e). This is an important difference compared to the quantum Hall states with φ k (y) = φ(y − k 2 ), where the momentum dependence of the wave function stabilizes ferromagnetism. Thus, the reconstructed modes can be mapped onto a model of one-dimensional (1D) electrons. For sufficiently low electron density, the ground state of the 1D electron gas in a translationally invariant HF calculation is always ferromagnetic such that the symmetry-breaking does not rely on the presence of the helical modes, but can be seen as a generic characteristic of the mean field solution at low electron density. An investigation of the overlap between the wave functions of the helical and reconstructed states indicates that the exchange coupling within the reconstructed states is a factor 20 larger than the exchange coupling between reconstructed and helical states. For this reason, the presence of the helical states does not affect the spin polarization of the reconstructed ones, and ferromagnetism of the reconstructed states can be studied on its own. 1D ferromagnetism and fluctuations: In 1D, fluctuations have a strong effect on the stability of ferromagnetism in the reconstructed edge states. For 1D electrons, the dispersion relation can be linearized around the Fermi points, which makes bosonization and a non-perturbative treatment of interactions possible, with a Hamiltonian given by [62] (2) Here, the free part of the Hamiltonian is given by with ν = ρ, σ denoting the charge or the spin part, u ν the excitation velocity, K ν the Luttinger parameter, and with fields obeying bosonic commutation relations . As can be seen from the previous equations, the Hamiltonian decomposes into a charge and a spin channel, and the terms with c 1 and c 2 are necessary to describe FM in the Luttinger liquid [63]. the resulting theory is similar to the Landau functional for the description of phase transition in classical systems. Table I shows the products and ratios of excitation velocities and Luttinger parameters for three different fillings of the reconstructed modes. The values were obtained using the interaction matrix elements and Fermi velocities from the self-consistent electrostatic solutions obtained in the previous section (for details see supplemental material [61]).
The Luttinger liquid fixed point described by the Hamiltonian Eq. (3) with c 1 = c 2 = 0 is characterized by a dynamical exponent z = 1 such that naively any kind of long-range order seems to be forbidden by virtue of a mapping to a classical 2D model at finite temperature with a continuous symmetry [65]. It has however been shown that u σ /K σ can become negative such that the formally irrelevant higher order terms of the fields and higher order gradients must be included to ensure stability [63]. The resulting effective theory as presented in Eq. (3) has dynamical exponent z = 2 such that a ferromagnetic ground state can in principle be stable to fluctuations, even in 1D [66].
In the above discussion, we have assumed that the bulk is magnetically inert, and that magnetism only forms in the the narrow one-dimsional strip occupied by the reconstructed states. At half filling, magnetism in the bulk is only possible when an excition condensate forms, which was predicted to happen in inverted bilayer quantum wells [67] but is unlikely to happen in the BHZ model.
In the case of a strongly screened and thus shortranged interaction with on-site interaction strength U , the transition into the ferromagnetic LL phase happens for U > 2π v F , using the criterion of a negative coefficient u σ /K σ calculated naively without considering renormalization effects. On the other hand, a restricted Hartree-Fock calculation yields the Stoner criterion [68] U > π v F , agreeing with the above estimate within a factor of two.
In the limit where excitations into higher transverse subbands have a negligible amplitude, and when the interaction range is less than the average interparticle spacing, the interaction is determined from on-site matrix elements of the interaction V , taken with regards to self-consistently determined wave functions. The Lieb-Mattis theorem [69] states that in this limit, the ground state can never be ferromagnetic, independent of the on-site interaction strength, or filling of the system. When applying the Lieb-Mattis theorem to the strictly one-dimensional reconstructed states (see Fig. 1 e), the prediction of ferromagnetism from Hartree-Fock and bosonization is questionable, at least in the limit of a well screened and short ranged interaction.
For the opposite limit of a long-range Coulomb interaction, both an unrestricted Hartree-Fock calculation and the bosonization approach [70] predict the formation of a charge modulated state known as a Wigner crystal [71]. The ground state of the Wigner crystal does not have ferromagnetic order, but instead, shows a tendency towards antiferromagnetic ordering. The Wigner crystal regime is realized if the screening length κ is much larger than the interparticle distance, that is κ π/(2k F ), with k F the Fermi momentum. The transition into the ordered state happens in 1D if [72] with m * denoting the effective mass, m 0 the bare electron mass, a the lattice constant, a 0 the Bohr radius, and r the dielectric constant. For a realistic screening length of κ = 30a [13], the critical Fermi momentum for ferromagnetism determined from HF is k F,FM · a π/20 ≈ 16k F,WC · a [61]. k F,FM is slightly increased as κ → ∞.
This HF analysis predicts the existence of a spurious transition between a Wigner crystal and ferromagnet in the limit of long range Coulomb interactions. In addition, the HF prediction of ferromagnetism for a local Hubbard interaction is also in disagreement with the Lieb-Mattis theorem [69], casting serious doubt on the stability of HF ferromagnetism for intermediate values of the screening length κ as well. Thus it appears there is no evidence for the occurrence of a generic FM phase. Even if one takes the HF solution for finite κ at face value, for a realistic choice of parameters it would only occupy a small fraction of Fermi momenta compared to the total extent of the BZ. In conclusion it seems unlikely that ferromagnetism will significantly compromise the robustness of edge conductance in topological insulators. Dynamical polarization: While we have argued above that static ferromagnetic order is unlikely to be realized in the reconstructed states, we now investigate the effect of dynamical spin polarization on the quantization of edge conductance. Going beyond HF, there are lowenergy spin exchange processes between the helical (h) and reconstructed states (r) described by the Hamiltonian: According to the discussion in the previous section we discard the possibility that intra-edge interactions lead to a phase transition, neglecting them for the moment, and consider the total Hamiltonian H = H r + H h + H ex where: Here Ψ r,α (x) annihilates an electron with spin α at position x, J is the spin exchange coupling, the spin operator with σ the vector of Pauli matrices, and S ± (x) = S x (x) ± iS y (x) the corresponding ladder operators. τ = ±1 describes the two chiralities of the reconstructed modes, and a summation over repeated spin indices is implied.
In first order perturbation theory, matrix elements for flip-flop scattering will generically vanish due to the inability to simultaneously satisfy energy and momentum conservation in the presence of different velocities for the two edge modes. In the presence of potential disorder however, flip-flop scattering is possible, and the difference in densities of up and down spins of the reconstructed modes can be characterized by a chemical potential difference ∆µ, leading to a net spin density in the z-direction: S z r = ∆µ/(2π v r ). Even without knowledge of the detailed dynamics of spin polarization, it is possible to see that in the asymptotic state, the spin density is fixed by ∆µ = eV , where all polarization processes cease.
The dominant effect of intra-edge interactions is to renormalize the density of states close to the Fermi level, and as a consequence the spin susceptibility. Within the previously discussed Luttinger liquid model, the spin density is modified to S z r = ∆µ(K σ /2π u σ ). Due to the renormalization factor K σ /u σ , the spin density of reconstructed modes can be significant, even for small transport voltages (see table I).
This spin polarization acts like an external Zeeman field on the helical states, and in the presence of Rashba disorder can give rise to backscattering of the helical states and a finite mean free path. However, it turns out that Rashba disorder not only gives rise to backscattering in the presence of an effective Zeeman field due to the spin polarization of the reconstructed states, but also strongly renormalizes this spin polarization. In order to discuss this effect, we describe random Rashba spin orbit coupling [28,41,43,73]

(SOC) by the Hamiltonian
We assume that the spatially random coupling has Gaussian correlations: a(x)a(x ) = A 0 · δ(x − x ). We find the backscattering mean free path of the helical states to be given by [37]: For the random Rashba disorder in Eq. (7) the wave functions for the reconstructed states can be determined analytically [61]. Scattering processes due to Rashba that do not preserve chirality are in principle possible. Their matrix element in Fourier space is proportional to (k + k ), where k and k are the momenta of initial and final states, respectively. For different chiralities, these momenta have the same magnitude but opposite sign such that their sum vanishes. Thus, the main contribution is due to forward scattering processes. We obtain for the approximate eigenstates of H = H r + H R ψ = e iθ0τ 0 σ y e iθxτ x σ y e iθyτ y σ y e iθzτ z σ y ψ 0 , where for example θ 0 (x) = −2k F,r x 0 dx a(x )/ v r and ψ 0 does not depend on the position variable (for details see supplemental material [61]). For the specific choice of ψ 0 = |+, ↑ , due to random Rashba SOC the expectation value of the z-component of spin is reduced by a factor 2 R /L, where R = 2 v 2 r /(A 0 k 2 F,r ) is the characteristic Rashba length and L the system size.
In the presence of both potential and Rashba disorder, the reconstructed states are no longer extended states, but become localized on a length scale dis = 2 v 2 r /U 0 and thus do not contribute to the edge conductance. Here, U 0 is the variance of the disorder potential. For finite U 0 , the expectation value of spin is reduced by a factor of 2 R / dis resulting in and an associated mean free path Stronger exchange coupling or higher transport voltage both increase the effective Zeeman field and thereby lower the mean free path. Calculating for a reasonable set of material/model parameters: v h = 2.2 · 10 5 m · s −1 , v r = 1.83 · 10 5 m · s −1 , J = 0.165 eV · a, A 0 = 4.5 · 10 −4 · meV 3 ·µm 3 , U 0 = 5·10 −2 meV 2 ·µm, and a transport voltage of eV = 0.5 meV, yields a mean free path of macroscopic length. The naively calculated value without the inclusion of the Rashba reduction factor is 50 µm for eV = 0.5 meV. Our analysis shows that the presence of random Rashba SOC substantially increases , supporting near quantized edge conductivity in QSHI's.
In the weak disorder limit, A 0 and U 0 are proportional to the impurity density. As long as the distance between impurities is larger than the range of the impurity potential, A 0 and to U 0 will be proportional to each other, since A 0 depends on the the gradient of the electric field and U 0 on the gradient of the impurity potential itself. Thus, A 0 ∝ U 0 [74] such that diverges in the weak disorder limit. In the opposite limit of strong disorder with dis k F,r 1, the reconstructed states can be seen as individually localized magnetic impurities. Their emerging magnetic moments can induce backscattering in the helical modes [47][48][49] and under the assumption that each spin acts as a independent scatterer we estimate a lower bound for the mean free path of 5 nm [61]. Conclusions: We have argued that the edge conductance in quantum spin Hall systems is remarkably robust against the emergence of additional edge states due to edge reconstruction. Quantum fluctuations mask the existence of a time reversal breaking ferromagnetic Hartree-Fock solution, giving rise to a strongly enhanced magnetic susceptibility instead. As a consequence, dynamical processes can lead to a spin polarization of the reconstructed edge, that in combination with random Rashba coupling, yields a nearly macroscopic finite mean free path on the order of hundreds of microns. We emphasize that disorder plays a crucial role here in hindering backscattering by parametrically reducing the effects of dynamical spin polarization. Only when disorder is strong enough to produce localized magnetic moments, do we predict stronger backscattering in line with experimental measurements on ballistic edge transport. Thus future applications of quantum spin Hall devices remain promising, provided that sufficiently clean samples outside the strong disorder regime can be fabricated. with m denoting the slope, j = 1, . . . , N y the lattice site, and Θ(x) the step function with Θ(0) = 0. The lattice site w separates the edge from the bulk. Further, we include the effects of ionic background with an ion density with n 0 denoting the ion density in the bulk (see Fig. S1). The interaction Hamiltonian is given by with ψ τ σ (r) annihilating an electron at position r with orbital index τ = (E1, H1) where we use the notation of Ref.
[S1], spin σ, and with U denoting the on-site interaction strength, κ the screening length, and r the relative dielectric constant. Note that the position vector has to be understood as r = (a · j x , a · j y ) T , with j x , and j y labelling lattice sites in the xand y-direction. We express the field operators as with φ nkσ (j y ) denoting the normalized eigenfunction for the transverse y-direction localized at unit cell j y , band index n, momentum k, and spin σ. α τ,n are the expansion coefficients when going from the orbital to the band basis, and L x = N x · a the system length in x-direction. Since we have N y lattice sites with an electron-and a hole-like state, we can have a total of 2N y different bands per spin such that n = 1, . . . , 2N y . In the following, we will use j y ≡ j, and k x ≡ k. Substituting Eq. (S5) into Eq. (S3), and using the orthogonality relation τ α * τ,n α τ,n = δ nn we obtain the Hamiltonian given in the main text Eq. (??) with the interaction matrix element given by and the effective one-dimensional interaction potential We use that c † nkσ c n k σ = f 0 (ε σ n (k)) · δ nn δ kk δ σσ with f 0 (ε) denoting the Fermi function, and ε σ n (k) the dispersion of band n with spin σ and momentum k and obtain the Hartree matrix element with the scaled electron density n e (j) defined as The matrix element for the electron-ion interaction is the same as in Eq. (S8), except that n e → −n ion . The total electrostatic term V ES = V H + V ion is thus given by For the Fock matrix element we obtain We selfconsistently determine the electronic occupation number of the different bands by minimizing the total energy obtained by summing the dispersion relation ε σ n (k) = ε 0,n (k) + ε σ F,n (k) over occupied states. The first term is the selfconsistent Hartree dispersion and the second one is the Fock contribution. The selfconsistent solution including the helical and reconstructed bands only is shown in the right panel of Figure S2.

MORE ON FERROMAGNETISM
In the intermediate regime where κ is still well above the interparticle spacing, but only a few times the interparticle distance, we can use restricted HF theory to formulate a criterion for the existence of the ferromagnetic solution on the Fermi momentum of the reconstructed states: with K 0 (x) denoting the modified Bessel function of the second kind, where u ≈ a/20 is a short-distance cutoff such that U = e 2 /(4π 0 r u). The Stoner criterion from the main text is reproduced at κ a/4.

LUTTINGER LIQUID ANALYSIS
The Luttinger liquid is described by the Hamiltonian Eq. (??) defined in the main text. It is characterized by the excitation velocity u ν , and the Luttinger parameter K ν , with with ν = ρ, σ for charge or spin. They are given by [S2] We obtain for the interaction constants with the subscript r denoting the reconstructed states, and k F,r their Fermi momentum. However, below a certain value of the electron density the ratio u σ /K σ can become negative [S3]. Yang [S3] argues that at this point, the standard LL description breaks down and higher order terms of the gradient and the field itself must be kept in the Hamiltonian to ensure stability. The additional terms are where c 1 and c 2 are positive constants determined by the interaction potential and the single-particle dispersion. Including these terms, it was shown that there is a ferromagnetic transition which contradicts the results from the standard LL theory [S3].    [S4], and the scalar disorder variance was inferred from [S5].