Kane-Mele Hubbard model on a zigzag ribbon: Stability of the topological edge states and
quantum phase transitions
Chung-Hou Chung,1,2Der-Hau Lee,1and Sung-Po Chao2,3
1Department of Electrophysics, National Chiao-Tung University, HsinChu, Taiwan, 300, Republic of China 2Physics Division, National Center for Theoretical Sciences, HsinChu, Taiwan, 300 Republic of China
3Department of Physics, National Tsing-Hua University, HsinChu, Taiwan, 300 Republic of China (Received 10 February 2014; revised manuscript received 13 May 2014; published 14 July 2014) We study the quantum phases and phase transitions of the Kane-Mele Hubbard (KMH) model on a zigzag ribbon of honeycomb lattice at a finite size via the weak-coupling renormalization group (RG) approach. In the noninteracting limit, the Kane-Mele (KM) model is known to support topological edge states where electrons show helical property with orientations of the spin and momentum being locked. The effective interedge hopping terms are generated due to finite-size effect. In the presence of an on-site Coulomb (Hubbard) interaction and the interedge hoppings, special focus is put on the stability of the topological edge states (TI phase) in the KMH model against (i) the charge and spin gaped (II) phase, (ii) the charge gaped but spin gapless (IC) phase, and (iii) the spin gaped but charge gapless (CI) phase depending on the number (even/odd) of the zigzag ribbons, doping level (electron filling factor) and the ratio of the Coulomb interaction to the interedge tunneling. We discuss different phase diagrams for even and odd numbers of zigzag ribbons. We find the TI-CI, II-IC, and II-CI quantum phase transitions are of the Kosterlitz-Thouless (KT) type. By computing various correlation functions, we further analyze the nature and leading instabilities of these phases. The relevance of our results for graphene is discussed.
DOI:10.1103/PhysRevB.90.035116 PACS number(s): 03.65.Vf, 71.27.+a, 73.20.−r I. INTRODUCTION
Recently, there has been growing interest in topological in-sulators (TIs) and superconductors which support gapless edge (surface) states while the bulk remains insulating [1,2]. These surface states come as a consequence of the spin-orbit (SO) couplings, and are protected by the time-reversal symmetry (TRS) [1,2]. The topological nature of TIs lies in the nontrivial topological Z2 invariant [3], while it becomes trivial for an ordinary band insulator (BI). The theoretical predictions [4–6] of TIs have been soon observed experimentally in various insulators with strong SO couplings [7]. In two-dimensional systems, these topological states have been predicted in the framework of the quantum spin Hall insulator (QSHI) [3,8–12], and have been realized experimentally soon after in HgTe/CdTe quantum well structures [4]. Unlike the integer quantum Hall insulator where the chiral (one propagating mode of electrons with a single spin species) edge states are generated by an external magnetic field which breaks TRS, the TRS preserving QSHI systems lead to helical edge states in the absence of a magnetic field in which propagation direction at one edge is opposite for opposite spins [12]. These one-dimensional helical edge state electrons are protected by TRS [3] and are free of spin-flip backscatterings [2]. As a result, they lead to perfect transmission in charge transport along the edge [13].
A simple theoretical model was first introduced by Hal-dane [8] and later proposed by Kane and Mele [3,9] (the KM model) to capture the helical edge states of QSHIs. The KM model was aimed to describe edge states in graphene. Though the SO coupling in graphene is expected to be too small to observe the edge states, the KM model is regarded as a generic model for 2D TIs. The existence of the helical edge states in KM model has been well studied. Recently, more attention has been put on the stability, exotic quantum phases and
phase transitions of the helical edge states and possible exotic quantum phases in the correlated Kane-Mele Hubbard [14–21] model upon including the on-site Coulomb repulsions (the Hubbard U > 0 term) in the KM model. In a pioneering work by Meng et al. in Refs. [14,15] via quantum Monte Carlo (QMC) and dynamical mean-field approaches, the helical edge states are stable up to a finite Hubbard interaction, and a gaped spin-liquid phase was predicted in the phase diagram of the KM Hubbard model at half filling for small to intermediate range of U . The authors in Ref. [19] have studied the effects of long-range Coulomb interactions on the edge states of a finite-sized zigzag KM ribbon. The 1D Luttinger liquid physics with power-law correlations for the helical edge states in the presence of a finite on-site electron-electron interaction has been addressed in the framework of the KM Hubbard model analytically via bosonization in Refs. [22,23] and numerically via QMC in Ref. [18]. Meanwhile, the doping effect on the KM Hubbard model was addressed in Ref. [24] where the spin liquid phase was argued to become a superconducting state.
In this paper, we present a theoretical analysis on the KM Hubbard model at half-filling and away from half-filling from a different perspective: we analyze the model on a finite-sized zigzag ribbon (where the helical edge states have been realized numerically in the tight-binding KM model [17]) with a ribbon width L= (N − 1)b (N being the number of zigzag chain in a ribbon and b is defined in Fig.1) in the weak-coupling (weak on-site Coulomb U ) limit via perturbative renormalization group (RG) combined with the bosonization approaches. Note that at a general level, the on-site Hubbard interaction U term considered here can take either positive (repulsive) or negative (attractive) values although the repulsive Hubbard interactions are more likely to be realized in QSHIs. Note also that one can alternatively study the model on an armchair ribbon, which was suggested to support edge states in graphene (equivalent to the KM model without SO coupling) [25]. We
FIG. 1. (Color online) Honeycomb lattice of a finite-sized zigzag ribbon of the tight-binding Kane-Mele model with the ribbon size N= 4 (N being the number of zigzag chains along x axis) along y axis. The honeycomb lattice consists of two interpenetrating triangular lattices denoted by sublattice A (dark circles) and sublattice B (open circles) with lattice vectors a1 and a2 (dashed arrows). The zigzag ribbon shows translational symmetry along x axis. The nearest-neighbor lattice vectors between nearest-neighbor A and B sites are denoted by ei=1,2,3 with a lattice constant a. The red (black) arrows within sublattice A(B) represent the directions of the next-nearest-neighbor hopping term λSOin the KM model (see text). The gray shaded region represents for the super unit cell of the zigzag ribbon, which repeats itself along x axis.
shall emphasize here the stability of the helical edge states against the combined short-ranged on-site Coulomb (Hubbard) interaction and finite-size effect of the zigzag KM ribbon, as well as possible other emerged quantum phases and phase transitions (QPTs) [26] among them.
The finite-size effect manifests itself in the structure of the energy spectrum and in an effective interedge tunneling term. We further find that these behaviors for even number of zigzag KM ribbons (N= even) are different from those for N = odd. For N = even, a finite energy gap is found at half-filling where the Fermi energy is at the Dirac point ka= π. This small gap is due to breaking of the sublattice translational invariance at the boundaries, and can be explained in terms of an effective finite single-particle interedge tunneling, which decays exponentially with increasing L. Away form half-filling, the energy dispersion becomes gapless at the Fermi level. For N= odd, however, the energy spectrum is gapless and the single-particle interedge tunneling vanishes for both half-filling and away from half-filling. Nevertheless, for both
N = even and N = odd, two-particle processes, effective
interedge two-particle spin-flip and interedge umklapp (two-particle backscattering) terms, are generated via second-order interedge hoppings.
Our stability analysis of the KMH ribbon is summarized as follows. For N = even, the energy gap at half-filled at the Dirac point gives rise to a charge and spin gaped (insulating) (II) phase [22]; at a generic filling, however, the two-particle processes when combined with the effect of the Hubbard
U term lead to the instabilities of the helical edge states towards a charge gapless but spin gaped (CI) phase [22] in the RG analysis via the Kosterlitz-Thouless type of quantum phase transitions. When L→ ∞, the interedge hopping term vanishes, the TI phase at half-filling is unstable against the charge gaped but spin gapless (IC) phase [22] for arbitrary
U >0, while it is stable away from half-filling. For N= odd, the single-particle interedge tunneling is absent, while the combined two-particle interedge hoppings and the on-site Coulomb interactions make the TI unstable for any finite U or interedge tunneling. As a result, the TI phase moves towards CI or IC or II phase depending on the ratio of Coulomb interaction and the interedge tunneling. The phase transitions for II-IC and II-CI are of the KT type.
By computing various correlation functions, we further analyze the instabilities of the helical edge states, the CI and IC phases towards the charge-density-wave (CDW), spin-density-wave (SDW) as well as the singlet (SS) and triplet (TT) superconducting states.
The remaining parts of the paper is organized as follows. In Sec.II, the Kane-Mele Hubbard at a finite size is introduced. The model is re-expressed in terms of the scalar and vector current operators. In Sec.III, the stability of the helical edge states is addressed via weak-coupling RG analysis. We also address the nature of the quantum phase transitions between the TI and other quantum phases. We conclude in Sec.IV.
II. MODEL HAMILTONIAN
A. The noninteracting Kane-Mele zigzag ribbon Before studying the interacting Kane-Mele Hubbard model, we summarize the main results for the noninteracting Kane-Mele (KM) model on a zigzag ribbon of honeycomb lattice given by the following Hamiltonian [3]:
HKM= −t ij,σ c†iσcj σ + iλSO ij,σ νijciσ† s zc j σ + H.c., (1) wherei,j and i,j refer to the nearest-neighbor (NN) and next-nearest-neighbor (NNN) sites, respectively. The NN and NNN lattice vectors for the honeycomb lattice are denoted respectively by ei=1,2,3and ai=1,2[17]: e1 = ¯a(0,1), e2 = ¯a/2( √ 3,−1), e3= ¯a/2(− √ 3,−1), a1 = ¯a/2( √ 3,3), a2 = ¯a/2(− √ 3,3) (2)
with ¯a being the lattice constant between nearest-neighbor
Aand B. The spin-orbit coupling term is represented by the imaginary NNN hopping λSOterm within the same sublattice where νij = 1 for i,j ∈ A (red counterclockwise arrows in
Fig.1) and νij = −1 for i,j ∈ B (blue clockwise arrows in
Fig.1). In the absence of the SO coupling, the KM model on zigzag ribbon reduces to the tight-binding Hamiltonian of a 2D zigzag graphene nanoribbon (ZGNR) [27], which shows two inequivalent Dirac points located at k≡ kx= ±2π3a with
kxbeing momentum along x axis with a≡
√
3 ¯a. Meanwhile, there exists a zero-energy flat band extended in the interval of 2π/3 ka 4π/3, known to correspond to the edge state of ZGNR [27,28]. It has been shown that the magnitudes of the edge state wave functions decay exponentially with distance away from the two edges, and the edge states are completely localized at the edges for ka= π [29,30].
In the presence of SO coupling, the KM Hamiltonian
HKM for a finite-sized zigzag ribbon (see Fig. 1) on honey-comb lattice supports helical edge states R,↑(↓)1(2),L,↓(↑)1(2) with
topological nature [3,17]. Here, R,↑(↓)1(2) stands for the wave function of the right-moving edge state electron with spin up (spin down) along the edge 1 (2), respectively. The indices 1 and 2 refer to the top and bottom edge, respectively. Similarly,
L,↓(↑)1(2) stands for the wave function of the left-moving edge state electron with spin down (spin up) along the edge 1 (2), respectively. The helical nature of these topological edge states manifest itself in the lock-in between the electron spin configuration and the direction of its momentum.
In the limit of large ribbon size N 1, the electron operator
cσi(x) near the edge is decomposed approximately in terms of these well-localized edge states as
c1(2)↑(↓)(x)≈ R,↑(↓)1(2)(x)e
ikFx, c↓(↑)
1(2)(x)≈ L,↓(↑)1(2)(x)e−ikF
x.
(3) The Hamiltonian of the edge Hedgeis therefore given by
Hedge= −ivF
dx(R,†↑1∂xR,↑1− L,†↓1∂xL,↓1
+ R,†↓2∂xR,↓2− L,†↑2∂xL,↑2) with vFbeing the Fermi velocity.
At a finite system size, however, the edge state electron wave functions acquire an additional functional dependence on y axis [c↑(↓)1(2)(x,y)] and are found to extend over a finite range in bulk via diagonalizing the tight-binding KM ribbon. The Hamiltonian of the edge states in this case is given by
Hedge= vF
dk
dyk[ ¯R,†↑1(k,y) ¯R,↑1(k,y) − ¯L,†↓1(k,y) ¯L,↓1(k,y) + ¯R,†↓2(k,y) ¯ ↓ R,2(k,y)− ¯ †↑ L,2(k,y) ¯ ↑ L,2(k,y)], (4) where ¯R/L,↑(↓)1(2)(k,y) are the edge state electron operators for a KM ribbon at a given momentum k and y obtained via partially Fourier transforming c↑(↓)1(2)(x,y) along the x axis:
¯
R,↑(↓)1(2)(k,y)=
dxe−ikxc↑(↓)1(2)(x,y),
(5) ¯
L,↓(↑)1(2)(k,y)=
dxe−ikxc↓(↑)1(2)(x,y).
Note that ¯R/L,↑(↓)1(2)(k,y) can be obtained numerically as the eigenstates of the Dirac dispersed helical edge states via diagonalizing the finite-sized zigzag KM ribbon. As shown in Figs.2and3, we numerically diagonalize the KM model at N = even (N = 4,16) and N = odd (N = 5,15) zigzag ribbon [17,31]. Two pairs of Dirac dispersed edge states ( ¯R,↑(↓)1(2), ¯L,↓(↑)1(2)) emerge in the energy spectrum of a finite-sized KM zigzag ribbon, and they tend to intersect at the Dirac points ka= ±π. However, at the Dirac points, a finite energy gap is developed for N = even, while no gap is seen for all
N = odd (see Fig.2). We shall focus on this even-odd effect in more details below.
Similar to the case for ZGNR, for 2π/3 ka 4π/3, we find the square magnitude of the two degenerate edge state eigenfunctions|(y)|2= | ¯
L/R,i(k,y)|2(except for N =
even and ka= ±π) show a symmetrical exponential decay
0 π 2π ka -2 -1 0 1 2 E/t 0 π 2π ka -2 -1 0 1 2 E/t (a) (b)
FIG. 2. (Color online) Energy spectrum of the finite-sized Kane-Mele model on a zigzag ribbon for (a) N= 4 and (b) 5 of honeycomb lattice. Here, we set t= 1 and λSO/t= 0.2.
from one edge to the other with respect to the ribbon center (y= L/2) from both edges into the bulk as a function of the distance to the corresponding edge. Here, y measures the distance to the edge along y axis and y= 0 corresponds to the first (top) zigzag chain. Also, to simplify the discussions, we use an integer index y/b+ 1 = Ni= 1,2, . . . ,N with
y= (Ni− 1)b for labeling the Nith zigzag chain along y axis
for a ribbon with N zigzag chains; y= 2b corresponds to the position of the third (Ni= 3) zigzag chain. As shown in
Figs.4(b)and5, the decay of these edge states is well fitted by the following exponential form:
| ¯L/R,i(k,y)|2 ∝ e−βy/b, (6)
where β is the decay constant depends on the momentum k. For N= even and at the Dirac point ka = π, we find the right and left moving edge states get hybridized so that the square magnitudes|(y)|2= | ¯
hy,i(y)|2of the two degenerate edge
states are maximized on both edges [see Fig.4(a)]. Note that we find, via eigenvector analysis of our numerical results through
0 π 2π ka -2 -1 0 1 2 E/t 0 π 2π ka -2 -1 0 1 2 E/t (a) (b)
FIG. 3. (Color online) Energy spectrum of the finite-sized Kane-Mele model on a zigzag ribbon for (a) N= 16 and (b) 15 of honeycomb lattice. Here, we set t= 1 and λSO/t= 0.2.
0 2 4 6 8 10121416 (y / b) + 1 10-18 10-15 10-12 10-9 10-6 10-3 100 |Ψ|2 0 2 4 6 8 10121416 (y / b) + 1 10-10 10-8 10-6 10-4 10-2 100 |Ψ|2 (a) (b)
FIG. 4. (Color online) The square magnitude of the edge state wave function ||2 of the KM zigzag ribbon at half-filling as a function of y/b+ 1 (defined in text) for N = 14 and (a) for ka = π and (b) for ka= π ± 0.2π. Here, ||2(blue circles and red squares) represent for the square magnitude of the two edge state wave functions, which are degenerate eigenstates at the corresponding wave vector k. In (a), the two hybridized degenerate edge state wave functions = hyb,i=1,2 (red and blue symbols) lead to the same square magnitude,|hyb,1|2= |hyb,2|2, in (b), we make the following identifications: (y)= R,↑1(blue) and (y)= L,↑2(red). The solid lines are guides to the eyes in (a), and in (b) they are fits to the exponential form in Eq. (6). We set λSO/t= 0.1.
exact diagonalization of the finite-sized KM ribbon, that these distinct two hybridized edge state wave functions: ¯hy,1(y)=
¯
hy,2(y) show the same magnitudes,| ¯hy,1(y)| = | ¯hy,2(y)|. Numerically, the values of|(y)|2as a function of y for a given edge state are obtained approximately by summing over the
0 2 4 6 8 10121416 (y / b) + 1 10-20 10-16 10-12 10-8 10-4 100 |Ψ|2 0 2 4 6 8 10121416 (y / b) + 1 10-10 10-8 10-6 10-4 10-2 100 |Ψ|2 (a) (b)
FIG. 5. (Color online) The square magnitude of the edge state wave function ||2 of the KM zigzag ribbon at half-filling as a function of y/b+ 1 for (a) N = 15 and ka = π and (b) for N= 15 and ka = π ± 0.2π. Here, ||2(blue circles and red squares) represents for the square magnitude of the two edge state wave functions, which are degenerate eigenstates at the corresponding wave vector k. The solid lines are fits to the exponential form in Eq. (6). We set λSO/t= 0.1. Note that in (a) ||2is shown for only even values of y/b (see text).
square of the matrix elements of the corresponding edge-state eigenvector contributed from both sublattices:
|(y)|2= |
A(y)|2+ |B(y+ ¯a/2)|2.
We also find that the square magnitude|(y)|2at ka= π for N= even [see Fig. 4(a)] oscillate along y axis. Similar oscillations are found for N= odd but not shown in Fig.5(a)as the values of|(y)|2for N = odd near edges are vanishingly small and go beyond the logarithmic scale shown there. This oscillatory behavior agrees qualitatively with that shown in Ref. [19]. The reason why the energy dispersion at the Dirac point for N = even is gapped, while it is gapless for
N= odd is due to the change of interference patterns of
the overlap in edge state wave functions between different sublattice structures in these two cases. For N = even, the hard-wall boundary condition breaks one of the sublattice discrete translational invariance along the y direction (see Fig. 1), which opens up a gap in the energy spectrum by introducing a finite interedge hopping [19]. For N = odd, however, this symmetry is preserved and the t⊥term vanishes by destructive quantum interference. These results are further confirmed numerically via Eq. (9) below based on eigenvector analysis of the zigzag KM ribbon.
Based on our numerical results, the edge states are much more localized at the Dirac point ka= ±π: β(k = π/a) > 1 compared to that at other values of k. For 2π/3 < ka < π , however, the edge state wave functions extend over a finite region in the bulk [see Fig.4(b)]. In both cases, a weak but finite overlap between edge and bulk electron wave functions is expected to be present in the zigzag KM ribbon, which generates an effective interedge hoping t⊥term approximately as (see Fig.6and Sec.II B):
Ht⊥ = t⊥ σ=↑,↓ dxc†σ1 c2σ+ H.c. ≈ t⊥ dxe2ikFx(R,†↑1L,↑2+ R,†↓2L,↓1)+ H.c. (7) with x= na and n = ±1, ±2, . . . . The value of t⊥in Eq. (7) can be estimated numerically via diagonalizing the finite-sized KM ribbon: Ht⊥ = t⊥ σ=↑,↓ dx dyc†σ1 (x,y)c2σ(x,y)+ H.c. ≈ t⊥ dy[ ¯R,†↑1(kF,y) ¯L,↑2(kF,y) + ¯R,†↓2(kF,y) ¯L,↓1(kF,y)]+ H.c. (8) ⊥ t t⊥ 1 2 bulk
FIG. 6. (Color online) Schematic diagram for the interedge hop-ping term t⊥(red or blue dashed line).
The Ht⊥ turns out to be important in our RG analysis on
the stability of the helical edge states (see below). The magnitude of t⊥can be estimated via the overlap integral [32] of the opposite edge state wave functions through exact diagonalization of the tight-binding KM model at a finite-sized ribbon [see Eq. (8)] [33]:
t⊥≈ t
L
0
dy[ ¯R,∗↑1(y) ¯L,↑2(y)+ ¯L,∗↓1(y) ¯R,↓2(y)+ c.c.], (9) where we have dropped the kFdependence in ¯L/R,ασ (kF,y) in
Eq. (9). At half-filling, kFa = ±π, hence e2ikFx = 1 and Ht⊥
can in general survive. However, N = even and N = odd lead to different results in this case as explained below.
For N = even, breaking of the sublattice translational invariance at the boundaries results in a finite t⊥. This leads to opening up a gap in the excitation spectrum at the Dirac point when combining Eqs. ((4)) and (7):
(k− π/a) ≈ ±
vF2(k− π/a)2+ (/2)2 (10) with = 2t⊥. We numerically analyzed the gap as shown in Fig.7. The existence of a finite t⊥not only agrees with the energy gap at the Dirac point, it also explains the hybridization of the left and right moving edge states that we found in numerics as the eigenstates of the edge states in the presence of t⊥ are linear combinations of left and right moving edge states. It is clear from Fig.7(a)that the magnitude of the gap decreases with increasing the ribbon size L. In fact, it shows an exponential decay [see Fig.7(b)]:
≈ 0e−αL (11)
with α being the decay constant.
Note that the decay of the small gap was found to be power-law fashion in Ref. [19] by a different (analytical) approach based on the analytical eigenstates for KM model
14π/15 π 16π/15 ka -0.2 -0.1 0 0.1 0.2 / t 4-ZGNR 6-ZGNR 8-ZGNR 10-ZGNR 12-ZGNR 4 6 8 10 12 N 10-10 10-8 10-6 10-4 10-2 100 Δ / t λ/t = 0.1 λ/t = 0.2 λ/t = 0.3 λ/t = 0.4 λ/t = 0.5 (a) (b) ε
FIG. 7. (Color online) (a) Energy spectrum vs momentum k of the topological edge states of the finite-sized (N zigzag chains) Kane-Mele model on a zigzag ribbon of honeycomb lattice near the Dirac point ka= π for different ribbon sizes. Here, we set t = 1, λSO/t= 0.5. (b) Energy gap at the Dirac point as a function of N for different values of λSO.
on 2D honeycomb lattice. With increasing λSO, we find the magnitude of increases with increasing λSO, which comes as a result of the increase in bulk band gap SO. We show in Sec.IVthat this gaped phase corresponds to the charge and spin insulating (or II) phase. In the limit of infinite ribbon width
L→ ∞, the gap vanishes and the gapless Dirac spectrum is
recovered. However, for N= odd, the sublattice translational symmetry at boundaries leads to cancellations in the overlap integral Eq. (9) between sublattices A and B.
At a generic filling away from half-filled, the oscillatory phase factor e2ikFx in t
⊥ term results in cancellations upon averaging over x and Ht⊥ hence vanishes. As shown below,
we also numerically confirmed this result via Eq. (9). Though the Ht⊥ term survives only for N= even and at half-filling,
as shown below, additional two-particle scattering terms are generated via second-order interedge tunnelings, which play an important role in all above-mentioned cases in our stability analysis of the helical edge states in KMH ribbon.
B. The Kane-Mele Hubbard model on a zigzag ribbon Based on the above results for the noninteracting KM model on a finite-sized zigzag ribbon, we now perform an analytical analysis via perturbative RG approach on the weakly interacting KM model (the KM Hubbard model) by including a weak on-site Hubbard U term in HKM. Upon including the on-site Hubbard U term, the Hamiltonian of the Kane-Mele-Hubbard (KMH) model reads
HKMH= HKM+ HU, HU = U dx dy[n↑(x,y)n↓(x,y)], (12)
nσ(x,y)= c†σ(x,y)cσ(x,y).
To simplify our calculations, we consider HKMapproximately as three different contributions: (i) the well-localized edge state Hedge, (ii) the insulating bulk states Hb, and (iii) a weak
coupling between edge and the bulk states Htdue to the
finite-size effect:
HKM≈ Hedge+ Hb+ Ht, (13)
where the edge part Hedgeis defined in Eq. (4), the bulk part
Hbof HKMis given by
Hb=
k,α=↑,↓
HKMcbα(k),c†,αb (k), (14)
and the edge-bulk overlap term Htreads
Ht = t dx[e−ikFxR,†↑1cb,↑1(x)+ eikFx†↓ L,1c ↓ b,1(x) + e−ikFx†↓ R,2c↓b,2(x)+ e ikFx†↑ L,2cb,↑2(x)], (15) where t∼ O(t,λSO). where cσb,1(2)(x) stands for the bulk electron operators near edge 1(2). We further sim-plify the Hubbard U term HU in Eq. (12), and
contributions as HU = HU,e+ HU,b, HU,e= U dx i=1,2 [n↑i(x)n↓i(x)], (16) HU,b= U dx dy[n↑b(x,y)n↓b(x,y)].
Here, i= 1(2) refers to the top (bottom) edge, cbα(k) is
the electron destruction operator in the bulk. Also, the Hb
term, representing the KM model of the bulk electrons, shows an energy dispersion Eb(k) with an energy gap so∼
6√3λSO[17]. For the periodic 2D KM model, Eb(k) has been
shown to be (see Ref. [17])
Eb(k)= ± |gk|2+ γk2, gk= t 3+ 2 cos(√3ky)+ 4 cos( √ 3ky/2) cos(3kx/2), γk= λSO[− sin( √ 3ky)+ 2 cos(3kx/2) sin( √ 3ky/2)]. (17) To simplify our analysis, we assume here the bulk bands are well-separated by the bulk gap SOin the presence of a finite spin-orbit coupling λSO, and|U| λSO. We rewrite the on-site Hubbard U term along the edges, HU,e, by the current operators
defined below for the ease of renormalization group analysis in the bosonization language [34,35]:
JLρ(R)=
i=1,2
JLρ(R),i, JL,ρ1(2)= L,†↓(↑)1(2)L,↓(↑)1(2),
(18)
JR,ρ1(2)= R,†↑(↓)1(2)R,↑(↓)1(2), JLa(R)=x,y,z= L†α(R)σαβa Lβ(R).
In this bases, HU,eis written as
HU,e= Hρ+ Hσz, Hρ = gρ dxJLρJRρ, (19) Hσz= gσz dx JLzJRz,
where JR/Lρ is the U(1) scalar current operator and JLz(R)=
1 2[
†↑
L,2(R,1)L,↑2(R,1)− L,†↓1(R,2)L,↓1(R,2)] is the z component of the SU(2) vector current operator JLa=x,y,z(R) . gρand gσztake the
following bare (initial) values in the context of renormalization group analysis: gρ(μ0)≡ gρ,0= U/2, gσz(μ0)≡ gz,σ0= −2U
with μ0being the bandwidth of the tight-binding KM model. We now turn our attention to Ht term in Eq. (13).
Integrating out the bulk electrons cα
b in Eqs. (14) and (15),
an effective interedge tunneling term Ht⊥as shown in Eq. (7)
is generated where t⊥∼ Dbulk(t)2/SO with Dbulk being the average electron density of states in the bulk. The estimation for t⊥ here can be compared to that in Eq. (9) via numerical diagonalization of the KM ribbon. Note that the interedge hopping t⊥(or the bulk gap SO) is enhanced with increasing spin-orbit coupling λSO: t⊥∝ (t)2/SO ∝ λ2SO/SO∝ λSO [see Fig.7(b)]. Apart from Ht⊥, the linear term in t⊥, for both
half-filling and away from half-filling cases, Ht⊥ term will generate through the second order perturbation theory [36] the following two two-particle scattering terms which turn out to
(a) (b) 1 2 bulk 1 2 bulk gum gσ
FIG. 8. (Color online) Schematic diagrams for (a) the interedge umklapp gum(red and blue arrows) and (b) the interedge spin-flip g⊥σ
processes.
be important in the stability analysis of topological edge states: ˜ Ht⊥ = Hum+ Hσ⊥, Hum = gum dx ei4kFx R,†↑1R,†↓2L,↑2L,↓1 +1 2( †↑ R,1(x)R,†↑1(x+ a)L,↑2(x)L,↑2(x+ a) + R,†↓2(x)R,†↓2(x+ a)L,↓1(x)L,↓1(x+ a)) + H.c. , Hσ⊥ = g⊥σ dx(JL+JR−+ H.c.), (20)
where Hum and Hσ⊥represent for the interedge umklapp and
interedge spin-flip terms, respectively (see Fig. 8), and the transverse components of the SU(2) vector current operators
JL/R+ , JL/R− are defined as
JL+(R)≡ JLx(R)+ i JLy(R)= L,†↑2(R,1)L,↓1(R,2),
(21)
JL−(R)≡ JLx(R)− i JLy(R)= L,†↓1(R,2)L,↑2(R,1).
Similar to Eq. (8), the bare couplings for Hσ⊥ and Hum,
gum(μ0)≡ gum0 and gσ⊥(μ0)≡ gσ⊥,0can be estimated
numeri-cally as gum0 ≈ t 4 L 0 dy
R,∗↑1(y)R,∗↓2(y)L,↑2(y)L,↓1(y) +1
2( ∗↑
R,1(y)R,∗↑1(y)L,↑2(y)L,↑2(y) + R,∗↓2(y) ∗↓ R,2(y) ↓ L,1(y) ↓ L,1(y))+ c.c. , gσ⊥,0≈ t 2 y 0
dy[L,∗↑2(y)L,↓1(y)R,∗↓2(y)R,↑1(y)+ c.c.]. (22) Note that the interedge umklapp term Hum depends sensi-tively on the electron filling factor. At half-filling, ei4kFx = 1,
Hum therefore in general survives. For N = even, we find −g0
um= g⊥,0σ = t⊥2/tvia the energy gap at the Dirac point.
For N= odd, by substituting the edge state wave functions that we numerically obtained based on the tight-binding KM
2 4 6 8 10 12 14 16 N 10-10 10-8 10-6 10-4 10-2 100 t⊥ N=odd t /
FIG. 9. (Color online) The exponential decay of ¯t⊥as a function of odd number of zigzag chains N .
ribbon into Eq. (22), we find−g0um= g⊥,0σ ≡ ¯t⊥2/t, where
¯t2 ⊥≈ t2
L
0
dy|R,↑1(y)|2|L,↑2(y)|2. (23) Note that though the first-order interedge hopping t⊥ term for N = odd vanishes due to sublattice symmetry mentioned above, the second order interedge hopping processes with bare couplings defined in Eq. (23) can, in general, survive due to constructive quantum interference between edge state wave functions in the overlap integral. This result has been confirmed numerically via Eq. (24) based on eigenvector analysis of the KM ribbon. We further find numerically that ¯t⊥ shows an exponential decay with increasing the ribbon width
L, similar to the case for N = even: ¯t⊥∝ e−γkL
(24) with γkbeing the decay constant (see Fig.9). Note that at
half-filling, γk=π/a 1 (or t⊥/t 1) due to the well-localized
edge states.
When the system is away from half-filling, however, the oscillatory factor ei4kFx in H
um leads to cancellations upon summing over x, and therefore Humterm vanishes completely. Nevertheless, Hσ⊥term still survive: g⊥,0σ ≡ ¯t2
⊥/t.
Note that similar two-particle scattering processes Hσ⊥and
Hum terms have been considered in Ref. [22] in the context of the tunneling between helical edge states in a quantum point contact (QPC) as well as in Ref. [37]. However, the authors in Ref. [22] studied the effect of interedge single-and two-particle scattering processes on the helical edge states for a fixed electron-electron interactions (or Luttinger parameter K), while in Ref. [37], the authors did not specify the origins of these two-particle scattering terms. By contrast, the two-particle scatterings we consider here come as a result of second-order interedge tunnelings. Furthermore, we treat the combined effects of the interedge two-particle scatterings Hσ⊥,
Humcontributed from the interedge hopping Ht⊥as well as Hρ,
Hz
σ terms via on-site Hubbard U term in the weak-coupling
limit on equal-footing.
Combining Eqs. (7) and (19)–(21), the effective Hamilto-nian of two weakly coupled helical edge states is therefore
given by
Hedgeeff = Hedge+ Ht⊥+ HU,e+ ˜Ht⊥
= Hedge+ Ht⊥+ Hσ⊥+ H z
σ+ Hρ+ Hum, (25) where Hedge can be re-expressed in terms of the scalar and vector current operators, similar to that for an one-dimensional noninteracting electrons at half-filling [34,35]:
Hedge= dx π 2v c F i=1,2
JL,iρ JL,iρ + JR,iρ JR,iρ
+2π 3 v s F( JL· JL+ JR· JR) (26) with the bare values for the Fermi velocities in the charge and spin sectors given by vc
F = v s
F = vF. Note that our
effective Hamiltonian for the edges Eq. (25) describes two weakly coupled helical Luttinger liquids. In particular, Hedge, describing two noninteracting helical edge states, exhibits
U(1)× SU(2) symmetry; while as the combined transverse and z component of the vector current operator product
Hσ⊥+ Hσz term, describing the couplings between (Hσ⊥) and within (Hσz) the two edges, breaks the SU(2) spin rotational
symmetry down to Z2 symmetry as gum= gσ⊥ in general
[see Eqs. (20) and (21)]. Our effective model for the weakly coupled helical Luttinger liquids Hedge can be characterized as a one-dimensional fermionic Hubbard model with SU(2) spin-anisotropic interactions [17,34,35]. The breaking of the SU(2) symmetry of the model comes as a result of the Hubbard
Uor interedge hopping term at the edges [see Eq. (12)].
III. RG ANALYSIS AND PHASE DIAGRAM OF THE KMH MODEL
We now analyze Eq. (25) via renormalization group approach to understand the stability of the edge states in the presence of Hubbard interactions. Note that the Hamil-tonian (25) is closely related to the spin anisotropic Hubbard model for one-dimensional electrons where electron-electron interactions break the SU(2) symmetry [34,35]. Following the similar RG analysis to Refs. [34,35], we may separate the four couplings (gρ,gum,gσ⊥,g
z
σ) into two pairs belonging to the spin
sector (gz
σ,gσ⊥) and the charge sector (gum,gρ), respectively.
Under RG transformations, these couplings exhibit the prop-erty of spin-charge separation, i.e., the renormalization of the couplings in the spin and charge sectors will remain in its own sector. We shall also analyze the single-particle interedge hopping Ht⊥term under RG. Below, we separately discuss the
RG scaling equations for the half-filled and for a generic filling away from half-filling for both N= even and N = odd.
We summarize our results in Table I. For width of the ribbon L→ ∞, the interedge hopping term vanishes. The relevant perturbations are given by Eq. (19) only. In this two-dimensional limit, the TI phase at half-filling is unstable against the charge gaped but spin gapless (IC) phase [22] for arbitrary U > 0, while it is stable away from half-filling. Details for this case are discussed in Ref. [22].
TABLE I. Summary of RG analysis on the Kane-Mele Hubbard model on a zizag ribbon. “?” refers that the phase transition type is not specified in this work.
Category Phase Transition Critical Point/Line Phase transition type RG flow figure
N= even, at half-filling TI↔II U= t⊥= 0 ? Fig.10
N= even, away from half-filling TI↔CI t2
⊥/t= 2U KT Fig.12
TI↔CI U= ¯t⊥= 0 ?
CI↔II ¯t2
⊥/t= −U/2 KT
N=odd, at half-filling II↔IC ¯t2
⊥/t= 2U KT Fig.15
IC↔TI U= ¯t⊥= 0 ?
TI↔II U= ¯t⊥= 0 ?
N= odd, away from half-filling TI↔CI ¯t2
⊥/t= 2U KT Fig.17
A. N= even 1. At half-filling
As shown previously, at half-filling (kFa = ±π), the KM
model for a finite-sized zigzag ribbon induces a finite interedge hopping term, t⊥= 0. It can be shown that under RG transformation [34], Ht⊥in Eq. (7) is a relevant operator with
scaling dimension [t⊥]= −1. Hence, the RG scaling equation reads [34]
dt⊥
dln μ = −t⊥, (27)
where μ is the running cutoff in energy. Under RG transfor-mation, the running cutoff scale μ is lowered from μ0>0 to zero. It is clear that t⊥flows to a strong coupling fixed point,
t⊥(μ= 0) = ∞. As a result, both gσ⊥and gumbecome relevant under RG as their magnitudes are proportional to t⊥2. When the two-particle spin-flip processes gσ⊥ term becomes relevant, a
spin gap is opening up, while a charge gap develops when the two-particle backscattering gum term becomes relevant. Therefore the t⊥→ ∞ fixed point corresponds to the charge and spin gaped (or charge and spin insulating II) phase (see Fig.10).
2. Away from half-filling
We now proceed to address the case of finite doping away from half-filling, kFa = π. In this case, the interedge
hopping term Ht⊥ and umklapp term Hum vanish due to the oscillatory exponential factors e2ikFx and e4ikFx, respectively
(see Sec.II). The RG scaling equations for both finite-sized and
t2⊥ / t
TI
U / t
II
charge gaped, spin gaped
FIG. 10. Quantum phase diagram of the Kane-Mele Hubbard model at half-filling for N= even as a function of U/t and t2
⊥/t. The helical topological edge states (TI phase) is stable only at U= t⊥= 0 (dark circle). For a finite ribbon size, t⊥= 0, the system flows to a charge and spin gaped (charge and spin insulating or II) phase.
infinite-sized ribbons are reduced to [34,35]
dgρ
dln μ = 0, (28)
in the charge sector with g0
ρ= U and dg⊥σ dln μ = −g ⊥ σg z σ, dgz σ dln μ = −(g ⊥ σ)2, (29)
in the spin sector with (gσz,0,g⊥,0σ )= (−2U,t⊥2/t). Due to the
absence of Ht⊥and Humterms in this case, the coupling gρis
marginal under RG up to the second order in gum[see Eq. (32)], which is at the level of accuracy in our RG analysis for all other couplings.
Via Eq. (28), it is clear that the system will not develop a charge gap under RG as gρdoes not diverge: gρ(μ)= gρ0 1.
The RG flows in the spin sector, however, suggest that the topological edge states may undergo the Kosterlitz-Thouless transition upon increasing t⊥to a charge gapless but spin gaped (CI) phase characterized by the following fixed point:
CI : gσ⊥,0+ g z,0 σ >0, g z σ(μ→ 0), gσ⊥(μ→ 0) → ∞, gρ(μ→ 0) = 0, gum(μ→ 0) = gum0 1. (30) The TI-CI phase boundary is set by the separatrix g⊥σ + g
z σ = 0
(or when t2
⊥/t = 2U, see Fig.11). The helical edge states are therefore stable for t2
⊥/t <2U , while it is unstable against the CI phase for U <t⊥2
2t. Combing RG flows in both charge and spin sectors, this spin gaped phase corresponds to the charge conducting but spin insulating (or CI) phase (see Fig.12).
B. N= odd 1. At half-filling
At half-filling, kFa= ±π and t⊥= 0, all the four couplings
(gρ,gum,gσ⊥,gσz) exist in general under RG transformations.
Their initial (bare) couplings at μ= μ0 are given by: (gum0 ,gρ0)= (−¯t⊥2/t,U/2), (gσz,0,g⊥,0σ )= (−2U,¯t⊥2/t). The RG
scaling equations in this case can be casted in a spin-charge separated form [34,35] and are readily obtained via the operator product expansion (OPE) for the current algebra in the one-dimensional Hubbard model with broken SU(2) symmetry (see, for example, Appendix in Chap. 17 of Ref. [34]):
dgρ dln μ = −g 2 um, dgum dln μ = −gumgρ, (31)
g g TI σ z σ CI t 2 (−2U, t / )
FIG. 11. (Color online) The RG flows of the Kosterlitz-Touless type for the spin sector (g⊥σ,gz
σ) of the zigzag Kane-Mele Hubbard
ribbon for N = even away from half-filling. The black circle stands for the initial (bare) couplings. The arrows indicate the directions of the RG flows upon decreasing the curt-off scale μ from μ0. The red line represents a line of fixed points in the TI phase, the TI-CI phase boundary is defined by the separatrix line (thick black arrow). Note that the coupling gρdoes not flow under RG in this case [see
Eq. (34)].
in the charge sector and
dgσ⊥ dln μ = −g ⊥ σg z σ, dgσz dln μ = −(g ⊥ σ)2, (32)
in the spin sector.
As shown in Figs. 13 and 14, the generic RG flows of Eqs. (31) and (32) are of the Kosterlitz-Thouless (KT) type. In the charge sector, the RG flows for gum and gρ with the
bare couplings (g0
um,gρ0)= (−¯t⊥2/t,U/2) are always towards
either the strong-coupling charge and spin gaped II phase for −2¯t2
⊥/t < U < 12¯t 2
⊥/t or towards the charge conducting and spin insulating CI phase for U <−2¯t2
⊥/t. Similarly, in the spin sector, the TI phase is unstable against either the II phase
t2 t2 U / t / t KT CI TI 0 / t = 2 U
FIG. 12. (Color online) Quantum phase diagram of the Kane-Mele Hubbard model away from half-filling for N= even as functions of t2
⊥/t and U/t. The helical topological edge states (TI phase) are unstable towards the charge conducting and spin insulating CI phase for t2
⊥/t >2U . The TI-CI quantum phase transition set by the boundary t2
⊥/t= 2U is of the Kosterlitz-Thouless (KT) type (red dashed arrows). gρ −g um 2 _ t (U/2, t / ) TI CI II
FIG. 13. (Color online) The RG flows of the Kosterlitz-Touless type for the charge sector (gρ,gum) of the zigzag Kane-Mele Hubbard ribbon at half-filling for N= odd. The black circle stands for the initial (bare) couplings at (g0
ρ,−g
0
um)= (U/2,¯t⊥2/t). The arrows indicate the directions of the RG flows upon decreasing the curt-off scale μ from μ0. The red line represents a line of fixed points in the CI phase, the CI-II boundary is defined by the separatrix line (thick black arrow) and its quantum transition is of the Kosterlitz-Thouless (KT) type. The topological TI phase is stable only at the origin U= 0 = ¯t⊥.
for−gz,0
σ < gσ⊥,0(i.e., ¯t⊥2/t >2U ) or against a charge gaped
but spin gapless IC phase for−gz,0
σ > gσ⊥,0(i.e., ¯t⊥2/t <2U )
(see Fig. 14). Therefore the TI phase is unstable against any infinitesimal U= 0 and ¯t⊥= 0. The II-IC and II-CI quantum phase transitions are of the KT type. Combining the RG flows for both spin and charge sectors, we obtain the global phase diagram shown in Fig.15for N = odd and at
g gσ z σ _ t 2 (−2U, t / ) TI IC II
FIG. 14. (Color online) The RG flows of the Kosterlitz-Touless type for the spin sector (gσ⊥,gz
σ) of the zigzag Kane-Mele Hubbard
ribbon for N= odd at half-filling. The black circle stands for the initial (bare) couplings at (gz,0
σ ,g⊥,0σ )= (−2U,¯t
2
⊥/t). The arrows indicate the directions of the RG flows upon decreasing the curt-off scale μ from μ0. The red line represents a line of fixed points in the IC phase, the II-IC phase boundary is defined by the separatrix line (thick black arrow) and its quantum transition is of the Kosterlitz-Thouless (KT) type. The topological TI phase is stable only at the origin U= 0 = ¯t⊥.
t2 t2 _ /t = 2 U t2 _ /t = − U /2 KT U / t / t KT IC _ CI II TI
FIG. 15. (Color online) Quantum phase diagram of the zigzag Kane-Mele Hubbard ribbon for N= odd at half-filling as a function of U/t and ¯t2
⊥/t. The helical topological edge states (TI) are unstable against any U= 0 or ¯t⊥= 0, and towards the IC, CI, and II phases for U > ¯t2
⊥/(2t), U <−2¯t⊥2/t and−2¯t⊥2/t < U < 12¯t 2
⊥/t, respectively. The II-IC and II-CI phase transitions are of the Kosterlitz-Thouless (KT) type (red dashed arrows).
half-filling: II :−2¯t⊥2/t < U < 1 2¯t 2 ⊥/t, gρ(μ→ 0), gum(μ→ 0) → ∞, ¯t2 ⊥/t >2U, gσz(μ→ 0), gσ⊥(μ→ 0) → ∞; IC : U > 1 2¯t 2 ⊥/t >0, gσz(μ→ 0), g⊥σ(μ→ 0) → 0, (33) gρ(μ→ 0), gum(μ→ 0) → ∞; CI : U <−2¯t⊥2/t <0, gσz(μ→ 0), g⊥σ(μ→ 0) → ∞, gρ(μ→ 0), gum(μ→ 0) 1.
The topological edge states (TI) are unstable against the charge and spin insulating II phase for−2¯t2
⊥/t < U < 12¯t⊥2/t, against the charge insulating abd spin conducting IC phase for
U > 1
2¯t 2
⊥/t >0, and against the charge conducting and spin insulating CI phase for U <−2¯t⊥2/t <0. Therefore the TI phase is unstable for any U = 0 or ¯t⊥= 0. The II-IC and II-CI quantum phase transitions are of the KT type. Our results on the stability of the TI phase for KM Hubbard model on a zigzag ribbon are different from those in Ref. [22] through bosonizing the infinite-sized helical Luttinger liquid at a fixed interaction strength set by the Luttinger parameter K =
1− U 2π vF 1+ U 2π vF . There, they showed that TI is stable for 1/2 < K < 2. The difference lies in the fact that the interedge tunneling t⊥ arising from the finite-size effect plays an important role here while it was absent in Ref. [22].
2. Away from half-filling
We now proceed to address the case of finite doping away from half-filling, kFa = π. In this case, the umklapp term Hum vanishes as mentioned in Sec. II. The RG scaling equations
g g TI σ z σ CI 2 _ t (−2U, t / )
FIG. 16. (Color online) The RG flows of the Kosterlitz-Touless type for the spin sector (gσ⊥,gz
σ) of the zigzag Kane-Mele Hubbard
ribbon away from half-filling for N= odd. The black circle stands for the initial (bare) couplings. The arrows indicate the directions of the RG flows upon decreasing the curt-off scale μ from μ0. The red line represents a line of fixed points in the TI phase, the TI-CI phase boundary is defined by the separatrix line (thick black arrow). Note that the coupling gρ does not flow under RG in this case [see
Eq. (34)]. reduce to
dgρ
dln μ = 0, (34)
in the charge sector with g0
ρ= U and dg⊥σ dln μ = −g ⊥ σg z σ, dgz σ dln μ = −(g ⊥ σ)2, (35)
in the spin sector with (gz,0
σ ,g⊥,0σ )= (−2U,¯t⊥2/t). Note that
via the same argument that leads to Eq. (29), the coupling gρ
here is marginal under RG up to the second order in gum. Via Eq. (34), it is clear that the system will not develop a charge gap under RG as gρdoes not diverge: gρ(μ)= gρ0 1.
The RG flows in the spin sector [see Eq. (35)], however,
t2 t2 U / t / t KT CI TI 0 _ _ / t = 2 U
FIG. 17. (Color online) Quantum phase diagram of the zigzag Kane-Mele Hubbard ribbon away from half-filling for N= odd as functions of ¯t2
⊥/t and U/t. The helical topological edge states (TI phase) are unstable towards the charge conducting and spin insulating CI phase for ¯t2
⊥/t >2U . The TI-CI quantum phase transition set by the boundary ¯t2
⊥/t= 2U is of the Kosterlitz-Thouless (KT) type (red dashed arrows).
suggest that the topological edge states may undergo the Kosterlitz-Thouless transition upon increasing ¯t⊥ to a spin gaped phase. Combing RG flows in both charge and spin sectors, this spin gaped phase corresponds to the charge conducting but spin insulating (or CI) phase (see Fig.16). The TI-CI phase boundary is set by the separatrix gσ⊥+ gz
σ= 0
(or when ¯t⊥2/t= 2U, see Fig.16). The helical edge states are therefore stable for ¯t⊥2/t <2U , while it is unstable against the CI phase for U < ¯t⊥2
2t (see Fig.17).
IV. INSTABILITIES, ORDERINGS, AND CORRELATION FUNCTIONS OF THE KANE-MELE-HUBBARD MODEL
We now investigate further the nature of the TI, CI, IC, and II phases. In particular, we focus on instabil-ities towards various orderings and correlation functions in these phases. Various correlation functions with spe-cific orderings can be defined for this purpose: (i) the charge-density-waveOCDW correlation, (ii) the spin-density-wave OaSDW=x,y,z correlation, and (iii) the singlet OSS and triplet OaTS=x,y,z superconducting pairing operators, where [35] OCDW= R,†↑1(x)L,↑2(x)+ R,†↓2(x)L,↓1(x), Ox SDW= R,†↑1(x) ↓ L,1(x)+ †↓ R,2(x) ↑ L,2(x), Oy SDW= −i[R,†↑1(x)L,↓1(x)− R,†↓2(x)L,↑2(x)], Oz SDW= †↑ R,1(x) ↑ L,2(x)− †↓ R,2(x) ↓ L,1(x), (36) OSS= R,†↑1(x) †↓ L,1(x)+ †↑ L,2(x) †↓ R,2(x), Ox TS= R,†↑1(x)L,†↑2(x)+ L,†↓1(x)R,†↓2(x), Oy TS= −i[R,†↑1(x)L,†↑2(x)− L,†↓1(x)R,†↓2(x)], Oz TS= R,†↑1(x)L,†↓1(x)− L,†↑2(x)R,†↓2(x).
Note that some of the operators defined above involve helical electrons on both edges, different from those defined for a standard Luttinger liquid in one-dimensional interacting elec-trons where all elecelec-trons are along the same one-dimensional wire [34,35]. To investigate the above correlation functions, it is useful to bosonize the Hamiltonian Eq. (25) as [22]
Hedgeeff = dx α=c,s vα 2 Kα(∂xα)2+ 1 Kα (∂xα)2 + t⊥ 2π a0 cos(√2π c+ 2kFx) cos( √ 4π s) + g⊥σ (2π a0)2 cos(2√2π s)+ 1 8π gz σ (2π a0)2 ×(∂xs− ∂xs)+ gum (2π a0)2 cos(2√2π c+ 4kFx) +1 4 gρ (2π a0)2 α=c,s ((∂xα)2− (∂xα)2) ,
where via bosonization formulas [22,34,35],
Lσ = 1 √ 2π a0 ησe−i √ 4π φLσ , Rσ = 1 √ 2π a0 ησei √ 4π φRσ, and the bosonic fields defined as
σ = φLσ+ φRσ , σ = φLσ− φRσ , c(s)= 1 √ 2(φ↑± φ↓) , c(s)= 1 √ 2(↑± ↓). with ησbeing the Klein factor and a0being the short-distance cutoff. In terms of these boson fields, the correlation functions mentioned above are given by
OCDW= e−2ikFx π a0 e −i√2π ccos(√2π s), Ox SDW= e−2ikFx π a0 e −i√2π ccos(√2π s), Oy SDW= − e−2ikFx π a0 e −i√2π c sin(√2π s), Oz SDW= i e−2ikFx π a0 e−i √ 2π csin(√2π s), (37) OSS= 1 π a0e i√2π ccos(√2π s), Ox TS= 1 π a0 ei√2π ccos(√2π s), Oy TS= − 1 π a0e i√2π csin(√2π s), Oz TS= 1 π a0 ei√2π csin(√2π s).
Based on the phase diagram via weak-coupling RG and the bosonized form of the Hamiltonian, we analyze below the instabilities and the behaviors of various correlation functions for (i) the charge and spin gapless (TI) topological edge states, (ii) the CI phase, (iii) the IC phase, and (iv) the II phase.
1. The topological edge states (TI) phase
In the gapless topological edge states—the charge and spin conducting state—various correlation functions can be computed via correlation functions of the boson fields, given by O†CDW(0)OCDW(r) ∼ e−2ikFx 1 r Kc+Ks ∼ e−2ikFx 1 r 1/K+K , O†xSDW(0)O x SDW(r) ∼ e−2ikFx 1 r Kc+1/Ks ∼ e−2ikFx 1 r 2K ,
OSDW†y (0)O y SDW(r) ∼ e−2ikFx 1 r Kc+1/Ks ∼ e−2ikFx 1 r 2K , (38) OSDW†z (0)OSDWz (r)∼ e−2ikFx 1 r Kc+Ks ∼ e−2ikFx 1 r 1/K+K , OSS† (0)OSS(r) ∼ 1 r 1/Kc+Ks ∼ 1 r 2/K , O†xTS(0)O x TS(r) ∼ 1 r 1/Kc+1/Ks ∼ 1 r K+1/K , O†yTS(0)O y TS(r) ∼ 1 r 1/Kc+1/Ks ∼ 1 r 1/K+K , O†zTS(0)O z TS(r) ∼ 1 r 1/Kc+Ks ∼ 1 r 1/(2K)
with Kc= K and Ks = 1/K in the helical Luttinger
liq-uid [22]. Note that in the conventional spinful Luttinger liquids where Ks= 1, the above correlation functions get modified
accordingly [23,34,35].
2. The CI phase
Now, we analyze instability and correlation functions in the charge conducting and spin insulating (CI) phase. As shown in Eqs. (30) and ((37)), gσ⊥,gσz → ∞ while gρ,gum→ 0 in this phase. In the bosonized form of the Hamiltonian, this implies that s is pinned to a constant value [22,34,35]:
s ∼ nπ/
√
8π . As a result, its conjugate variable s is
disordered and exhibit exponentially decaying correlation functions [34,35]. The corresponding leading correlation functions have the following power-law behaviors:
O†CDW(0)OCDW(r) ∼ 1 r Kc ∼ 1 r K , (39) O†SS(0)OSS(r) ∼ 1 r 1/Kc ∼ 1 r 1/K .
Note that due to the disordered nature of the s field, the
SDW as well as the TS orderings vanish: O†x,y,zSDW Ox,y,zSDW → 0,O†x,y,zTS O
x,y,z
TS → 0. Therefore we find the leading instabil-ities of the CI phase are towards the CDW and superconduc-tivity (SC). For repulsive interactions K < 1 (or U > 0) that we consider here, the CDW order is dominating over the SC order as CDW correlators decay more slowly than that for SC orders. However, for attractive interactions K > 1 (or U < 0), it is the SC order that dominates the CI phase.
3. The IC phase
We now analyze the instability of the charge insulating but spin conducting (IC) phase. It is clear from Eq. (37) that
c field is pinned to a constant value in this phase: c∼
nπ/√8π . The correlation functions for the CDW and SDW
orderings are given by
O†CDW(0)OCDW(r) ∼ 1 r Ks ∼ 1 r 1/K , O†xSDW(0)O x SDW(r) ∼ 1 r 1/Ks ∼ 1 r K , (40) O†ySDW(0)O y SDW(r) ∼ 1 r 1/Ks ∼ 1 r K , O†zSDW(0)O z SDW(r) ∼ 1 r Ks ∼ 1 r 1/K .
On the other hand, due to the pinning of the c field,
its conjugate field c is completely disordered. Hence
the SS and TS orderings are suppressed: OSS† OSS → 0, O†x,y,zTS O
x,y,z
TS → 0. For repulsive Hubbard term U > 0 (or
K <1), the SDW orderings along x− and y− directions
are the leading instabilities of this phase as their correlation functions decay more slowly compared to the others. The system shows quasi-long-ranged magnetic order. This phase shares similarities to the Mott insulating phase in the sense that interactions lead to a metal-insulator transition and at the same time to a state with magnetic order. In fact, this phase corresponds to the SDW phase found in the mean-field approach of the KM Hubbard in Ref. [17]. For the attractive Hubbard model U < 0 (or K > 1), however, the leading instabilities go towards the CDW and SDW along the z axis.
4. The II phase
Finally, we analyze the charge and spin insulating II phase. This phase occurs for a finite-sized ribbon at half-filling where all the couplings—the interedge hopping term t⊥, the umklapp term gum, scalar density-density interaction gρ, the
two-particle spin scattering terms gσ⊥,z–become relevant under
RG, t⊥,gum,ρ,g⊥,zσ → ∞. From the bosonized Hamiltonian
Eq. (37), this phase requires the pinning of both c and
s fields at c,s ≈ nπ/
√
2π , leading to exponential decay of all the correlation functions associated with the orderings in Eq. (36) except for the CDW ordering with a constant correlator. Whether or not the II phase found here is related to the gaped, charge-gaped (similar to II phase) spin-liquid phase found numerically via QMC in Refs. [14,15] or furthermore to the Anderson’s resonant-valence-bond (RVB) spin liquid need further investigations.
V. DISCUSSIONS AND CONCLUSIONS
Before we conclude, we would like to make a remark on the possible realization of our system in experiments. Though graphene has been originally proposed to be a candidate for QSHIs [3,8,9], its negligible SO coupling makes it unpractical. Nevertheless, there has been proposals based on density functional theory and tight-binding simulations to significantly increase the intrinsic spin-orbit coupling of the KM type in graphene by doping heavy adatoms, such as indium or thallium [38]. The KM Hubbard model and our results here are therefore relevant for these adatom-doped graphene where
a weak Coulomb repulsive interaction is expected and the doping level is controllable via applying gate voltages.
In summary, we have studied the stability of the helical edge states and quantum phases and phase transitions of the Kane-Mele Hubbard (KMH) model on a finite-sized zigzag ribbon of honeycomb lattice. We first focused on the finite-size effect of the Kane-Mele (KM) zigzag ribbon in the absence of the on-site Hubbard interaction. We reproduced in the energy excitation spectrum the well-known Dirac-dispersed topological edge states. In additions, due to the finite ribbon size, we have shown that a finite interedge hopping between two edge states exist, which falls off exponentially with increasing ribbon width. This interedge hopping term generates via second-order perturbation two important two-particle scatterings: the interedge spin-flip and the interedge backscattering (or the umklapp) terms. These three terms lead to instabilities of the topological edge states. We further analyzed the instabilities of the topological edge states, as well as possible quantum phases and phase transitions upon including a weak on-site repulsive Hubbard interaction on the zigzag KM ribbon. Via perturbative RG approach we found the combined effects from the interedge hopping and the on-site Coulomb interactions lead to the instabilities of the topological
edge states (TI phase) against (i) the charge and spin insulating II phase, (ii) the charge insulating but spin conducting IC phase, and (iii) the charge conducting but spin insulating CI phase, depending on N= even/odd, the electron density (filling factor), and on the ratio of the Coulomb interaction
Uand the interedge tunneling t⊥, U/t⊥. Via RG analysis we found that the quantum phase transitions for TI-CI, II-IC, and II-CI are of the Kosterlitz-Thouless type. Via bosonization approach, we furthermore investigated the instabilities towards new orderings, including the CDW, SDW and superconducting orders by computing correlation functions of these orderings in the helical edge states, as well as in the CI, IC, and II phases. Our theoretical results can serve as a basis to investigate further both theoretically and experimentally correlation effects or Mott physics in interacting topological insulators.
ACKNOWLEDGMENTS
We acknowledge M. Cazalilla and C. Y. Mou for helpful discussions. This work is supported by the NSC Grant Nos. 98-2918-I-009-06 and 98-2112-M-009-010-MY3, the NCTU-CTS, the MOE-ATU program, the NCTS of Taiwan, Republic of China.
[1] M. Z. Hasan and C. L. Kane, Rev. Mod. Phys. 82, 3045 (2010).
[2] X. L. Qi and S. C. Zhang, Rev. Mod. Phys. 83, 1057 (2011).
[3] C. L. Kane and E. J. Mele,Phys. Rev. Lett. 95,146802(2005). [4] B. A. Bernevig, T. L. Hughes, and S. C. Zhang,Science 314,
1757(2006).
[5] L. Fu and C. L. Kane,Phys. Rev. B 76,045302(2007). [6] H. Zhang, C.-X. Liu, X.-L. Qi, X. Dai, Z. Fang, and S. C. Zhang,
Nat. Phys. 5,438(2009).
[7] M. K¨onig, S. Wiedmann, C. Br¨une, A. Ruth, H. Buhmann, L. W. Molenkamp, X.-L. Qi, and S. C. Zhang,Science 318, 766(2007); D. Hsieh, D. Qian, L. Wray, Y. S. HOr, R. J. Cava, and M. Z. Hasan,Nature (London) 452,970(2008); D. Hsieh, Y. Xia, L. Wray, D. Qian, A. Pal, J. H. Dil, F. Meier, J. Osterwalder, G. Bihlmayer, C. L. Kane, Y. S. Hor, R. J. Cava, and M. Z. Hasan,Science 323,919(2009); Y. Xia, D. Qian, D. Hsieh, L. Wray, A. Pal, H. Lin, A. Bansil, D. Grauer, Y. S. Hor, R. J. Cava, and M. Z. Hasan,Nat. Phys. 5,398(2009); Y. L. Chen, J. G. Analytis, J.-H. Chu, Z. K. Liu, S.-K. Mo, X. L. Qi, H. J. Zhang, D. H. Lu, 1 X. Dai, Z. Fang, S. C. Zhang, I. R. Fisher, Z. Hussain, and Z.-X. Shen,Science 325,178(2009); P. Roushan, J. Seo, C. V. Parker, Y. S. Hor, D. Hsieh, D. Qian, A. Richardella, M. Z. Hasan, R. J. Cava, and A. Yazdani,Nature (London) 460,1106(2009); D. Hsieh, Y. Xia, D. Qian, L. Wray, J. H. Dil, F. Meier, L. Patthey, J. Osterwalder, A. V. Fedorov, H. Lin, A. Bansil, D. Grauer, Y. S. Hor, R. J. Cava, and M. Z. Hasan,ibid. 460,1101(2009).
[8] F. D. M. Haldane,Phys. Rev. Lett. 61,2015(1988).
[9] C. L. Kane and E. J. Mele,Phys. Rev. Lett. 95,226801(2005). [10] B. A. Bernevig and S. C. Zhang,Phys. Rev. Lett. 96,106802
(2006).
[11] L. Fu, C. L. Kane, and E. J. Mele,Phys. Rev. Lett. 98,106803 (2007).
[12] C. Wu, B. A. Bernevig, and S. C. Zhang,Phys. Rev. Lett. 96, 106401(2006); C. Xu and J. E. Moore,Phys. Rev. B 73,045322 (2006).
[13] X.-L. Qi and S. C. Zhang,Phys. Today 63,33(2010).
[14] Z. Y. Meng, T. C. Lang, S. Wessel, F. F. Assaad, and A. Muramatsu,Nature (London) 464,847(2010).
[15] M. Hohenadler, Z. Y. Meng, T. C. Lang, S. Wessel, A. Muramatsu, and F. F. Assaad,Phys. Rev. B 85,115132(2012). [16] Z. Y. Meng, H.-H. Hung, and T. C. Lang,Mod. Phys. Lett. B
28,1430001(2014).
[17] S. Rachel and K. Le Hur,Phys. Rev. B 82,075106(2010). [18] M. Hohenadler and F. F. Assaad,Phys. Rev. B 85, 081106
(2012).
[19] M. Zarea, C. B¨usser, and N. Sandler,Phys. Rev. Lett. 101, 196804(2008).
[20] D. Zheng, G. M. Zhang, and Congjun Wu,Phys. Rev. B 84, 205121(2011).
[21] S-L. Yu, X. C. Xie, and J.-X. Li,Phys. Rev. Lett. 107,010401 (2011).
[22] Jeffrey C. Y. Teo, and C. L. Kane,Phys. Rev. B 79, 235321 (2009).
[23] B. Braunecker, C. Bena, and P. Simon,Phys. Rev. B 85,035136 (2012).
[24] Jun Wen, Mehdi Kargarian, Abolhassan Vaezi, and Gregory A. Fiete,Phys. Rev. B 84,235149(2011).
[25] P. A. Maksimov, A. V. Rozhkov, and A. O. Sboychakov,Phys. Rev. B 88,245421(2013).
[26] S. Sachdev, Quantum Phase Transitions (Cambridge University Press, Cambridge, UK, 2000); S. L. Sondhi, S. M. Girvin, J. P. Carini, and D. Shahar,Rev. Mod. Phys. 69,315(1997).
[27] K. Nakada, M. Fujita, G. Dresselhaus, and M. S. Dresselhaus, Phys. Rev. B 54,17954(1996).
[28] M. Fujita et al.,J. Phys. Soc. Jpn. 65,1920(1996). [29] L. Brey and H. A. Fertig,Phys. Rev. B 73,235411(2006). [30] Y. W. Son, M. L. Cohen, and S. G. Louie,Phys. Rev. Lett. 97,
216803(2006).
[31] V. N. Do and T. H. Pham,Adv. Nat. Sci.: Nanosci. Nanotechnol. 1,033001(2010).
[32] N. W. Ashcroft and D. N. Mermin, Solid State Physics paperback edition (Thomson Press, India, 2003).
[33] D. H. Lee and C. H. Chung (unpublished).
[34] A. O. Gogolin, A. A. Nersesyan, and A. M. Tsvelik, Bosoniza-tion and Strongly Correlated Systems (Cambridge University Press, Cambridge, 1998).
[35] T. Giamarchi, Quantum Physics in One Dimension (Oxford University Press, Oxford, 2004); J. von Delft and H. Schoeller, Ann. Phys. 7,225(1998).
[36] The derivation of Eq. (20) involves going back from second quantization to the first quantization representation in the edge state model Hamiltonian (5) and treats interedge hopping Eq. (8) as the perturbation term Ht⊥. Energy shift in the second-order
perturbation givesH2
t⊥/H0 ∼ t⊥2/tas the hopping along the same edge is determined by scale t. Writing first quantization results back to second quantization gives Eq. (21).
[37] Y. Tada, R. Peters, M. Oshikawa, A. Koga, N. Kawakami, and S. Fujimoto,Phys. Rev. B 85,165138(2012).
[38] C. Weeks, J. Hu, J. Alicea, M. Franz, and R. Wu,Phys. Rev. X 1,021001(2011).
