RKKY to Kondo crossover in Helical Edge of a Topological Insulator

Two spatially separated magnetic impurities coupled to itinerant electrons give rise to a dynamically generated exchange (RKKY) inter-impurity interaction that competes with the individual Kondo screening of the impurities. It has been recently shown by Yevtushenko and Yudson (2018), that the RKKY interaction and the RKKY vs. Kondo competition become nontrivial on helical edges of two-dimensional topological insulators where there is lock-in relation between the electron spin and its direction of motion. Kondo screening always takes over and dominates at large inter-impurity distances and it can also dominate all the way to short distances if the Kondo coupling is sufficiently large and anisotropic. In the present paper, we study the Kondo-RKKY competition in detail on a qualitative and quantitative level. For this we employ the numerically exact numerical renormalization group (NRG) for a broad parameter scan of two Kondo coupled impurities vs. magnetic anisotropy, impurity distance and temperature, and comment on the role of finite bandwidth. We give a pedagogical introduction on the the setup of the two-impurity setting within the NRG in the helical context. Overall we establish a plain crossover from RKKY to Kondo with increasing impurity distance which permits an intuitive physical picture by simply comparing length scales set by the Kondo screening cloud, the thermal length scale vs. the impurity distance.


I. INTRODUCTION
Physics of magnetic impurities (MIs) coupled to helical electrons on one-dimensional (1D) edges of two dimensional (2D) time-reversal invariant topological insulators (TIs) attracted the attention of researchers soon after the experimental discovery of the TIs [2][3][4][5].This interest resulted from a search of possible backscattering mechanisms which could make the virtually protected helical conductance sub-ballistic in relatively long samples [6][7][8][9].Since the helical electrons possess a lock-in relation between the spin projection on the quantization axis and the direction of propagation (the so-called chirality), their backscattering is expected to involve some nontrivial spin processes, e.g. the spin flip.The MI can provide such an inelastic backscattering of the individual helical electrons.However, the helical conductance can be suppressed only if the spin conservation on the edge is violated, see papers [10,11] and references therein.If the edge is spinconserving and the MI does not break the spin U(1) symmetry, it backscatters the helical electrons but cannot influence the dc conductance [12].The anisotropic MI is able to suppress the helical conductance only if it breaks the spin conservation and is not Kondo screened [13][14][15].The latter requires either the temperature being larger than the Kondo temperature, T > T K , or a large value of the MI spin, S > 1/2.This points out the importance of understanding the Kondo effect in TIs which is substantially different from that in usual (non-helical) 1D wires in the presence of the electron-electron interaction or the magnetic anisotropy of the XXZ type [16].
Here as well as throughout this paper, T K ≡ T (1) K represents * weichselbaum@bnl.gov the Kondo scale of a single, possibly anisotropic spin-half impurity coupled to a helical edge [16].Its value is identical to the plain non-helical Kondo model given that that helicity in the non-interacting bath is irrelevant from the point of view of a single impurity.
The Kondo effect can be suppressed by the indirect exchange MI interaction, the Ruderman-Kittel-Kasuya-Yosida (RKKY) interaction [17], if the helical edge is coupled to a dense array of the MIs.According to the simple picture of the Doniach criterion [18], the "winner of the RKKY-Kondo competition" can be found by comparing T K with the characteristic RKKY energy scale E R .The latter means the energy gap which opens after the RKKY correlations lift a degeneracy in the energy of the uncorrelated MIs.If T → 0 but E R ≫ T K , the RKKY correlations overwhelm the Kondo screening which may lead to many nontrivial effects, including Anderson localization of the helical electrons caused by the random magnetic anisotropy [19,20], and magnetically correlated phases [21][22][23][24].The RKKY-induced magnetic order in helical 1D systems is nontrivial and qualitatively different from that in their non-helical counterparts, cf.Refs.[25,26] and references therein.
Since T K does not depend on the MI density while E R typically decays with increasing the inter-impurity distance, R -see Eq. (20) below, one can surmise that there is a characteristic distance defined by the equality T K ≃ E R (x c ), which separates the RKKY-and Kondo-dominated phases, R < x c and R > x c respectively.This conclusion is inspired by a misleading analogy with the physics of non-helical wires [27][28][29][30][31].One of us (in collaboration with V. I. Yudson [1]) have recently considered two MIs coupled to the helical edge and shown that, if either the electron-electron interaction or the magnetic XXY anisotropy or both are strong, x c shrinks and the Kondo effect overwhelms the RKKY in-arXiv:2208.02275v3[cond-mat.str-el] 1 Sep 2023 teraction over all macroscopic inter-impurity distances.This unexpected conclusion has been drawn based on phenomenological arguments and on the analytical consideration of limiting cases.This theory reveals many quantitative features but is far from being complete.In particular, it cannot predict whether the above mentioned phases with finite and vanishing x c are separated by a crossover or by a phase transition and whether the MIs remain somehow correlated even in the Kondo-dominated phase.
In the present paper, we expand and complete the theory of Ref. [1].We present an analytical theory of the RKKY correlation between the two MIs coupled to the helical electrons but, as the main working tool, we have chosen the numerical renormalization group (NRG; [32][33][34]).This powerful and well-established method has allowed us to answer the aforementioned open questions.In particular we will show that 1) the different phases are separated by the crossover, and 2) the MIs are uncorrelated, i.e. independently screened, in the Kondo dominated phase.
The paper is organized as follows: We introduce the model in Sec.II, followed by analytical considerations in Sec.III.The remainder of the paper then is dedicated to a detailed analysis and discussion of the model based on the NRG in Sec.IV, concluded by summary and outlook.In appendix App.A, we give a detailed pedagogical derivation of how the helical 2-impurity system is setup and mapped into the standard NRG machinery.In particular, this highlights the possibility to map the system onto a Wilson ladder, with complex coefficients only within the coupling to the impurity.Furthermore, we included in the appendix a brief reminder on the poor-mans scaling of the anisotropic Kondo model, as well as a plain second order perturbative derivation of the RKKY effective Hamiltonian and RKKY energy that is complementary to Sec.III. 1. Helical edge mode of the one-dimensional edge in a 2D topological insulator, e.g., as it occurs in the Kane-Mele model [35].Focusing on low energies, the dispersion vs. momentum k is given by ε kσ = σvk with spin σ ∈ {↑, ↓} ≡ {+1, −1} indicated by the blue (red) line, respectively, where v denotes Fermi velocity.We assume a finite half-bandwidth D, throughout.This gives rise to an effective Brillouin zone with lattice spacing a [cf.Eq. ( 2)].

II. THE MODEL A. Hamiltonian of the helical edge
We study spins coupled to a helical edge mode in a 2D topological insulator.This may be approached in various ways.The edge mode can be simulated (i) by fully modeling an underlying 2D lattice model in real space, such as the Kane-Mele model with the Dresselhaus spin-orbit interaction [35].This naturally introduces a cutoff in terms of bandwidth of the helical edge mode which is located inside the gap of the continuum of extended bulk states.Together with interacting correlated impurities, this may be simulated numerically, for example, using the density matrix renormalization group (DMRG, [36,37]).There the non-interacting 2D lattice without the impurities can be conveniently mapped to an effective 1D impurity setting via Lanczos tridiagonalization [38].However, this approach bears significant overhead in terms of the precise choice of the underlying 2D lattice model and its parameters, the details of which are considered irrelevant for the low-energy physics.Conversely, (ii) the model can be considerably simplified by focusing on a single pure effective 1D helical edge described in energy-momentum space at low energies as depicted in Fig. 1 and described by Eq. ( 1) below.The latter representation of the ballistic non-interacting edge modes may be (iii) exactly transformed into a 1D real-space lattice realization which, however, involves long-range hoppings [cf.Sec.A 2 and in particular Eq.(A12) for more on this].Working in energy-momentum space, instead, is appealing from an analytical point of view [1], but is also perfectly well-suited for the numerical treatment via the NRG.With the additional goal to scan many orders of energy scales with Kondo physics in mind, the energy-momentum representation is preferred over real-space lattice descriptions.Therefore we follow approach (ii) throughout this work.
The standard model Hamiltonian for a single helical edge mode can be written in discrete form in momentum space as follows [cf.Fig. 1], where ĉkσ are the annihilation operators of the helical fermions with spin σ ∈ {↑, ↓} ≡ {+1, −1}, and τ α with α ∈ x, y, z the Pauli matrices.The fermions are described by a linearized dispersion relation with the Fermi velocity v. Below we consider a finite half-bandwidth D (UV cutoff) such that ε kσ ∈ [−D, D].This is required in the numerical context yet and also for the sake of regularization.We assume that the TI edge is oriented along the x-axis, with xy being the TI plane, and the z-direction is the quantization axis for the spins of the helical edge modes.This convention corresponds to the experimentally relevant situation where the quantization axis is often fixed being perpendicular to the TI plane [6][7][8][9].The helical system in Eq. ( 1) respects time-reversal symmetry (TRS), in that ε kσ = ε −k,−σ .The crossing of the spin selective dispersions, i.e., the Dirac point in Fig. 1 is chosen for simplicity at k 0 = 0 without restricting the case.The precise choice of k 0 is irrelevant for our purposes, as it can be absorbed into the definition of the basis states, and hence can be gauged away.[39] By assuming a finite half-bandwidth D and having translational invariance by working in momentum space, this directly implies an effective Brillouin zone (BZ) with momentum range |k| ≤ k max , see Fig. 1, with a discontinuous dispersion across the Brillouin zone (BZ) boundary.Conversely, this defines an effective lattice constant, a = π/k max , via the one-particle dispersion, having |ε( π a )| = D, i.e., a ≡ πv D . (2) Below, we use the local density of states ρ 0 = aρ 1D as experienced by an impurity where ρ 1D = 1/2πv is the constant one-particle density of states of a 1D system with linear dispersion relation [40] (thus units are [ρ 0 ] = 1/energy, whereas [ρ 1D ] = 1/(energy • distance)).Using Eq. ( 2), the local density of states becomes ρ 0 = 1 2D , which is consistent with standard NRG conventions.For more on the effects of finite bandwidth and an effective 1D lattice defined by Eq. ( 2), see Sec.A 2 [cf.Eq. (A11)]).
In numerical simulations, furthermore, adhering to standard NRG conventions, we choose the unit of energy D := 1.By also setting the unit of distance a := 1, this fixes the velocity to v = 1/π.We also set ℏ = k B = 1.
The last expression in Eq. ( 1) provides a more compact notation using the spinor ĉk ≡ (ĉ k↑ ; ĉk↓ ), where the semicolon denotes a column vector.This explicitly shows that the SU(2) spin symmetry is broken in the helical setting.It reduces to the abelian U(1) symmetry with preserved component of the total spin S tot z [this is in contrast to a chiral system where both spins move in the same direction [40,41], which thus preserves SU(2) spin symmetry].In the helical case, we thus only explore the combination of abelian symmetries U(1) charge ⊗ U(1) spin [see Sec.A 5 for further comments on symmetries].
The Hamiltonian (1) can be written in real-space in the standard continuous form where Ψ ≡ ( ΨR ; ΨL ) is the spinor constructed from the slow helical fields of the right and left moving electrons ΨR,L .For the free electrons in the helical edge mode, i.e., the bath to which the spin impurities are coupled, spins up and down and, respectively, chirality (the direction of propagation) of the right/left movers can be denoted in more general fashion by the index in particular, ΨR,L ≡ Ψ↑,↓ .When used as variable in equations below, the particular meaning of this index is always clear in context.The Hamiltonian Eq. ( 3) may be defined on a finite-length system with periodic boundary conditions (BCs) and hence discrete momenta, or in the thermodynamic limit with a continuous energy-momentum space.Sec.II A represents the simplest effective model that describes a single helical edge, for example, in HgTe/CdTe quantum-well heterostructures that possess axial and inversion symmetry around the growth axis.

B. Coupling between the helical fermions and MIs
Let us introduce two MI spins Ŝα η separated at distance x and located symmetrically around the origin at positions for left and right impurity, respectively [to be differentiated from right / left movers denoted by σ ∈ {R, L} in Eq. ( 4)].
By working in energy-momentum space, the bath operators at the location of the impurity are obtained via Fourier transform [e.g., cf.App.A].With this, one can focus both, analytically and numerically, on the two impurities being located along a single edge without having to worry about periodic boundary conditions.This is possible, when the actual (2D) sample is always considered much larger than the impurity distance x.
Thus without restricting the case, one is free to think of the two impurities as being symmetrically located at ±x/2 along a straight edge around some arbitrary but fixed origin.The exchange interaction between the helical electrons and these two MIs is described by the Hamiltonian Here j 0 ≡ ρ 0 J and j z ≡ ρ 0 J z are the constant dimensionless exchange couplings, such that, for example, 2πv j 0 = aJ with J the coupling strength of the impurity in units of energy and [aσ α η ] = 1 dimensionless as typically used within the NRG, having with τ ± ≡ 1 2 (τ x ± iτ y ).The exchange interaction (6) may be anisotropic, having J ̸ = J z , while it always conserves the z-projection of the total (electron and MIs) spin.
In the next section, we adhere to the standard notations of the literature devoted to the analytical study of Kondo impurities coupled to helical electrons, where the Kondo coupling is measured in units of the Fermi velocity, J = (2πv)j, and similarly for Jz .

III. ANALYTICAL THEORY OF THE RKKY REGIME
Let us briefly review the RKKY theory for the helical edge mode coupled to two MIs [1] using field theoretical machinery.

A. Non-interacting fermions
To describe degrees of freedom of the MIs, one should approximately integrate out the helical fermionic edge modes.A natural way to do this is to exploit the formalism of functional integrals.As a result, one arrives at an effective action for the impurity spins: Here S 0 and S int are the action of the free electron system, i.e., the bath, and of the electron-impurity interaction (6), respectively; Z 0 is the statistical sum of the electron system without spin impurities.The functional integration is performed over fermionic (Grassmann) variables.The spin degrees of freedom in the action are described by one of the known approaches (e.g.coherent spin representation or a set of Majorana Grassmann variables, etc., see the book [42]).A particular choice of the spin variable is not important for our current purposes.
Let us start from the simplest case of the non-interacting fermions.The electron action S 0 in the Matsubara representation reads with σ = ±1 as in Eq. ( 4).The Matsubara Green's functions of the helical electrons in the momentum-frequency and space-frequency representation are given by The combined action reads where Calculating the Gaussian integral over the Grassmann variables, we obtain the contribution to the spin action, This expression is formally exact and, being properly regularized, describes all effects of the electron coupling to the spin impurities.Following the standard RKKY scheme, let us now focus on the weak coupling regime by restricting ourselves to terms up to second order in J in the action.This yields where X j ≡ (τ j , x j ).As argued below, other second order combinations do not contribute.After Fourier transform to the Matsubara frequencies, the first term in (15) takes the form with x ≡ x 2 − x 1 > 0 being the inter-impurity distance.The most interesting is the low temperature regime, x ≪ L T ≡ v/T , where the summation over frequencies in (16b) can be replaced by the integration over dω 2πT , resulting in Expressions similar to Eqs. ( 16) are governed also by the second term in (15).Combining all terms together, we obtain We are interested in slow motion of the MI spins with characteristic frequencies being much smaller than the inverse timeof-flight of the electron between the MIs, |Ω n | x v ≪ 1.In this case, e − x|Ωn | v ≃ 1 and, returning back to the imaginary time, we arrive at the expression This is the action of a system described by an effective RKKYlike Hamiltonian of the MI spins: with the RKKY energy scale This anisotropic spin coupling is ferromagnetic for all distances.Note that the RKKY coupling in a normal metal also starts out with a ferromagnetic sign at short distances [43].The first two expressions in Eq. (20b) are valid in the wide-band limit.In the presence of a finite but large bandwidth, this RKKY scale can also be rewritten as in the last expression where the RKKY coupling in units of the bandwidth (D) is simply given by the dimensionless j 2 0 [Eq.( 6)] p q X V 7 j Z F 5 e q z V 8 e W f z R 9 c i X / y y 9 z o f f T D r 7 / k D P w p D P 8 I W B F F 8 G H e A j K J w h M n A X 9 X e 8 c H H D 1 X v 9 t P Z r P / Z O r O 6 4 A q Y p M Z M i F / B t K E a B J O 8 9 Z L a 8 I q y K 5 r z i Y W K F t x M m 1 U w L X 5 h m T n O S m 2 X A r x i 7 0 4 0 t D B m W a S 2 s 6 C w M P 9 q H f k / b V J D d j R t h K p q 4 I q t D 8 p q i a H E X c p 4 L j R n I J c W U K a F v S t m C 2 q z B f s X 3 i q E 0 Z D E Z G j f P g o D 3 w 8 t i C N y S I 5 + h 3 A e D M h w E L y x a Z y g d f X Q L n q O 9 h F B M T p G r 9 E Z G i O G v q D v z o a z 6 X x z k b v l b q 9 b X e f n z D P 0 V 7 k 7 P w C X B b e 3 < / l a t e x i t > |ui = |""i < l a t e x i t s h a 1 _ b a s e 6 4 = " / e 3 t x 1 J q e D + p 9 k W M g q B 4 N S C + d 6 n J X Y n w q L S m q 4 r i W V g 1 L I S z G C n q d G 5 O D 6 0 / l 1 1 / T Q K 0 O a F d Y / g 3 S u P p y Y i t y 5 S Z 7 6 Z C 5 w 7 J 5 6 M / E 5 r 1 d h d t y f K l N W C E Y u F m W V p l j Q W V V 0 q C x I 1 B N P h L T K / 5 X K s b B C o i 9 0 U U K 7 x W P e 8 r e 3 w y Z j o S d x x I / 4 8 f 8 S z p o N 3 m o 0 P / o 2 3 p E F 1 s l r s k / e E E 5 i 0 i E f y C n p E k l u y T f y g / w M v g b f g 1 / B 7 0 V 0 K b i f e U U e I f j 7 D y K M p c Y = < / l a t e x i t > |di = |##i < l a t e x i t s h a 1 _ b a s e 6 4 = " t / N z m w 0 l p Q e S p h u P 0 9 O v c P z 4 D 6 1 R h v u O 0 h E E u x k Z l S g r 0 0 r D e v h g l V p i x B r p H a a I h w w u a j I p z I 6 w t z p 9 R q 8 Y T f E z X h v U d 1 m B R G L K I e t K M 4 t 1 4 T n g n C j u U N 9 g C O 9 3 t 5 M u P u + 7 0 c F i / 9 X f J K g e D U g v n + p y V O J g J i 0 p q u K w l l Y N S y F M x h r 6 n R u T g B r P F h p f 0 k 1 d G N C u s P w b p Q v 1 7 Y i Z y 5 6 Z 5 6 p O 5 w I n 7 1 5 u L / / P 6 F W b t w U y Z s k I w c v l Q V m m K B Z 3 X R U f K g k Q 9 9 U R I q / x f q Z w I K y T 6 U p c l d F o 8 5 i 2 / e y d s M h Z 6 E k d 8 l 7 e f S j h q N n i r 0 f z m 2 9 g n S 6 y T j 2 S b f C a c x K R L D s g h 6 R F J f p J r 8 p v 8 C a 6 C X 8 F N c L u M r g S P M 1 v k G Y L 7 B 2 1 o p 4 M = < / l a t e x i t > E FIG. 2. Eigenstates and eigenlevels of the RKKY Hamiltonian (20).Red lines and arrows show various decay channels of the singlet state (dashed lines) and up/down states (solid lines).These include spin-flips (outer paths) or phase flip (center path), as probed by the transverse ⟨S + L ∥S − R ⟩ω etc. and the longitudinal ⟨S z L ∥S z R ⟩ω dynamical spin correlation function, and which in the helical context reflect backward and forward scattering, respectively.divided by the distance of the impurities in units of the lattice spacing [Eq.( 2)].The wideband limit in the analytical approach therefore implies two assumptions: based on the second order approach used to derive Eq. (20b), this implies (i) j 0 ≡ ρ 0 J ≪ 1, i.e., J ≪ D. Yet via Eq.(2) [see also Eq. (D7)], the wideband limit also implies (ii) x ≫ a.
The ground state of this MI spin Hamiltonian is the triplet state with S z =0 (cf.Fig. 2): Hence despite the ferromagnetic coupling in Eq. (20a), the spins are antialigned.Yet by not being in a singlet state, they are still also exposed to Kondo screening.The forward-scattering (∼J z ) parts of the electron-MI coupling do not give a contribution ∼ S z 1 S z 2 to the secondorder effective spin action.This is because the contribution of the S z 1 S z 2 -term to Eqs.(15-16b) contains the product of two Green's functions of the same chirality, G σ (x, iω n )G σ (−x, iω n ), which vanishes in Eq.(16b) due to Eq. (11).Cross-scattering of the type S z S ± is also absent in the effective second-order spin action in Eq. ( 18) because, in equilibrium, the electrons of different chiralities are not correlated.Namely, the averaging of the corresponding combinations of fermion operators contain three operators of the same chirality, e.g.
)⟩ which vanishes due to the property ⟨Ψ L (X 1 )Ψ † R (X 2 )⟩ = 0.In addition to the ground state |t⟩, the effective spin Hamiltonian (20) possesses three other eigenstates with higher energies (cf.Fig. 2): the remaining triplet states denoted as the degenerate doublet ('up' and 'down') states, and the singlet state, The excitation energies of the doublet and the singlet relative to the ground state are E R and 2E R , respectively (cf.Fig. 2).
There is an important physical difference between the ground and all excited states of the effective spin Hamiltonian.The first one is stable and corresponds to the ground state of the total many-particle (electrons + spins) system projected onto the spin sector.To be specific, in the wideband limit D → ∞ (i.e., j 0 → 0) and zero temperature, the pair of impurities live in an exact non-correlated product ground state with the helical edge channel, i.e., |g⟩ ≡ |t⟩ ⊗ |0⟩ edge .In contrast, the excited states of the spins remain connected to the many-electron "reservoir" and as such cannot be simply written like product states for true eigenstates of the entire system.Instead, they represent resonant states, in the sense that they give rise to (narrow) resonances in the dynamical spin response functions.
By using the Fermi Golden rule, one can estimate decay rates of the excited spin states: the decay rate of the states |u⟩ and |d⟩ is ∝ J2 0 E R while the state |s⟩ has the parametrically smaller decay rate J2 z J4 0 E R (see [44] for more details).

B. Taking into account non-perturbative effects of electron interactions or finite Jz
Non-perturbative effects of the electron interactions and of a finite J z on the indirect exchange interaction of two MIs attached to the helical edge can be described by using the bosonized theory.The standard free (without impurities) action reads as: The single bosonic field ϕ describes both the spin and chiral degrees of freedom [16,21].K and u are the Luttinger parameter and speed of the helical plasmons, respectively.K incorporates effects of the electron interaction.The boson-spin exchange interaction is described by actions: × S + j e −2iϕ + c.c. .
Here, subscripts fs and bs stand for forward-/backward scattering, and ξ is the lattice constant which is usually needed to make the bosonized theory regular.By using the Emery-Kivelson gauge transformation [45], one can completely reduce the effect of J z to changing the dimension of the backscattering, which can be described by the effective dependence of K on J z : Below, we assume that Jz is included in K.
The bosonic theory is not quadratic and a functional integral over the bosonic fields D{ϕ}e −(S b +S bs ) (28) cannot be calculated exactly, even formally.However, the effective action of the indirect spin interaction can be obtained for small J by calculating the integral over ϕ as the first cumulant.This is similar to the renormalization group (RG) treatment of the sine-Gordon theory [40]: one Taylor expands the exponential in Eq.( 28) in S bs up to the second order, calculates the integral over ϕ and re-exponentiates the answer: where ⟨A⟩ S b ≡ D{ϕ} Ae −S b .Using the well-known expression for the bosonic correlation function we arrive at the perturbative (in J) expression for the effective spin action This theory is non-local in time and, thus, takes into account retardation effects which are beyond the Hamiltonian formulation of the RKKY theory.The action can be reduced to a local one if ξ ≪ |x 1 − x 2 | ≪ L T and 1/2 < K ≤ 1.In this case, Eq. ( 30) reduces to These expressions were analyzed in Ref. [1].They coincide with answers derived in the previous section for the noninteracting case, K = 1.
The application of the bosonized theory has several advantages.Firstly, it takes into account non-perturbative effects of the electron interaction and of J z .Besides, it is straightforward to go beyond the quadratic approximation in J and derive renormalization of this Kondo coupling constant [16,21].

IV. NRG ANALYSIS OF THE 2HKM
We proceed by presenting the NRG results for the 2HKM in this section, before we explain in detail how we setup the two-impurity helical setup within the NRG framework in the subsequent Sec. A. An exhaustive NRG parameter scan for the 2HKM at J = 0.1 is shown in Fig. 3.All panels have the anisotropy parameter J z /J on their vertical axis whereas the horizontal axis shows energy in various forms [temperature [first column, i.e., Figs.3(:,1)], NRG energy scale [second column Figs.3(:,2)], or frequency (remainder of columns, Figs.3(:,3-5))].Throughout, energy is decreasing towards the right, as motivated by the NRG approach (second column), where large energy scales come first, followed by a zoom into exponentially small energy scales towards the right.Each row labelled by a letter shows data for a fixed impurity distance x as indicated in the left panel.This distance is always chosen integer, i.e., on the grid (A11), and increases exponentially towards lower panels, with the value specified in the left panels.
The inverse time to travel in between the impurities naturally gives rise to its own energy scale which may be interpreted as the 'Thouless energy' of a quantum dot confined by the magnetic impurities.This energy scale in itself is independent of any impurity properties other than their distance [see also Eq. ( A9)].Now since the impurity distance explicitly enters the construction of the Wilson ladder (A34), there is always a clear qualitative change in the NRG finite size spectra, aka.energy flow diagram at the energy scale E x , as visualized in a condensed graphic way in Fig. 3(:,2), i.e., the second column.There one observes a distinct change (whitish to gray transition) right at E x (vertical dotted line).This 'curtain' is opened, i.e., moves towards the right as the impurity distance is increased from the top to the bottom panels.As such this 'unveils' the underlying Kondo physics.At largest distance shown in the bottom row of panels, Fig. 3(h,:), E R has already dropped below the smallest Hence in this case the impurities are fully Kondo screened individually, and RKKY physics plays no role any longer, and hence is absent for all J z ≥ 0.
At larger energies above the coherence scale, E ≳ E x , i.e., to the left of the vertical dotted black line, the system is described by effectively independent impurities.Since information cannot travel faster than v between the impurities, the impurities do not yet 'see' each other at energy scales E > E x .In this sense, the finite size spectra in the NRG look identical in Fig. 3(:,2) for large energies E ≳ E x , i.e., to the left of vertical dotted line across all panels in Figs.3(:,2).
For low energies, Fig. 3 shows that the smallest energy scale from the point of the impurity is max(E R , T K ).Therefore E R serves as a low-energy cutoff.For example, the energy flow diagram is converged below E R [uniform gray area to the right of E R (blue vertical marker) in Figs.3(:,2)].On explicit physical grounds this is replicated in terms of static inter-impurity correlations for T < E x in Figs.3(:,1), or by having no particular structure in the spectral density for |ω| < E R in the right panels.For this reason the brown line which depicts the Kondo scale, is only shown above E R and hence terminated at the blue vertical marker line.For J z above this crossing point, Kondo screening sets in and eventually fully dominates.This is seen in the brightening of the dark red inter-impurity spin correlations in Fig. 3(:,1) towards white (no correlations) when increasing J z (vertical direction).(for precise prefactors, see [46]).This encodes cumulative effective thermal weights in three arbitrary but fixed low-energy symmetry sectors into red-green-blue (RGB) colors for all odd Wilson shells n (to avoid even/odd effects) based on an effective inverse temperature βn ≡ 4ωn.The symmetry sectors chosen for this were q ≡ (Q; S tot z ) ∈ (0; 0, 1, −1) with Q the total charge relative to half-filling.The remainder of the columns shows dynamical spin-spin correlation functions ⟨ Ŝη ∥ Ŝη ′ † ⟩ω as indicated at the top of each column at zero temperature (T = 1.8 10 −11 ).The solid-dotted lines shows the respective derived inverse static susceptibility T ηη ′ S = 1/4χ ηη ′ S [cf.Eq. ( 34)].Blue (red) indicates positive (negative) value, respectively.The number at the top right of each panel indicates the maximum absolute value the spectral data A of the broadened NRG spectral data, and hence gives an impression of numerical range.Blue (red) shading in all panels except for the second column indicates positive (negative) values, respectively.The dashed lines in (a,c;3) and (a,c;5) represent exponential fits as indicated with (a5) of the maximum of the spectral data for J/Jz ≤ 0.5.NRG parameters (e.g., see [46] for detailed definitions): Λ = 4, truncation energy Etrunc = 8; z-averaged over nz = 4, with log-Gauss broadening σ = 0.3 of the spectral data after z-averaging.2)], the dimensionless Kondo coupling in Fig. 3 is small, having j 0 ≡ ρ 0 J = 0.05 ≪ 1, and thus Ex ER = 1 πj 2 0 = 127.3.Therefore the bare RKKY regime spans about two orders of magnitude in energy scale.It competes with Kondo physics for the case T K (J, J z ) < E x , i.e., when the Kondo scale is small enough that the screening clouds of the two impurities overlap.We refer to this intersect as the intermediate energy regime.

The bare RKKY regime concerns the energy window
This intermediate regime becomes visible as the lighterblue shaded area in Fig. 3(a-f,2) in between the two vertical makers below the brown solid line which represents the single-impurity Kondo scale T K .The latter scale T K is cut off at E R , also visually so by terminating its line towards the right at the blue vertical marker.Hence the intermediate regime is present (i) if there is a brown line segment in between the two vertical markers (black dashed and blue), in which case (ii) the intermediate regime occurs below it.For example, with no brown line segment in the bare RKKY regime in the lowest two rows in Fig. 3, T K > E x dominates the low-energy regime throughout, and the system consists of two individually Kondo screened impurities.However, once T K < E x , the buildup of the Kondo screening cloud is affected by the presence of the other impurity.This leads to characteristic changes in the NRG energy flow diagrams, as seen in Fig. 3(:,2).

A. Spin-spin correlation functions and low-energy scales
Data for dynamical spin-spin correlation function A(ω) [as defined in Eq. (A53) in App.A 4] is shown as color-plots in the right panels of Fig. 3.The frequency of the maximum for fixed parameters sets the relevant low-energy scale for the spins.It is well traced by the inverse static susceptibility (solid-dotted lines) as derived from the dynamical susceptibility with S ∈ {S z , S ± }, and the normalization convention for In case that η = η ′ , i.e., inter-impurity correlations, only a single label may be shown, e.g., (the η label may be skipped altogether then, since the impurities are considered identical, hence by symmetry, e.g., . Strictly speaking, the interpretation as an energy scale is only justified for 'diagonal' correlations, (here intra-impurity η = η ′ ), as this guarantees positive spectral data and hence a respective positive energy scale.But it is useful to also include η ̸ = η ′ here for the sake of the argument and presentation.We will also refer to χ ηη ′ S± which involves a spin-flip, as the transverse susceptibility, and χ ηη ′ Sz which involves a phase flip, as the longitudinal susceptibility.The choice of the prefactor (1/4) is motivated by the standard definition of the Kondo temperature in the plain single-impurity Kondo model, T NRG K ≡ 1 4χ0 [33,47].For isolated RKKY impurities with the energy spectrum as in Fig. 2, the spectral functions have δ-peaks at energies ω = ±E R for the transverse, and ±2E R for the longitudinal correlation function, thus resulting in T ηη ′ ;0 S± ≡ E R /2 and T ηη ′ ;0 Sz ≡ sgn(ηη ′ )E R , respectively [note the normalization convention of the spin operators as indicated with Eq. ( 35)].Up to a sign, these are independent of the choice of η and η ′ , i.e., they become the same for inter-and intra-impurity susceptibilities.The latter is a direct consequence of the RKKY low-energy regime where the pair of impurities, even though spatially separated, act like a nearly-decouled microscopic unit governed by the RKKY Hamiltonian.In the following we will thus scale energies in the numerical data by the smaller energy scale (where subscript 'S' denotes 'spin'), We prefer T S over E R , since E R only represents a lowest order estimate, whereas T S includes the full many-body aspects of the problem and is thus also self-contained and thus consistent within the NRG.
The energy scales T ηη ′ S in Eq. ( 34) capture the low-energy scale of the impurity spins, as seen, e.g., in the lowest panels Figs.3(f-h,3), [with the impurity operators S as well as their location η and η ′ specified at the top of each of the right columns Figs.3(:,3-5)].There the inverse susceptibility T S± from the intra-impurity spin-spin correlation (soliddotted line) [Figs.3(:,3)] follows closely the analytical Kondo scale T K (J, J z ) (brown line) up to a constant prefactor of order one.This generally holds for ω > E x , i.e., to the left of the black dashed line in Figs.3(:,3).
The two right-most columns of Fig. 3, in contrast, show inter-impurity correlations.These can only be due to RKKY interactions, and hence diminish with increasing impurity distance.Once this distance exceeds the size scale of the Kondo screening cloud, i.e., T K > E x , the inter-impurity susceptibility becomes much smaller than the onsite susceptibility, such that its inverse is orders of magnitude larger in energy [e.g.compare solid-dotted line in Figs.3(gh;4,5) to T K (brown line)].Its physical interpretation is that inter-impurity correlations start to play a role relative to Kondo correlations only once the latter are sufficiently suppressed, e.g., by a large temperature scale T ≳ T LR S± ≫ T K .In the intermediate regime, where the low-energy physics is cut off by RKKY, the inter-and intra-spin correlations start to look identical when applying the same operators [e.g., compare Fig. 3(a,3) to Fig. 3(a,4)].This holds quantitatively as also seen by the overall scale (see maximum spectral weight A indicated with the panels).This can be understood based on the RKKY impurity state |t⟩ in Eq. ( 21) that (nearly) decouples as a product state from the bath channel [cf.App.A 7], thus having which holds for both, η = η ′ (intra-impurity) as well as η ̸ = η ′ (inter-impurity), while bearing in mind that the equal-time correlator above is identical to the integrated spectral data over frequency [cf.sum rules].
The longitudinal inter-impurity correlations are shown in the last column of Fig. 3. Deep in the RKKY, the dominant spectral weight of this non-diagonal dynamical correlation functions is expected to be negative, and half the absolute value, in agreement with the red shading (which indicates negative) and overall scale of the spectral data, e.g., in Fig. 3(a5) [see a detailed analysis of the effects of finite bandwidth on the precise value of the l.h.s. in Eq. (38) in App.A 7, and in particular Fig. 11 therein].Via spectral sum rule, the frequency integrated data yields the data in Fig. 3(:,1) at given temperature (in the present case, at the lowest temperature shown).Consistently, this appears in a deep red in the RKKY state, indicating that the spin-orientations of the two impurities are antiferromagnetically (AF) correlated.However, when crossing over into the Kondo regime, the weakened inter-impurity longitudinal correlations can even change sign and turn ferromagnetic (FM) [e.g., blue soliddotted lines in Figs.3(:,5) in the inverse susceptibility; this necessarily has to originate from corresponding positive spectral data, but its weight is too small, though, so it is not visible in shading in the spectral data as shown].The weak ferromagnetic correlation can also be partly seen in the static equaltime spin correlation [e.g., the very faint blue hue in the intermediate regime just left of the red area in Fig. 3(b,1) or (e,1), and hence at significantly elevated temperatures].The sign change towards weak ferromagnetic correlation observed in the longitudinal data is not systematic, though.For example, not all inverse susceptibilities in Figs.3(:,5) show a sign change.Furthermore, deep in the Kondo regime, the system can feature weak ferromagnetic (FM) inter-impurity spin correlations and together with a sign change depending on the frequency [e.g., faint red and blue shadings in Figs.3(g-h,5) just above J z = 0 vs. ω, which is partly also visibly replicated vs. temperature e.g., in Figs.3(g,1)].For energies below E x this should be well resolved by NRG, and not related to coherence affects averaged out by NRG above E x due to coarse graining [cf.discussion following Eq.(A9)].The appearance of weak AF as well as FM correlations across the impurities may be related to the fact, that for finite J z and finite bandwidth, subleading terms can generate an effective small longitudinal inter-impurity interaction J z ŜL z ŜR z which due to its oscillatory behavior vs. distance may be ferromagnetic as well as antiferromagnetic [cf.App.A 7, and also the discussion around Eq. (D11)].Consequently then, z-averaging of the spectral data may lead to apparent non-systematic behavior in the longitudinal correlations as a numerical artifact.In the present case, however, we do not dwell on this any longer as this concerns subleading effects.

B. RKKY scale in spectral data
The inverse transverse susceptibility deep in the RKKY regime resembles a straight line on a semilog plot [Fig.3(a, [3][4]].This also reflects the behavior of the maximum in the actual spectral data.Generally, from its very definition via the Kramers-Kronig integral relations, the inverse susceptibility is also sensitive to the precise line shape of the spectral data.Yet by having the maxima and inverse susceptibility run in parallel for small J/J z , this suggests similar line shapes.Tracking and fitting the peak maxima in the spectral data by the exponential fit specified with Fig. 3(a5) for J/J z ≤ 0.5, we obtain ω * ≃ 2.78•10 −4 e 0.82Jz (blue dashed line).It is lower-bound at J z = 0 by the analytically obtained RKKY scale E R [Eq.(20); blue vertical marker in Fig. 3], having E R = (ρ 0 J) 2 /x = 2.50•10 −4 for x = 10 [Fig.3(a,:)].The difference of about 10% is due to finite bandwidth [cf.App.A 7].
As seen from the fits in Fig. 3(a;3,5) the low-energy scale in the intermediate regime diminishes exponentially with decreasing J z , as in T S ∼ ω * ≃ ae bJz with b > 0. This is qualitatively similar to the one-impurity Kondo temperature in the anisotropic 1-impurity case (cf.App.B) yet with different renormalized coefficients.For one, this shows that the RG / poor-man's scaling for a single impurity needs to be stopped at the RKKY scale.Moreover, the relative slope in the exponents of the low-energy scale as seen in the semilog plots in Figs.3(:,4) changes in the intermediate regime with increasing impurity distance.The slope 1/b of the peak vs. ω * is larger than for T K (brown line) in Fig. 3(a;3,4), about comparable in Fig. 3(e,4), and smaller in Fig. 3(f-h;4) where RKKY is absent.This clearly underlines the continuous crossover from RKKY to Kondo.
A similar exponential fit on the maximum of the spectral data was also obtained for the longitudinal spin correlations in Fig. 3(a5).The peak in the spectral data occurs at a larger a ′ = 5.86•10 −4 at J z = 0 when compared to the corresponding fit in Fig. 3(a3).The slopes b are comparable, though (b ′ = 0.80 vs. b = 0.82).Therefore the maximum in the longitudical spectral data ⟨ ŜL z ∥ ŜR z ⟩ ω is systematically shifted at small J z by a factor of 5.86/2.78= 2.11 towards larger energies as compared to the maximum in the spectral data that requires a spin-flip, ⟨ ŜL + ∥ ŜR − ⟩ ω .Thus in the RKKY regime at x = 10 [Figs.3(a,:)], and consistent with the explicit analytical expression for E R in Eq. ( 20) This is in agreement with the effective twoimpurity level-spectrum in Fig. 2, where with reference to the Lehmann representation for spectral data in Eq. (A53), one needs to pay an energy E R for a spin-flip, whereas one needs to pay an energy 2E R for a sign-flip (triplet |t⟩ to singlet |s⟩ transition).The deviation of about 5% from the expected factor of 2 in the fitted values is within the spectral resolution of NRG, and thus likely due to z-averaging and broadening of the spectral data.The parameter scans in Fig. 3 give important hints: they show that E R (blue vertical marker in Fig. 3) as obtained from a second order perturbative approach [cf.Eq. (20b)], does not always describe the RKKY low-energy scale correctly, as it can get renormalized by the presence of single-impurity Kondo correlations.While for x = 10 this well coincides with the J z = 0 low-energy scale [e.g., compare to solid lines with symbols in Fig. 3(a,4)], when increasing the distance, the peak in the NRG data shifts towards the left of the vertical blue marker, i.e., towards values that are larger than E R .That is, the fit value for a shown with Fig. 3(a;3,5) effectively increases relative to E R when exponentially increasing the distance x (no additional fits shown, though, as this is a qualitative argument).Thus even if at J z = 0 the oneimpurity Kondo scale may still be many orders of magnitude smaller than E R , the value of E R already gets weakly affected (likely acquires logarithmic corrections) due to the underlying Kondo correlations, even if T K ≪ E R .If at the same time it also holds E R ≪ D = 1, then integrating out the helical edge mode starting from the initial band edge D towards zero energy, e.g., in a poor-man's scaling sense, this can be expected to introduce an RG flow also for E R .This suggests that for a consistent interpretation of the NRG data with analytics, the NRG data should be scaled by T (η) Sz [cf.Eq. ( 34)], rather than the lowest order analytical estimate E R in Eq. (20b).
Moreover, the peaks seen in the spectral data are rather narrow, i.e., within the resolution limit of the presented NRG data.Based on the chosen NRG parameters Λ = 4 with n z = 4 z-shifts (e.g., see Sec.A 3, or also [46] for more detailed definitions), we expect a best possible relative spectral resolution δω/ω > Λ 1/2nz − 1 = 0.19, hence the value used for the broadening of σ = 0.3 on the discrete NRG raw data.Based on perturbative approach, however, one may suspect significantly narrower features as shown in Fig. 3.With this in mind, the spectral features seen in data in Fig. 3 are likely overbroadened.

C. More detailed spectral analysis (line shapes)
The RKKY 'resonances' in the spectral data are expected to be potentially very narrow.As it turns out, though, while sharp peaks (near delta spikes) can be found in the bare discrete spectral NRG data, the precise location in frequency of these peaks is sensitive on the z-shift in the logarithmic discretization of the NRG setup (e.g., see Sec.A 3 or also [46] for detailed definitions).Specifically, as the z-shifts can shift energies by a factor of √ Λ in the discrete setup, in the present case, also the respective 'response' of the system in terms of the precise location of a narrow spectral resonance for the 2HKM model can also vary within a factor of √ Λ. Therefore blind z-averaging of the NRG data leads to artificial broadening and somewhat irregular z-averaged data, as demonstrated in Fig. 4 or Fig. 6 for two different impurity distances, x = 10 and x = 1 000, respectively.Instead, when scaling the spectral data for each individual z-shift by its respective T Sα (z) and averaging that resulting data, peak shapes are significantly improved, and in particular also narrower (e.g., compare Fig. 4 → 5, or Fig. 6 → 7).
For the analysis here, the rather large Λ = 4 is useful to emphasize how discretization manifests itself in the spectral data which is more subtle here, as it also leads to shifts.A major motivation for the larger Λ = 4 is, of course, that this FIG. 5. Exactly the same bare data as in Fig. 4, except that for each z-shift the frequency scale is first scaled by T ηη ′ Sα /T ηη ′ Sα (z) (with T ηη ′ Sα ≡ ⟨T ηη ′ Sα (z)⟩z the z-averaged value, and α ∈ {z, ±} chosen, respectively, for each curve), and then combined as in Fig. 4. The zshift specific TS z (z) is also determined within the NRG in the same calculation as part of the post analysis [cf.Eq. ( 34)].Overall, this procedure leads to a significantly improved quality of peak shape that is significantly narrower, as compared to the spurious spread of spectral peaks in the right panels in Fig. 4. also results in faster NRG calculation or, conversely, in better converged NRG spectral data.By careful z-averaging, we can obtain good spectral resolution for Λ = 4, nevertheless.However, given that RKKY peaks can be expected to become much narrower in terms of width vs. location, this ultimately will be also challenging for smaller Λ in any case.Overall, the FIG. 6.
Same analysis as in Fig. 6 except for an increased impurity distance x = 1 000, leading to the reduced value for TS z = 3.0•10 −6 = 1.19 ER.The broadening was also increase to σ = 0.15.FIG. 7. Exactly the same bare data as in Fig. 6, except that for each z-shift the frequency scale is first scaled by T ηη ′ Sα /T ηη ′ Sα (z) (similar to Fig. 5).
present analysis in terms of Λ = 4 already clearly supports narrow features.More importantly, the location of the peaks, and hence the corresponding energy scale, can be considered reliable and significantly more accurate than the width, given also that these are directly related to inverse static susceptibilities.
The data in the right panels of Fig. 5 gives a good impression for the spectral data in the RKKY regime.As compared to Eq. ( 36) with T S = 0.50 E R , the actual data in Figs. 4  & 5 gives T S = 0.556 E R which is reasonably close.Furthermore, the peak location in the spectral data is expected at ω * = E R ≃ 2T S for the transverse spectral data, and at ω * = 2E R ≃ 4T S in the longitudinal data, again consistent with the data in Fig. 5.The peak width in the longitudinal data related to a phase flip when transitioning to the singlet state [cf.Fig. 2] is found to be comparable as in the transverse data [44].However, this width is limited by the broadening σ as seen with the individual data in the left panels.While some structure can be observed with Fig. 5(a) resolved by the spread with z-shifts, the data in Fig. 5(b) is more smooth this way, just showing the broadening σ = 0.10 used.This suggests that the data in (b) is still likely overbroadened by given σ = 0.10.Finally, as already argued with Fig. 3, the inter-and intra-impurity correlations are identical to each other deep in the RKKY regime.Here this is seen by having the dashed lines in the right panels of Fig. 5 on top of the solid lines [note the sign change, though, as indicated with legend in Fig. 5(d)], which again is rooted in the triplet ground state in Eq. (21).
Repeating the same analysis for the increased impurity distance x = 1 → 1 000, the data in Fig. 7 starts to show a qualitative change in the spectral line shape.While having increased the broadening to σ = 0.15, as reflected by the individual peaks in the left panels, the overall lineshape in Fig. 7(d) is still largely comparable.The interand intra-impurity correlations do lie nearly exactly on top of each other, thus also suggesting a clear RKKY low-energy state still.In particular, peak position are still at the expected ω * = 2T Sz or 4T Sz for the transverse or longitudinal data.The numerical value for T Sz , however, further deviates from the plain second-order perturbative RKKY scale, having T Sz = 3.0•10 −6 = 1.19 E R which by now is clearly unequal from the naive expected value of 0.5 E R (see also Fig. 8 for a more detailed analysis in this regard).Yet when scaling the data consistently fully within the NRG framework, the schematic picture in Fig. 2 still works well when substituting To conclude this section, we reemphasize that the NRG approach correctly describes the position of the peaks in the response functions as well as the integrated spectral weight, but yields resolution-limited information about the peak width and thus its overall shape.The latter unavoidably results from the coarse-graining in energy space intrinsic to the NRG approach.

D. RKKY vs. Kondo Energy scales
An analysis of the low-energy scales vs. impurity distance is shown in Fig. 8.The longitudinal inter-impurity susceptibilities are negative due to the antiferromagnetic RKKY correlations; note the sign with T LR Sz in the legend.In the crossover region when leaving the RKKY regime, this can become positive, though [cf.legend with Fig. 8(c)].As seen in all panels, the energy scale T S 'drifts' away from the bare RKKY scale E R (black dashed line) already many orders above the Kondo scale T K (indicated or specified with the upper panels).If the distance exceeds the temperature length scale v/T (left panels) or the inverse Kondo scale (right panels), the impurities effectively become independent of each other, such that RKKY is cutoff by max(T, T K ).In Fig. 8, the end of the RKKY regime is seen where inter-impurity correlations start to deviate significantly from their intra-impurity counter part (e.g., compare blue vs. yellow in the lower panels).By approaching the wide-band limit (or equivalently, by reducing J, as in right vs. left panels in Fig. 8), one can observe in the lower panels that T Sz ≃ 2(T S± ≡ T S ) is obeyed over a wide distance (energy) window.Also deviations from the relative factor of 2 diminish towards the wide-band limit.The shad- ing in the lower panels of Fig. 8 indicates variations (i.e., the standard deviation) due to by z-shifts in the logarithmic discretization.These also become smaller towards the wide-band limit [e.g., compare Fig. 8(d → b)].

V. SUMMARY AND OUTLOOK
We presented a broad parameter scan of the two-impurity helical Kondo model (2HKM) vs. impurity distance, coupling anisotropy, temperature, and finite bandwidth within the NRG framework.We emphasize that this setup substantially differs from a chiral edge mode, where both spins propagate in the same direction [for a recent NRG study on the latter, see e.g., Lotem et al. [41]; the chiral model is different from the helical edge mode discussed here, both in the setup as well as in the physics].With the NRG being non-perturbative in character, our presented results are reliable quantitatively in the full parameter regime, with the only real constraint being energy resolution in spectral data.We have established a plain crossover from RKKY to Kondo with increasing impurity distance or, conversely, increasing Kondo coupling or its anisotropy.The Kondo screening of the impurities individually by the helical edge mode is tuned continuously into an effective mutual 'RKKY' screening of the impurities themselves.In this sense, the Kondo renormalization flow is cutoff by the RKKY energy E R once it exceeds the Kondo scale.The RKKY occurs in the energy window determined by the inverse time scale required for the impurity to travel in between the impurities (Thouless energy E x ), and the actual RRKY energy E R ≪ E x .If the Kondo scale falls within this window one comes across a continuous crossover.
The low-energy effective RKKY Hamiltonian gives rise to narrow resonances with the helical edge.While the presented NRG analysis, the location and overall spectral weight is reliable, the spectral width, however, is likely (much) below the energy resolution of our NRG data deep within the RKKY far from any Kondo screening.Hence the precise linewidth e.g., in the dynamical spin-spin correlation data for the fully interacting model is left for future analytical and numerical studies.

Hybridization function
The helical 1D edge mode above is assumed to constitute a fermionic macroscopic bath that interacts equally with two impurities dησ symmetrically located at positions x η = η 2 x with η = ±1 [Eq.( 5)].The NRG approach first focuses on Anderson type impurities with explicit hybridization of the impurities with the bath.The switch to Kondo-type impurities can be taken as a subsequent step, e.g., via Schrieffer-Wolff transformation.Thus the starting point is the hybridization Hamiltonian, This defines the bath states f0ησ 2 ĉkσ that the impurities couple to, which permits the interpretation that these are the bath states at the respective location of the impurities.The states f0ησ are normalized, but not necessarily orthogonal yet, hence the tilde (cf.Sec.A 2).
In Eq. (A2), an electron with spin σ can hop on and off the corresponding helical branch, with the spin σ preserved in the process.The hopping amplitude V k is assumed independent and thus symmetric in the spin and impurity indices.For simplicity, the bath is also assumed featureless, characterized by two parameters only: a hybridization strength Γ and a finite half-bandwidth D. In the continuum limit, the hybridization function for each impurity individually is given by Γ(ε) = πρ ε |V ε | 2 ≡ Γ ϑ(D − |ε|) (aka.the default NRG box distribution), with ρ ε = ρ 0 ϑ(D − |ε|) the one-particle density of states, assuming constant ρ 0 = 1 2D without restricting case.The hybridization is cut off sharply in energy, as depicted in Fig. 1.The integrated (norm-squared) hybridization strength which yields the splitoff normalization factor in the second line in Eq. (A2).
From the perspective of the impurities, the full hybridization function becomes a 2 × 2 matrix with η ∈ {R, L} indexing the impurities [Eq.( 5)].This includes a non-zero off-diagonal hybridization for η ̸ = η ′ .In the wide-band limit D → ∞, the hybridization becomes where τ ≡ σx v is the time required for a particle to travel from one impurity to the other, since for example from Eq. (A3) for a particle to travel from L → R (thus with x ̸ = 0, also τ ̸ = 0), This gives a non-zero contribution only for τ > 0. Assuming v, x > 0, for example, this requires σ > 0. Indeed, in given helical setting, a spin up electron travels to the right, that is the particle needs to be created at the left impurity (η ′ = −1) which then can propagate to the right impurity (η = +1).For spin down, the situation is vice versa.
We emphasize that the spectral representation Γ(ω) of the hybridization function cannot be simply written as −Im∆(ω) in the present case because the matrix elements in Eq. (A3) are complex.Instead, the imaginary part needs to be taken from the propagator only, i.e., with − Im where the spectral function Γ ηη ′ (ω) ≡ [Γ(ω)] ηη ′ needs to be differentiated from the constant Γ.Correspondingly, in matrix notation, which is non-zero now for both off-diagonal entries η ̸ = η ′ for either spin σ.Only by taking the spectral data as above, this can be simply completed to the full hybridization function using Kramers-Kronig transform on the complex spectral data, i.e., by folding the above possibly complex spectral function with 1/(ω + − ε).Within the wide-band limit and in the absence of interactions, the Green's function for the impurities with onsite energies ε dη becomes, e.g., in the spin-up channel (σ = +1, and hence τ > 0), At particle hole symmetry with ε dη = 0 on obtains for ω = 0, where the diagonal entries reflect half-filling based on the Friedel sum rule.Due to the helicity, the diagonal terms in Eq. (A7) are just the Green's function of decoupled impurities, G ηη (ω) = 1 ω−ε dη +iΓ .The non-zero off-diagonal term maintains the same matrix position as in ∆(ω) in Eq. (A4a), consistent with the directedness of propagation.Similar to Eq. (A4a), the off-diagonal term also shows oscillatory behavior in energy with period This period is fixed by the energy scale E x = v/x [Eq.( 33)] set by the inverse time τ required for an electron to travel from one impurity to the other.In particular, the period δω is independent of the energy scale ω in G(ω).Given that the NRG discretizes logarithmically in energy space, assuming E x ≪ D = 1, there is no way these oscillations can be resolved to orders of magnitude higher in energy all the way up to D. On the other hand, on physical grounds at high energies, e.g., at temperatures T > E x , the impurities effectively no longer see each other.In this sense it appears reasonable to expect for most physical quantities such as spin-spin interactions that do not explicitly resolve one-particle phases of propagation.In the presence of relaxation processes due to interactions, these rapid oscillations likely average out at large frequencies even for small temperature, and thus become less important in detail.Yet this needs to be verified in practice on a case-by-case basis by tracking the stability of the data with respect to the level of course graining.

Normalization of bath states and effects of finite bandwidth
For either spin, the hybridization term Eq. (A2) in the Hamiltonian defines two bath states f0ησ with η ∈ {R, L}.These are normalized, but not strictly orthogonal, hence the tildes with the f operators in Eq. (A2) as a reminder.These bath states at the location of the two impurities have finite overlap that can be determined from the fermionic anticommutator relation, with |0⟩ the vacuum state, and a ≡ πv D as in Eq. ( 2), and therefore πx a = Dx v = |τ D|.Rewritten in matrix notation it can be diagonalized, with eigenvalue matrix with eigenvalues sorted as s 0η = 1+ηr 0 , with η ∈ {+1, −1}, i.e., with the larger eigenvalue s 0+ ≥ 1 coming first.The off-diagonal term in Eq. (A10b), r 0 ≡ a δ a (x), represents a δ-function of width a.In this sense, the UV cut-off D again directly gives rise to the lattice constant already introduced with Eq. ( 2).It has unit of length and resembles the resolution in real space based on the given one-particle density of states encoded into the box distribution.By having a finite bandwidth, this translates into a cutoff in spatial resolution, which in the present case naturally gives rise to a well-defined lattice spacing: the overlap in Eq. (A10) between the two bath states becomes exactly zero for when n ̸ = 0 which naturally suggests a discrete grid with lattice spacing a.As pointed out with Eq. ( 2), choosing D = 1 and a = 1 then fixes the velocity to v = 1/π.For n = 0, i.e., x → 0, the off-diagonal overlap in Eq. (A10a) becomes 1, and thus identical to the diagonal case η = η ′ .In this case the two locations η and η ′ approach the same 'site' and thus become identical.
As an aside, we note that the argument here namely that a finite bandwidth naturally gives rise to a discrete lattice spacing, can also be straightforwardly carried over to a plain spinless tight binding chain, with the minor difference that there due to the structured density of states, the Fermi velocity v differs from the value above by a factor of order 1.Furthermore, the Fourier transform of any quadratic Hamiltonian in momentum space yields its full (long-range) hopping structure in real space.Though a purely 1D lattice model is unable to describe the topologically nontrivial phase of the 2D topological insulator, the one-particle dispersion may nevertheless also be Fourier transformed for the isolated single helical edge as in Eq. (1).By starting from momentum space, however, in the present case this would mandate periodic (or infinite) BCs with spin-dependent complex hopping coefficients t nσ ≡ σt n , where with ε kσ ≡ σε k , for a chain of N sites, with distance x = x n = na > 0, having t 0 = 0.These hopping amplitudes are long range, decaying like 1/x, and by construction break inversion symmetry, t n = t * −n ̸ = t −n .The long-range hoppings cannot be eliminated by permitting deviations from the linear dispersion close to the band edge.For example, the altered dispersion εk = 2v a sin ka 2 which still obeys εk ≃ ε k = vk for small k, also results in long-range hopping.This is due to the fact that the one-particle dispersion is discontinuous across the boundary of the Brillouin zone.Given this complications for numerical lattice simulations, the starting point in energymomentum space in Eq. ( 1) is more natural and convenient.Nevertheless, from the above it is clear that the notion of lattice spacing and Brillouin zone are perfectly valid also for a single helical edge mode even if the one-particle dispersion is discontinuous across the boundary of the Brillouin zone.
The analytic structure of the overlap in Eq. (A10) is closely related to static fermionic correlations versus distance.For example, consider a filled helical Fermi sea for energies ε ∈ [−D, 0].Then at zero temperature with τ = σx v [Eq.(A4a)], the integral only includes the filled Fermi sea [whereas Eq. (A10a) integrated over the entire 'Brillouine zone'].By comparison with the corresponding Fermionic correlations for plain 1D tight-binding chain, this suggests the Fermi wave vector k f as half the extent of the filled Fermi sea in momentum space, i.e., given half-filling, The leading phase factor in Eq. (A13) has subtle consequences when computing charge correlations, and results in features that qualitatively differ from a plain tight binding chain.Assuming a continuous, i.e., non-discretized 1D edge mode in terms of lattice sites spaced by a, then based on Eq. (A10), the overlap diminishes to zero for x ≫ a.In particular, this includes the wide band limit where a → 0. However, for the sake of orthogonality of the fermionic states, the vanishing of the overlap at finite bandwidth can be simply guaranteed by adhering to the discrete grid in Eq. (A11) with lattice spacing a in complete analogy to a tight-binding chain.Hence, in order to avoid complications based on non-orthogonal f0ησ states, henceforth distances will be chosen on the grid (A11), i.e., with a := 1 having x ∈ Z.With this f0ησ → f0ησ become well-defined orthonormalized local bath modes at the location of the impurities, denoted by using hats now instead of tildes.
While finite bandwidth is physically meaningful when having particular 2D lattice models in mind, for a helical edge mode this cutoff is peculiar in that the helical branches merge with a continuum of bulk states.Therefore a sharp ultraviolet cutoff for an isolated 1D helical edge mode can have potentially artificial consequences.Lack of orthogonality of local bath modes discussed with Eq. (A10b) above is one example.The latter complication can be simply eliminated, though, by adhering to the effective discrete lattice in Eq. (A11).On a related footing, the hybridization function in Eq. (A4a), which is closely related to the dynamical one-particle propagation in between the impurities, reflects the directedness of motion via the step functions θ(τ ).This step function, however, is strict for infinite bandwidth only.For finite bandwidth it also contains, in particular, a non-zero oscillatory tail for τ < 0. That is, for |ω| ≪ D and τ < 0, rather than strictly being zero, the amplitude for this enhanced backscattering probability decays like 1 D e i(ω±D)τ ∼ 1 D e iωτ [see also Eq. (A6)] where the oscillatory behavior with phase Dτ = πx a is similar to Eq. (A10b), thus having e iτ D = ±1 on the grid Eq. (A11).However, this backscattering probability decays with increasing bandwidth D, which eventually enforces strict directionality.Yet for the above reasons finite bandwidth can generate a weak subleading contribution J z S z L S z R to the RKKY Hamiltonian (20) in the helical system [cf.App.D].

Coarse graining
For the sake of a numerical treatment, the continuum of the bath needs to be discretized.Here we use the NRG which, by construction, always discretizes in energy-momentum space.This allows us to target a single edge mode with plain linear dispersion.To be specific, the NRG coarse-grains on a logarithmic grid D/Λ −(n+z) in energy space with Λ > 1 (typically Λ ≳ 2) a dimensionless discretization parameter, n ∈ N, and z ∈ [0, 1[ a plain 'z-shift' of the logarithmic discretization [49,50].
Consider therefore some arbitrary but fixed energy interval I l ≡ [ϵ l , ϵ l+1 ] of width ∆ϵ l ≡ ϵ l+1 −ϵ l > 0 and average energy εl ≡ 1 2 (ϵ l +ϵ l+1 ) within the continuum of the bandwidth.Here l > 0 will refer to energy intervals at positive energies, with energy increasing with increasing l (this is contrary to the NRG, hence l ∼ N − n with N the number of levels with ϵ l > 0).Since the helical mode in Fig. 1 is symmetric around ϵ = 0, the coarse graining for positive and negative energies is also chosen symmetrically around ϵ = 0, having ϵ −l = −ϵ l such that l < 0 corresponds to negative energies.The index l = 0 is generally considered excluded here, as it is typically used to refer to the entire bandwidth, e.g., as with r 0 in Eq. (A10b).Having ϵ −l = −ϵ l , the index l thus resembles momentum, in that the simultaneous inversion of momentum and energy for a given spin flavor leaves the Hamiltonian of the edge mode invariant.
The energy εl is differentiated here from ε l (note the different font) with the latter eventually used for the effective level position for the full interval l, typically having ε l ≃ εl similar but not exactly the same [50].When coarse graining the bath, the integral for the hybridization is split up into intervals, with {ĉ εσ , ĉ † εσ }=δ σσ ′ δ(ε−ε ′ ), such that [ĉ εσ ] = energy −1/2 .The coarse-grained discrete and thus dimensionless bath modes clησ are normalized and orthogonal with respect to spin σ, but not yet with respect to the impurity location η, as emphasized by using tildes.The states clησ for each individual interval need to be orthonormalized in any case even if the distance x is chosen on the grid in Eq. (A11).This orthonormalization has to occur prior to the subsequent mapping of the so-called star geometry in Eq. (A15) between the impurity and the bath states to an effective one-dimensional (1D) chain geometry, as the latter requires properly orthonormalized Fermionic levels.
Orthonormalization of the pair of bath states within each interval l for given spin can be achieved starting from their overlap which is again simply related to the fermionic anticommutator similar to Eq. (A10), Equation ( A32) is now the starting point for the Lanczos block-tridiagonalization of the bath, given the initial pair of orthonormal states f0 for the zeroth site with real coefficients in the symmetric / antisymmetric basis ĉl of the bath.The bath itself is represented in diagonal form just by the energies of ĉl , and hence is clearly also real.The resulting Wilson chain then consists of two intercoupled chains, with α n , β n ∈ R 2×2 and β 0 reserved for the coupling to f−1 ≡ U 0 d, i.e., the impurities.By choosing the symmetric / antisymmetric basis above, incidentally, it turns out on numerical grounds, that within numerical double-precision accuracy, there are no Creutz-couplings [51], i.e., in between fnη and fn+1,−η .That is, β n turns out diagonal in η, and the Wilson chain becomes a pure ladder (to be referred to as Wilson ladder [52]).Similarly, with the bath setup being symmetric in energy and thus at half-filling, there are also no onsite energies along the Wilson sites.All α n are thus purely offdiagonal, and encode the rung couplings of the Wilson ladder.In summary, the block entries α n and β n in Eq. (A33) have the structure, As argued with Eq. (A18), the coupling to the impurities is dominated at low energies by the symmetric channel.Therefore the coupling strengths along the symmetric and antisymmetric channel can differ, as quantified by the relative difference δ n .This splitting scales like δ n ∝ 1/x ≪ 1, and therefore is small for early Wilson shells close to the impurity for large impurity distance x ≫ 1.The case n = 0 is special, as it refers to the coupling of the impurities to the f0 states.In the spirit of symmetrically coupled Anderson impurities where there is actual hybridization of the impurities with the helical bath, it holds that δ 0 = 0 if x = x n is chosen on the discrete lattice Eq. (A11).In this case, the f0 states themselves are orthogonal already, and so there are strictly no impurity cross-coupling due to non-orthogonality of the f0 states.
The fact that the block-tridiagonalization can proceed in real arithmetic has several advantages, in practice.The switch to the symmetric / antisymmetric basis prevents certain numerical errors from piling up during the blocktridiagonalization that are related to slowly drifting complex phases.Also by dealing with real arithmetic, phase conventions on basis states only refer to signs.The pair of states within every tridiagonalization step represent symmetric / antisymmetric states, and hence are already orthogonal, by construction, so there is no immediate explicit need for a Schmid decomposition within the pair of states fn .
The total combined hybridization strength of the two impurities for given spin(-up) is [cf.Eq. (A15)] where the prefactor of 2 derives from the two impurity channels.This total combined hybridization strength is the same, on average, for the impurities f−1 in the symmetric / antisymmetric basis since l tr(T l T † l ) In the symmetric basis, however, the contribution of the individual channels can differ from each other, as reflected also in δ 0 ̸ = 0 symmetry, which is absent in a helical system.Rather, it is reduced to a discrete Z 2 particle-hole symmetry [cf.Sec.A 5].This is specific to the effective model of the helical edge used here, though, stressing that such a particle-hole symmetry is absent in 2D time-reversal invariant topological insulators [53].The complex phases with a 0 and b 1 close to the impurity are important and hence cannot be gauged away, since e.g.within a fixed spin flavor, the helical Hamiltonian is not time-reversal invariant [if derived from a real-space lattice the Hamiltonian necessarily would have to include complex spinorbit coupling terms; in the diagonal eigenbasis, as in Eq. ( 1), any quadratic Hamiltonian becomes real of course].
The coupling a 0 is fully determined by ⟨0| f0+ Ĥbath f † 0− |0⟩ with |0⟩ the vacuum state, and hence can be expressed analytically, The precise value for a 0 is thus sensitive on the discretization, as seen with the second line.In the continuum limit, ε l → ε may be pulled inside the integral, which yields the last expression having used εe iτ ε = (−i d dτ )e iτ ε .While r 0 = 0 on the grid (A11), the derivative is generally non-zero with alternating signs and decaying like 1/τ ∼ 1/x.Therefore a 0 is non-zero for any x = x n .That is, the two impurities weakly see each other right away upon a single application of Ĥbath .On intuitive grounds, one can compare this to a tightbinding chain with long-range hoppings that decay like 1/x [cf.Eq. (A12)] which also gives rise to a similar behavior.

d. Decoupling of anti-symmetric sector at low energies
The antisymmetric sector gets suppressed once energies are much below ε ≪ E x ≡ v x [Eq.(A9)], or more precisely in the limit |τ ∆ϵ l ≪ 1|, This occurs in practice at very low energies in the NRG context.With r l → 1 − , while s l+ approaches the finite value of 2 − , the smaller eigenvalue s l− → 0 + approaches zero as 1 − 1 − which is numerically ill-conditioned.Hence form a numerical perspective, it is computed via an expansion around small ξ ≡ τ ∆ϵ l /2 = x∆k l /2, i.e., from Eq. (A17b), When setting s l− strictly to zero below some threshold s l− < 10 −16 , the block tridiagionalization eventually switches over to a hopping amplitude that decays twice as fast, because the antisymmetric levels that actually couple have been exhausted.On the other hand, keeping the asymptotic dependence in Eq. (A49) down to the lowest energies considered, the hopping amplitudes along the Wilson ladder always decay like ω n ∼ Λ −n/2 .That is, the antisymmetric channel x [Eq.(A9)] one may think of the bath states coupled to their respective impurity.Therefore the physics represents two independent copies of one-impurity problems down to energies Ex.Around the energy scale of Ex the impurities start to coherently interact with each other (if the impurities are, for example, already Kondo screened above the energy Ex, then the two impurities remain decoupled down to zero energy).The 2-impurity physics takes place at energies ε ≲ Ex where the relevant description of the bath effectively changes from the local (R, L) to the symmetrized (+, −) representation.The antisymmetric channel (η = −) starts to smoothly decouple, but stays in the system as a passive spectating bath space.The symmetric channel (η = +) is the one that remains fully coupled to the impurities.While the representation of the bath in the main text is always in the symmetric/antisymmetric configuration [cf.Eq. (A42), except for a final rotation of f0 back to the local representation].If the same unitary mixing of modes were to be applied for subsequent Wilson shells still, the property of having two independent copies of bath modes remains intact down to Ex. remains in the system, throughout.The reason for this is that the decoupling occurs smoothly.So once the energy scale (or more precisely, the energy resolution) drops below v/x, the antisymmetric channel does not decouple in an instant, and so it stays in the system, as schematically depicted in Fig. 9.At energies much below v/x, however, one can show in practice that the Wilson chain, indeed, switches over to two fully decoupled chains [cf.Fig. 10(b)].While the symmetric sector which remains coupled to the impurities, shows a smooth decay of the hopping amplitudes, in the antisymmetric channel the hopping amplitudes along its corresponding leg in the Wilson ladder becomes increasingly alternating (lower legs in Fig. 10): namely the paired up antisymmetric levels at energies ±ε l .They form strong bonds along the Wilson chain at zero energy, where bonding and antibonding states reveal the original ±ε l states in the star geometry.
Since the impurity distance directly enters the coarse graining in the NRG, there will always be a qualitative change in the NRG energy flow diagram around the energy scale v/x [cf.Figs. 9 and 10].But this change in the representation of the bath can become irrelevant for static or dynamic properties from the perspective of the impurity.In this sense, the 'apparent' energy scale strongly visible in the standard NRG energy flow diagrams at the energy scale v/x may be irrelevant for the impurity.Nevertheless, this may leave minor artificial fea- .They are all real-valued, where black (red) color shows positive (negative) hopping amplitudes, respectively.Around the energy 1/x (shell n ∼ 40), the structure of the Wilson ladder changes, as qualitatively already argued in Fig. 9.The upper leg corresponds to the symmetric (η = +), and the lower leg to the antisymmetric channel (η = −).The rescaled Creutz couplings have amplitudes below double precision accuracy, hence are absent.(b) Same as in (a), except that starting from the position 'blocktrafo' towards the right a nearest-rung shell-mixing numerically determined block transformation was performed on top of (a).This shows that at low energies the system can be exponentially decoupled towards later shells into two independent channels.tures (wiggles) in the temperature or frequency dependence of physical properties around the energy scale of v/x.This effect is expected to be more pronounced for coarser discretization (larger Λ), but to diminish for smaller Λ.
e. Block-tridiagonal structure for opposite spin Switching the spin has the same effect as changing the sign of energy or spatial inversion [cf.discussion after Eq. (A18)].This is also explicitly reflected in the the variable τ ≡ σx v [Eq.(A4a)] that appears in much of the above treatment.Therefore by construction of the starting vectors for the blocktridiagonalization of the spin-up channel above, from the point of view of the impurity, spin-down couples to an identical Wilson ladder of its own, except that the local f0 modes couple to the impurity levels in reverse order.In this sense, the backtransformation to the local representation of the impurity in Eqs.(A47) is useful as a prior step.Then the Pauli matrix σ x below accounts for the reversed order, while everything else remains the same for spin down as for spin up.The effect of the cross terms on the impurity spectral function, for example, are therefore, This makes intuitive sense, since motion in between the two impurities occurs in opposite directions for different spins.

Dynamical correlation functions
Within the NRG approach, correlation functions are computed by evaluating the Lehmann representation [33,49] which can be carried out exactly at arbitrary temperature T within the full density matrix approach (fdm-NRG [34]).To be specific, a retarded dynamical correlation function for two local impurity operators B and Ĉ, with ϑ(t) the Heaviside step function, B(t) = e i Ĥt Be −i Ĥt , and ⟨. ..⟩T the thermal average, is computed via its spectral Lehmann representation with a and b complete sets of many-body eigenstates, having Ĥ|a⟩ = E a |a⟩ with E ab ≡ E b −E a , thermal weights ρ a = 1 Z e −βEa , β ≡ 1 T , and Z ≡ a e −βEa the partition function.The iterative diagonalization within the NRG approach explicitly generates the full many-body state space Ĥ|a⟩ ≃ E a |a⟩ above.In practice, this yields a tractable symmetry-respecting, eigenstate decomposition of the entire impurity Hamiltonian on the Wilson chain [54] which while approximate, is well-controlled and complete.
Standard fermionic and bosonic correlation functions have an (anti)commutator in their Green's function, with the commutator (s= − 1) for bosonic correlations such as spin-spin correlation functions, or the anticommutator (s= + 1) for fermionic correlation functions such as the fermionic local density of states.Equivalently, in frequency space G R BC (ω) = G BC (ω) + s G C † B † (−ω), where by construction at zero temperature A BC (ω < 0) = 0 since E ab ≥ 0 in Eq. (A53), such that the two contributions in the (anti)commutator to the full correlation function are separated in frequency space, since they exclusively contribute to positive or negative frequencies only.
The important point with Eq. (A53), as already also pointed out with the hybridization function in Sec.A 1, is that for offdiagonal correlations B ̸ = Ĉ, the spectral data cannot only be negative, but fully complex, i.e., A BC (ω) ∈ C, if the matrix elements themselves are complex.In this sense, one cannot simply write the spectral data A(ω) as − 1 π Im G(ω).Instead, by only taking − 1 π Im (. ..) of the propagator 1 ω + −E ab , the full Green's function in frequency space can be simply obtained by standard means, i.e., via Kramers-Kronig transformation, or by folding with the propagator, By construction, with Eq. (A53) one still also has full access to the well-known simple spectral sum rules, which may be complex for B ̸ = C † in the present helical setting.This is relevant, for example, when computing the oneparticle correlation function across the impurity ⟨ dσL ∥ d † σR ⟩ ω [e.g., see Eqs. (A6) or (A7) for the non-interacting case].

Symmetries
The global symmetries of the effective 1D model of the helical edge in Sec.II also manifest themselves in the structure of correlation functions.Aside from the symmetry U(1) charge ⊗ U(1) spin already discussed when introducing the helical model system in Sec.II and also explicitly exploited in our numerical simulations, the isolated helical edge can support further symmetries that are actually absent in the original full-fledged 2D topological system.The original topological aspect is reflected here in the fact that the isolated helical edge exists as a valid physical model in the first place.With this in mind, the isolated helical edge with two Kondo impurities (2HKM) located symmetrically around the origin also preserves • Z 2 time reversal symmetry (Z TRS ): momentum together with spin reversal is preserved by the helical edge.The local impurity-bath Kondo interaction is also spin-reversal symmetric, with the impurities themselves located symmetrically around the origin [cf.Eq. ( 5)].Therefore k → −k together with the reversal of the impurities, η → −η, leaves the local impurity-bath interaction invariant.Now the combined operation k → −k and η → −η is equivalent to spatial inversion.Hence, for our model setup with an isolated helical edge, TRS can be translated into spatial inversion with simultaneous global spin reversal.
• Z 2 particle-hole symmetry (Z p/h 2 ): the helical channel in Fig. 1 was chosen such that for every level at ε kσ > 0 there is a level at ε −kσ < 0, having ε kσ = −ε −kσ .By having half-filling, this converts a particle to a hole or vice versa.Furthermore, by construction, the Kondo interaction of the impurities with the helical channel are also particle-hole symmetric.
The Hamiltonian in Sec.II has two impurities located symmetrically around the origin [cf.Eq. ( 5)].From an Anderson impurity point of view with explicit hybridization as in Eq. (A2), this is the only point where complex numbers enter the total Hamiltonian.From a numerical perspective, the Fermionic operators ĉkσ and dησ can be encoded by real matrix elements, such that the only complex entry is the phase e ikx η 2 .This also holds when switching from Anderson-type hybridization to Kondo spin interactions when represented in terms of Ŝη ± or Ŝη z , as their matrix elements are also real.With this perspective, it holds that complex conjugation of the Hamiltonian, Ĥ → Ĥ * , is equivalent to reversing the locations of the impurities dησ → d−η,σ .Denoting the latter by R I , with R I H|a⟩ = H * R I |a⟩ it holds that if |a⟩ is an eigenstate of H to eigenvalue E a , then so is |a ′ ⟩ ≡ [R I |a⟩] * .Hence, with K denoting complex conjugation, R I ≡ R I K is an anti-unitary symmetry of the system, with |a ′ ⟩ = R I |a⟩ also an eigenstate of the Hamiltonian with the same eigenenergy.
As a consequence of the above symmetries, for example, the spectral data of spin-spin correlations as in Eq. (A56) is real, after all.Based on the Lehmann representation in Eq. (A53), one encounters the matrix elements For the case that η = η ′ , i.e., intra-impurity spin correlations, this product of matrix elements can be combined into |⟨a|..|b⟩| 2 , which is real.For the case of inter-impurity spin correlations, η = −η ′ , taking the complex conjugate and inserting Now since a and b can be chosen to also be simultaneous eigenstates of R I , i.e., having a = a ′ and b = b ′ with the same eigenvalue w.r.t.R I , the Lehmann sum of the spin-spin spectral functions, while having complex matrix elements, yields a purely real result.In practice, a and b are not eigenstates of R I .Yet by explicitly resorting to complete many-body basis sets within the NRG [54] evaluated within fmd-NRG [34], the spectral data has an imaginary contribution with relative strength comparable to numerical noise based on double precision accuracy, and hence can be ignored.

From Anderson to (anisotropic) Kondo type model
The entire discussion above assumed Anderson-type impurities that hybridize with the helical edge mode.Now if the local Coulomb interaction U with each impurity is large, charge fluctuations get frozen out.Therefore in the lowenergy regime, charge fluctuations can be integrated out locally with each impurity via Schrieffer Wolff transformation.Assuming half-filling of each impurity with a single magnetic spin-half moment, second order perturbation theory based on Eq. (A2) yields the Kondo-type interaction f † 0ησ f0ησ ′ the spin operator with respect to the bath site at the location of impurity η.Therefore a single impurity interacting with in a helical edge mode is identical to a regular 1-impurity Anderson or Kondo model without a helical character.The dynamically generated low-energy Kondo scale 47] for the 1-impurity problem with ρ 0 the one-particle density of states around the Fermi level, however, also represents a relevant energy scale for the 2-impurity problem.Assuming the impurity distance on the grid (A11), the local bath operators f0ησ are already properly orthonormalized, such that there are no issues with crosstalk between the impurities.
For spin-independent hybridization between the Anderson impurities and the bath, the resulting Kondo coupling in Eq. (A59) is SU(2) spin symmetric.While bearing in mind that the helical bath mode itself already has the SU(2) symmetry broken, the spin symmetry can also be broken at the level of the Kondo Hamiltonian, giving rise to anisotropic local spin-spin interactions.Assuming a single preferred direction (z) with J x = J y ≡ J ̸ = J z , J z > J (J z < J describes an easy-axis (easy-plane) regime, respectively.The global U(1) spin symmetry thus remains preserved.

Effects of finite bandwidth with Kondo interaction
When starting from the Anderson model, one may scale the local Coulomb interactions properly to infinity in relation to other parameters, with the result that also in the numerical setting, one effectively arrives at the Kondo model [55].For the Anderson model the effective bandwidth relevant for the impurities is cut off by U if U < D. However, by taking U ≫ D when transitioning towards the Kondo model, bandwidth keeps playing a considerable role, and the universal wide-band limit is approached rather slowly.
Here for the anisotropic 2HKM in the present case, if the Kondo scale T K (J, J z ) ≪ E R ≪ D = 1 for each impurity in-dividually is just several orders of magnitude smaller than all other energy scales, effects of finite-bandwidth are still considerable.In order to reach the wide-band scaling limit (here J ≲ 0.02), the single-impurity Kondo scale is already many many orders of magnitude lower than the band width D = 1 (T K ≲ 10 −40 ), as demonstrated in Fig. 11, and also consistent with the literature on the single-impurity case [55].In the present case, we see that similar arguments also carry over to the RKKY regime, even if Kondo physics is irrelevant (in the sense that it sets in at much smaller energy scales, even much smaller than any temperature of practical interest).
In the RKKY regime and in the wide-band limit, the impurities are expected to be well decoupled from the bath and described by In practice, however, one sees substantial deviations from this expectation value up to nearly 20% for J = 0.1 in the NRG data even using J z = 0 as shown in Fig. 11(a).Hence these deviations must derive from higher-order processes that go beyond second order PT [see App.D].The coupling to the bath remains finite in the low-energy regime, thus inducing fluctuations in the impurity spin configuration.These deviations can be reduced systematically by lowering the Kondo coupling J, e.g., having a deviation already below 1% for J ≲ 0.02.As shown in Fig. 11(a) the static expectation value ⟨S z L S z R ⟩ approaches −1/4 in a polynomial fashion as J is lowered, down to the smallest J considered.Therefore, indeed, the deviations seen in Fig. 11(a) are clearly due to finite bandwidth.This demonstrates that finite bandwidth does play an observable role in the Kondo setting [55] when comparing numerical to analytical results if the latter strictly assumed the wideband limit.If one considers the wide-band scaling limit reached within deviations in observables of about 10 −3 , this suggests approximately J ≲ 0.02 in Fig. 11(a).On the single impurity level, this already corresponds to astronomically small Kondo scales T K ≪ 10 −40 in Fig. 11(b) consistent with earlier NRG studies [55].From a physics point of view, however, we do not expect that the observed minor variations change the overall physical picture in any significant qualitative manner.It follows from the above that where ±a 2 with a ≥ 0 is some constant of integration.If the starting point has |z 0 | > |x 0 | (easy-axis), then the positive sign is chosen for a 2 , whereas the regime |z 0 | < |x 0 | (easy-plane) has the negative sign.The contours described by Eq. (B3) exactly reflect the RG paths of the anisotropic Kondo (parabolas or hyperbolas separated by |x 0 | = |z 0 |, as shown by the gray lines in Fig. 12).It simply also follows from Eq. (B2) that for x → 0, (x, z) stops flowing.The model is physically equivalent for x 0 → −x 0 (vertical flip in Fig. 12) as this can be absorbed into a gauge transformation of the spin basis.
For z 0 > 0 and |x 0 | > |z 0 |, the anisotropic Kondo model flows to the isotropic strong coupling regime.That is integrating out the bath starting from large D needs to be stopped at some D * > 0 where x and z diverge, which thus defines the dynamically generated low-energy Kondo scale T K := D * .
In the easy-axis Kondo regime (J z > |J x |), with z ≥ a, the RG differential equations diverge at The Kondo scale vanishes to exponentially small energy scales for J x → 0, i.e., a → z − 0 .In the easy-plane regime (|J x | > |J z |), with |x| ≥ a, the RG differential equations diverge at The single-impurity Kondo scale from the previous section may be compared to the RKKY energy E R [cf.Eq. ( 20b)].The latter decays with distance, hence is largest for short distances, with the shortest possible distance given by one 'lattice spacing' [cf.Eq. ( 2)], hence E R ≤ (ρ 0 J) 2 in the adopted units [Sec.II A].With this in mind, Fig. 12 also includes contours of constant T K (J, J z ) − (ρ 0 J) 2 ≤ 0 (green lines), with the thick line representing T K (J, J z ) = (ρ 0 J) 2 .Therefore the parameter regime to the right of the thick green contour always has T K > E R for any x ≥ 1, meaning that in this regime Kondo is always dominant over RKKY.This regime has been predicted in Ref. [1].In particular, this includes the isotropic Kondo for J ≳ 0.25, corresponding to a dimensionless coupling j 0 ≡ ρ 0 J = 0.125.While still clearly below 1., this is on the upper end of what may be considered acceptable on physical grounds for the Kondo model: Note that the poor-man's scaling approach for the Kondo model only includes second order virtual processes by integrating out the bath which assumes (and thus is justified only if) J 2 D ≪ D, i.e., j 0 ≪ 0.5.Conversely, note that coming from an Anderson model, one has j eff 0 ≃ 4Γ πU with hybridization strength Γ and onsite interaction U .So in order for the Kondo model to be justified in the first place, the onsite interaction U needs to be large enough so that charge fluctuations can be integrated out via virtual second order processes (cotunneling) in order to obtain a pure effective spin Hamiltonian.This clearly requires U ≫ Γ, and thus also j 0 ≪ 1 on physical grounds.
Conversely, to the left of the thick green line in Fig. 12, RKKY can dominate in the low-energy regime if a pair of impurities is just brought close enough to each other.In particular, RKKY also occurs within the entire white region to the left where T K = 0.With E R > T K = 0 then, the impurity distance x can be taken to infinity while still seeing RKKY in the low-energy regime (note that since RKKY is second order, with E R ∝ J 2 /x, the sign of the Kondo couplings is irrelevant).
Appendix D: Second order perturbation processes entering RKKY Complimentary to the field theoretic approach in the main text, the RKKY processes can also be analyzed from the point of view second order perturbation theory (PT).Ultimately, this generates the same effective Hamiltonian as in the field theoretic approach, but purely within the Hamiltonian setting, e.g.via the Feshbach formalism [58,59].
Let the full Hamiltonian of two Kondo-coupled impurity spins η ∈ {R, L} ≡ {+1, −1} at location ηx 2 , i.e., at mutual distance x, interacting with a one-dimensional helical bath encoded by the Fermionic operators ĉkσ be given by with a ∈ x, y, z, and bath spins (with τ a the Pauli matrices), assuming a total of N levels k for each spin.Furthermore, the energies of the bath ε kσ = σvk are constrained to within a finite half-bandwidth |ε kσ | ≤ D. Taking J≡J x =J y and σ± η ≡ ŝ± η ≡ ŝx η ± iŝ y η , the Kondo interaction can be written from the point of view of the helical Fermions as, with the matrix notation in σ∈{↑, ↓}≡{1, 2} for the operators of each impurity [see also Eq. ( 13)], with anisotropy ∆≡ Jz J and τ ± = 1 2 (τ x ± iτ y ).The isotropic case ∆=1 just reduces to T η = Ŝ † η • τ .The effective Hamiltonian in between the impurities, i.e., their direct interaction, is generated via second order PT in the Kondo setting (4 th order in the hybridization for an Anderson model).Then with the shortcuts 1 ≡ (k 1 σ 1 ), etc., a typical second order process is given by with degrees of freedom such as impurity location η∈{L, R}≡{−1, 1}, or spin σ or momenta k to be summed over eventually.P 0 is a projector into the target low-energy regime of the bath without acting at the impurities.This is in the spirit of generating the low-energy effective Hamiltonian such as the Feshbach-Fano partitioning [58,59].There, by construction, one needs to start out of the low-energy regime, say by considering some infrared cutoff energy D * ≪ D. In the low-energy regime described by P 0 then, all bath levels with energy ε kσ < −D * are strictly occupied and all energy levels ε kσ > +D * are strictly empty.Now in Eq. (D2), the perturbation T on the r.h.s.generates a particle-hole (p/h) pair with particle energy ε 3 > 0 and hole energy ε 4 < 0. There-fore as long as either ε 3 > +D * or ε 4 < −D * (or both), this represents state space outside the low-energy regime P 0 .If D * ≪ D, the overwhelming number of processes will involve both, particle and hole at energy cost above D * .There in order to return to the low energy regime, exactly the same p/h pair must be destroyed again.This constrains the second processes above to δ 14 δ 23 , with the shortcut notation δ ij ≡ δ kikj δ σiσj .Actually, the processes δ 14 δ 23 can be included all the way down to zero energy, i.e., there is apriori no need for an infrared energy cutoff D * , so one can assume D * → 0 + as long as the sum is well-defined (which it turns out it is).With the bath projected via P onto the exactly filled helical Fermi sea, this generates the effective, purely local inter-impurity Hamiltonian (D3) having taken the continuum limit i ≡ (k i σ i ) → (ε i σ i ), with ρ 0 the constant one-particle density of states of the helical edge for a given spin at the Fermi level as in Eq. ( 2).The phase factor in the second line rewritten in terms of energies is already also specific to the helical edge.As an aside, we note that the energy denominator of I ηη ′ σ1σ2 can be rewritten into an integral over imaginary frequencies, which resembles a Matsubara summation at zero temperature.Via contour integral, the result is non-zero only if ε 1 and ε 2 have opposite sign as encoded into the Heaviside step function on the r.h.s.This already reflects the particle-hole nature of the excitations in Eq. (D3) where ε 1 < 0 and ε 2 > 0.
The case with opposite signs is included in Eq. (D3) via the overall sum.Assuming η̸ =η ′ , then the case η ↔ η ′ gives the additional contribution above, while also exchanging the labels 1 ↔ 2, such that the sign in the phase factor is properly restored, while also having [ T η , T η ′ ]=0 for η̸ =η ′ .Hence by summing over all second order processes, one can connect the second order perturbative approach here to the doublepropagator structure in Eq. ( 15) in the main text obtained from the field-theoretic approach.Collecting phase factors in k 1 and k 2 , the latter thus permits the interpretation that any second order process in the effective impurity Hamiltonian requires the free propagator of two particles shuttling back and forth in between the impurities (one particle needs to propagate the distance +x, and another the distance −x).Here in the helical setting, however, the directions that particles can move are constrained depending on their spin.This manifests itself in the overall structure of the resulting RKKY Hamiltonian [1].The integral I ηη ′ σ1σ2 in Eq. (D3) is generally well defined when both, ε 1 and ε 2 approach zero energy, irrespective of the phase factor since The imaginary part form i0 + does not contribute, as by the integral limits, it can only contribute at ε = ε ′ = 0, where by the double integral the integrated weight vanishes.When η=η ′ , the complex phases drop out, and with S 2 η ∝ 1 1 the integral in Eq. (D3) just generates a plain irrelevant shift in the global energy reference as estimated above (the generation of the single impurity Kondo couplings needs to consider finite D * and relax the condition δ 14 δ 23 above).Hence the following discussion focusses on η̸ =η ′ , i.e., η ′ = − η, in which case the phase factors in the enumerator are non-trivial and reflect the underlying helical physics.In contrast to the energy denominator, however, the phase factors carry momenta as arguments.Converting these into energies, in the helical setting with ε kσ = σvk, this involves signs depending on the spin as already indicated with the integral I ηη ′ σ1σ2 in Eq. (D3).Contributions that involve a spin flip (σ 2 = − σ 1 ), will have opposite relative sign of ε 1 vs. ε 2 in the phase factor in Eq. (D3) as compared to the energy denominator, whereas in the absence of a spin flip the relative sign is the same.To evaluate this integral including the phases, it is therefore convenient to change variables to symmetric and antisymmetric combinations With ε 1 < 0 and ε 2 > 0, by construction, 0 ≤ ε ≤ 2D and ε ′ ∈ [−ε, ε].The integral will converge with large D, such that the upper integral limit can be taken more loosely by deforming the integration area to ε ≲ D = 3D 2 , with the upper integral limit D assumed large, and eventually taken to infinity if well-defined.With this one obtains for the case including a spin flip, with τ ≡ with the dimensionless Kondo coupling strength j 0 .Overall, this generated direct impurity interaction is ferromagnetic and non-oscillatory with a plain decay with inverse distance x, having the RKKY energy scale as already encountered in Eq. (20) in the main text.The RKKY Hamiltonian is independent of the bandwidth D, and aside from the dimensionless Kondo coupling strength j 0 , only references the coherence scale E x ≡ 1 τ with τ = x/v the time required for a helical particle to travel from one impurity to the other.With the lattice spacing in Eq. ( 2), nevertheless, this energy scale may be rewritten as E x = D/π x/a .That is when measuring distance on the grid (A11), by definition, this involves a finite bandwidth, such that the bandwidth does appear in the definition of E x .In the continuum wide-band limit, however, the natural way to think about the coherence scale is E x = 1/τ without any reference to bandwidth.
In the absence of a spin flip, i.e., σ 1 = σ 2 the integral can be similarly evaluated, 2 , assuming that the integral is well-defined in the wideband-limit.Again the contribution from the imaginary part i0 + in the denominator vanishes, since dε δ(ε) εe i τ ε → 0, as already expected from Eq. (D5).By summing over spin or location, with τ ∝ ησ 2 , only the real part remains, ησ1 I η=−η ′ σ1=σ2 ≃ 4 τ sin τ D which seems to suggest 1/x behavior similar in magnitude to Eq. (D8).In the present case, however, the integral remains sensitive on the bandwidth D. Therefore the assumption for introducting D above does not hold, and the integral needs As apparent from the oscillating averaging structure in the second term, it also vanishes for large λ.In the asymptotic form for large λ, the leading term of the cosine integral Ci(λ) ∼ sin λ λ again drops out on the grid (A11).Therefore the subleading term Ci(λ) ∼ cos λ λ 2 becomes the dominant one for large λ.But with DCi(λ) ∼ D λ 2 ∼ 1 Dx 2 , this does not just decay faster over distance as compared to RKKY for a normal 1D metallic mode, but is also suppressed in the wideband limit.Therefore as expected the ZZ contribution to the RKKY Hamiltonian properly vanishes for D → ∞ even for finite J z .At finite bandwidth, however, there is a finite return probability, resulting in a small but finite S z L S z R interaction strength across the impurities.This may be considered acceptable on physical grounds, bearing in mind that for a true 2D model the spin-dependent directedness of motion only concerns the edge but not the bulk states gapped out to energies > D.
Finally, we point out that above integrals also appear in the theory of a normal metallic edge when one assumes the same linear dispersion and finite bandwidth as for the helical system here.However, by having additional branches of the opposite helicity in the dispersion for back-propagation, the same integrals appear and are summed over all spin interactions XX, YY, and ZZ.This way the RKKY interaction becomes isotropic for a normal metallic edge.The above argument that the leading oscillatory term vanishes is particular to the one-particle dispersion chosen here, and would also apply to the normal metallic edge with the dispersion indicated.This differs crucially from a system of free particles in a Fermi sea with a quadratic dispersion and a finite Fermi energy.In this case, the analytic structure of the respective integrals is different.As a consequence, this permits the leading 2k f oscillatory term ∝ 1/x to be present in a normal metallic 1D mode [43].

FIG. 3 .
FIG.3.NRG analysis of the 2HKM for J = 0.1 (i.e., j0 = 0.05) -Each row corresponds to a different integer impurity distance x as indicated with the left panel, which increases exponentially from top to bottom.The vertical axis is the same for all panels, namely Jz/J ≡ jz/j0, whereas the horizontal axis represents energy in various forms.Here we adopt the NRG energy scale point of view, throughout, which starts at large energies at the left and then proceeds towards exponentially smaller energy scales towards the right.Rows (columns) are labelled by letters (numbers), respectively, e.g., having (a1) for the upper left panel.The vertical black dotted [solid blue] marker replicated in all panels indicate the coherence scale Ex (33) [RKKY scale ER (20)], respectively, as labelled in panel (a1), having Ex/ER = 1/πj 2 0 = 127.3,throughout.Similarly, the brown solid curved line shows the analytical single-impurity Kondo temperature TK(J, Jz) [App.B] for reference.This curve TK is visually cutoff at ER, because TK is irrelevant at lower temperatures.The brown marker in panel (h2) shows the finite intercept TK(J = 0.1, Jz) = 4.31•10 −8 at Jz = 0.Each column shows the quantity indicated above the top panel: Panels (:,1) [first column] show the static (equal-time) inter-impurity correlations ⟨S L z S R z ⟩ vs. temperature and Jz.Panels (:,2) [second column] gives a visual impression of the changes along the NRG energy flow diagram vs. energy scale ωn ∼ Λ −n/2(for precise prefactors, see[46]).This encodes cumulative effective thermal weights in three arbitrary but fixed low-energy symmetry sectors into red-green-blue (RGB) colors for all odd Wilson shells n (to avoid even/odd effects) based on an effective inverse temperature βn ≡ 4ωn.The symmetry sectors chosen for this were q ≡ (Q; S tot z ) ∈ (0; 0, 1, −1) with Q the total charge relative to half-filling.The remainder of the columns shows dynamical spin-spin correlation functions ⟨ Ŝη ∥ Ŝη ′ † ⟩ω as indicated at the top of each column at zero temperature (T = 1.8 10 −11 ).The solid-dotted lines shows the respective derived inverse static susceptibility T ηη ′

FIG. 4 .
FIG. 4. NRG spectral data of the 2HKM for J = 0.1, Jz = 0, x = 10 for the dynamical correlations functions as partly indicated with the left axis [transverse in the upper panels, longitudinal in the lower panels; intra-impurity for the left panels, whereas in the right panels also includes the inter-impurity correlation for comparison; note the sign in the legend with panel (d)] -Left panels: individual spectral data for z-shifted logarithmic discretization (nz = 8 curves with z ∈ [0, 1[, having Λ = 4, broadening σ = 0.1), and the corresponding z-averaged data in the right panel.The horizontal and vertical axis are globally scaled by the z-averaged low-energy scale TS z = 1.4•10 −4 = 0.556 ER only, so with spectral sum rules in mind, the data shown is of order one.

FIG. 8 .
FIG.8.Low-energy scales vs. distance from inverse static susceptibilities [cf.Eq. (34)] as indicated in the legend of panel (a) for Jz = 0, throughout.Left panels (J=0.05)[right panels (J=0.10)]have temperature T < 10 −10 [the Kondo temperature TK ∼ 4•10 −8 ] as low-energy cutoff, respectively.The blue vertical markers translate the low-energy scale to distance, i.e., x0 ≡ j 2 0 /max(T, TK) which also corresponds to the crossing point of ER (black dashed) with max(T, TK) (horizontal marker).The lower panels redraw the data in the upper panels, but vertically scaled by TS.The shading indicates the standard deviation of the color matched energy scale due to z shifts in the NRG discretization.

FIG. 9 .
FIG. 9. Schematic representation of the 2-impurity Wilson setup of the bath giving rise to two intercoupled chains (for a quantitative example, see Fig.10).The system moves to small energies towards the right with the Wilson shells n having energy εn ∼ Λ −n/2 .At high energies ε ≫ Ex ≡ vx [Eq.(A9)] one may think of the bath states coupled to their respective impurity.Therefore the physics represents two independent copies of one-impurity problems down to energies Ex.Around the energy scale of Ex the impurities start to coherently interact with each other (if the impurities are, for example, already Kondo screened above the energy Ex, then the two impurities remain decoupled down to zero energy).The 2-impurity physics takes place at energies ε ≲ Ex where the relevant description of the bath effectively changes from the local (R, L) to the symmetrized (+, −) representation.The antisymmetric channel (η = −) starts to smoothly decouple, but stays in the system as a passive spectating bath space.The symmetric channel (η = +) is the one that remains fully coupled to the impurities.While the representation of the bath in the main text is always in the symmetric/antisymmetric configuration [cf.Eq. (A42), except for a final rotation of f0 back to the local representation].If the same unitary mixing of modes were to be applied for subsequent Wilson shells still, the property of having two independent copies of bath modes remains intact down to Ex.

FIG. 10 .
FIG.10.Typical Wilson ladder for the 2-impurity helical system as described by Eq. (A33) for an impurity distance x = 10 6 .The labels on top indicate the Wilson shell n (based on Λ = 2).The widths of the bonds are proportional to the hopping amplitudes rescaled by Λ −n/2 .They are all real-valued, where black (red) color shows positive (negative) hopping amplitudes, respectively.Around the energy 1/x (shell n ∼ 40), the structure of the Wilson ladder changes, as qualitatively already argued in Fig.9.The upper leg corresponds to the symmetric (η = +), and the lower leg to the antisymmetric channel (η = −).The rescaled Creutz couplings have amplitudes below double precision accuracy, hence are absent.(b) Same as in (a), except that starting from the position 'blocktrafo' towards the right a nearest-rung shell-mixing numerically determined block transformation was performed on top of (a).This shows that at low energies the system can be exponentially decoupled towards later shells into two independent channels.

2
A58) which is diagonal in η.When projecting the into the lowenergy Kondo regime of the Anderson model by scaling up local interactions, this leaves the representation of the bath untouched.Therefore from the NRG point of view the bath remains completely unaffected by whether one resorts to an Anderson-type or low-energy Kondo-type impurity setup.Based on the coarse-grained version in Eq. (A32) then, ĤK ∼ 2J η∈L,R Ŝη • Ŝ0η , (A59) with J the Kondo coupling, Ŝη ≡ Ŝdη the spin operator of impurity η, and Ŝ0η ≡ σσ ′ τ σσ ′ 2

FIG. 11 . 2 x
FIG. 11.Static ⟨S z L S z R ⟩ impurity correlations for TK ≪ ER ∼ J 2 x obtained by NRG towards the wide-band limit D ≫ J(= Jx = Jy) for Jz = 0. (a) Deviations of the computed ⟨S z L S z R ⟩ from the expected RKKY value of −1/4 on a loglog plot.Data is shown for various impurity distances x, with the polynomial fit obtained for x = 10 showing approximate quadratic behavior (red dashed line).(b) Comparison of energy scales for the same parameter range as in (a) where the temperature T ∼ 10 −9 used in NRG simulations in (a) is indicated by the horizontal gray line.The data for the Kondo scale TK was obtained via poor-man's scaling [cf.App.B].Red line shows a simple estimate for the RKKY energy at x = 10, demonstrating that (a) is deeply in the RKKY regime, i.e., ER ≫ TK, T .

FIG. 12 .
FIG. 12. Kondo scale of the anisotropic Kondo model from poor Man's scaling [Eqs.(B4) and (B5)].Jx < 0 is physically equivalent to Jx > 0, such that the lower half of the panel is a mirror image of the upper half.The thick black horizontal line (Jz < 0 for Jx = 0) represents a stable line of fixed points without a Kondo scale.The entire white region below the diagonal lines to the left flows towards it.The black dotted horizontal line (Jz > 0 for Jx = 0) represents a line of unstable fixed points that flow to strong coupling for any small but finite Jx (strictly at Jx = 0 one has TK = 0, represented by white color).The gray lines represent RG flow contours as in Eq. (B3).The Kondo temperature is defined by the starting point (x0, z0) on such a contour.The Kondo scale increases monotonically from left to right, and also with increasing |Jx|.The green contours show lines of constant TK (J, Jz) − (ρ0J) 2 ≤ 0. The thick green contour to the right describes TK (J, Jz) = (ρ0J) 2 = ER(x = 1) which represents the largest RKKY energy possible by putting two impurities right next to each other, with the closest distance being one 'lattice spacing' [cf.Eq. (20b)].

with
to be evaluated more carefully.The exact representation of the integral in the first line of Eq. (D10) yields with z ≡ τ ε, λ ≡ τ D, and λ ≡ | λ| = τ D, Ci() the cosine integral function.The first term vanishes on the grid (A11) having |λ| = τ D = πx/a a multiple of π.