Table Of Content

Phase Diagram and Entanglement of two interacting topological Kitaev chains Lo¨ıc Herviou,1,2 Christophe Mora,2 and Karyn Le Hur1 1Centre de Physique Théorique, Ećole Polytechnique, CNRS, Université Paris-Saclay, 91128 Palaiseau, France 2Laboratoire Pierre Aigrain, Ećole Normale Supérieure-PSL Research University, CNRS, Université Pierre et Marie Curie-Sorbonne Universités, Université Paris Diderot-Sorbonne Paris Cité, 24 rue Lhomond, 75231 Paris Cedex 05, France (Dated: April 5, 2016) A superconducting wire described by a p-wave pairing and a Kitaev Hamiltonian exhibits Ma- jorana fermions at its edges and is topologically protected by symmetry. We consider two Kitaev wires (chains) coupled by a Coulomb type interaction and study the complete phase diagram us- 6 ing analytical and numerical techniques. A topological superconducting phase with four Majorana 1 fermions occurs until moderate interactions between chains. For large interactions, both repulsive 0 and attractive, by analogy with the Hubbard model, we identify Mott phases with Ising type mag- 2 netic order. For repulsive interactions, the Ising antiferromagnetic order favors the occurrence of r orbital currents spontaneously breaking time-reversal symmetry. By strongly varying the chemi- p cal potentials of the two chains, quantum phase transitions towards fully polarized (empty or full) A fermionic chains occur. In the Kitaev model, the quantum critical point separating the topological superconducting phase and the polarized phase belongs to the universality class of the critical 3 Ising model in two dimensions. When increasing the Coulomb interaction between chains, then we identifyanadditionalphasecorrespondingtotwocriticalIsingtheories(ortwochainsofMajorana ] l fermions). We confirm the existence of such a phase from exact mappings and from the concept of e bipartite fluctuations. We show the existence of negative logarithmic corrections in the bipartite - r fluctuations, as a reminiscence of the quantum critical point in the Kitaev model. Other entangle- t mentprobessuchasbipartiteentropyandentanglementspectrumarealsousedtocharacterizethe s . phase diagram. The limit of large interactions can be reached in an equivalent setup of ultra-cold t a atoms and Josephson junctions. m - d I. INTRODUCTION model,inthesensethattheseinteractionswillbepresent n in most systems. This ladder is also a building step to- o wards building two-dimensional materials. c The discovery of topological effects in quantum solid The simplest model of superconducting wire is the well- [ stateshasrevolutionizedcondensedmatterstudies. From known and exactly solvable Kitaev’s wire2 : 2 the fractional Quantum Hall effect to topological super- v conductors, topological Hamiltonians have been a topi- L L−1 (cid:88) (cid:88) 8 cal subject of study during the last few years. One of HK{c}=−µ c†jcj + −t(c†jcj+1+c†j+1cj) 9 the main characteristics of such systems is the presence j=1 j=1 9 2 of anionic excitations, with fractional spin and charge. +∆(c†jc†j+1+cj+1cj). (1) 0 In particular, in the case of topological superconductors, . zero-energyMajoranamodesappearattheedgesofthese Here, c is a fermionic spinless annihilation operator, t is 1 systems or in vortices1,2. It is of particular interest in the hopping amplitude (it will serve as an energy scale 0 6 conjunction with the rise of quantum information3. The in the rest of the paper), j is a site index and ∆ is 1 topological nature of these modes, preventing any effect a Bardeen-Cooper-Schrieffer (BCS)11 p-wave supercon- : from small local interactions, protects them from deco- ducting pairing term generated by an interaction with a v herence and make them perfect candidates for the real- superconducting substrate. L is the number of sites in i X ization of quantum bits4. Indeed, several schemes have the wire. There has been several proposals and realiza- r been proposed to realize complete sets of quantum gates tions of this model, for example by coupling a semicon- a and memories, using superconducting wires with Majo- ductingnanowiretothebulkoftwo-orthree-dimensional ranafermionsateachextremities5–7. Theseschemesrely superconductors via a strong spin-orbit interaction and on controlled interactions between several of such wires. by applying a magnetic field to select one spin species The natural question concerns the effects of other un- in the wire7,12–15. Other implementations have been dis- controlled interactions that could arise due to the prox- cussed with ferromagnetic metallic chains16,17 and ultra- imity of these wires. Lately, numerous propositions on cold atoms18,19. Majorana fermions can also occur as a topological systems presenting solvable points have been result of purely intrinsic attractive interactions20. This made8–10. We present in this paper the general study of modelpresentsaZ topologicaldegeneracyinitsground 2 two topological superconducting wires in the presence of state, corresponding to a free Majorana fermion subsist- thesimplestinteraction,aCoulomb-likeinteractionmod- ing at each extremity of the wire. eled by an on-site repulsion à la Hubbard. It is a generic WeconsidertwointeractingKitaevchains,coupledvia 2 a Coulomb interaction: (cid:88) 1 1 10 H =g (n − )(n − ), (2) int j,1 2 j,2 2 MI−AF j g/t 5 nj,1/2 istheelectronnumberoperatorinthefirst/second h t wire at site j. Interpreting the chain index as a spin in- g n dex, then this can be identified as the well-known Hub- re 4MF t bard interaction, a staple of condensed matter physics n s 0 thoroughly studied for the last 40 years. This interac- tio c tiondoesnotbreakanyofthediscretesymmetriesofthe a er Polarized Polarized original problem. Indeed, while Kitaev’s model does not nt 5 conserve the number of fermions, it preserves fermionic I parityandhastime-reversalsymmetry(andparticle-hole MI F − symmetryforzerochemicalpotential). Introducingasec- 10 ond wire, though, has non trivial effect on the topology: 4 2 0 2 4 Chemical potential µ/t followingtheclassificationoftopologicalphasesproposed byFidkowskiandKitaev21,22,weknowthatthetopolog- FIG.1. (Coloronline)Sketchofthephasediagramofthein- icalnatureofasetofnidenticalwiresininteractionsde- teractingladderat∆=tobtainedwithanalytical(bosoniza- pendsonn. Inparticular,with2wires,weknowthatthe tion, exact mappings) and numerical methods (exact diag- system is only a symmetry protected topological phase onalization (ED), Density Matrix Renormalization Group (SPT). Hence, an arbitrarily small term breaking one of (DMRG)).Thechemicalpotentialsofthetwochainsaretaken theoriginalsymmetriesofthemodel,inourcasetimere- tobeequal. 4MF istheSPT gappedphasepresenting2Ma- versalsymmetry,liftsthedegeneracyofthegroundstate. jorana fermions at each of the ladder extremities. MI−AF A large variety of interactions can be considered. In and MI −F are two gapped Mott phases, either antiferro- particular, fine-tuned interacting terms have been added magneticorferromagnetic. Polarizedcorrespondstoatrivial to make the Kitaev ladder exactly solvable8,9. Supple- phasewithaquasi-emptyorquasi-fullladder. Inred,thegap- lessDCI phaseembodiesanextensionofthecriticalpointat mentary terms could be considered such as introducing g=0. It acquires an extension of order t as g goes to +∞. a hopping term between the two wires −t c†c +h.c or ⊥ 1 2 allowing for an orthogonal pairing term ∆ c†c† +h.c23. ⊥ 1 2 Whilethesetermscanbetunedtoobtainexactlysolvable points,theyalsobecomenegligibleforalargeenoughdis- lable interaction could be achieved between and inside tancedbetweenthewires. TheCoulombrepulsionscales the two wires. like1/d2. Thehoppingamplitude,scalingasexp(−d/χ), This problem can also be mapped in terms of two in- withχbeingacorrelationlength,isnegligibleford(cid:29)χ. teracting Ising spin-1 chains36: 2 Similarly, if d is larger than the coherence length of the  Cooper pair, one can safely ignore ∆⊥, as long as both (cid:88)2 (cid:88)L L(cid:88)−1(∆−t) these terms do not break the time-reversal symmetry, Hspin = − µσjz,w+ 2 σjx,wσjx+1,w i.e t and ∆ are real. Several interacting terms have ⊥ ⊥ w=1 j=1 j=1 been also considered in the case of one wire24–28. In this (cid:19) L work, we ignore the effect of intra-wire repulsive interac- −(∆+t)σy σy +g(cid:88)σz σz , (3) tions and assume that the Cooper channel dominates in 2 j,w j+1,w j,1 j,2 j=1 each wire. This gives a minimal model, which interpo- lates between Hubbard and Kitaev physics, and displays where w is a chain index and σx,y,z the Pauli ma- a competition between topological superconducting or- trices. This model will have the same phase dia- dering and Mott ordering. Reaching the large g limit gram as its fermionic counterpart, but different physi- could be eventually achieved experimentally by placing cal properties37. This representation favors another con- an insulating material between the wires, forming a ca- trolled experimental realization with cold atoms38,39 or pacitance between the two parallel wires. Coupling with usingJosephsonjunctionsaspseudotwo-levelsystems40, a bath could also allow to engineer such an interaction allowing to access the large g limit. Other Majorana- term (more complex interaction terms have been envi- Josephson models, similar to Ising models in transverse sionedrecentlyinRefs. 8and9). Notthattwosidegates fieldshavebeenproposed, seeforexampleRef. 41. Such couldbeusedtoscreenouttheinteractionsalongthetwo systemswouldallowustoreachthehighcouplinglimits, wires. Other interaction effects through charging energy andconsequentlytoprobeeasilythemoreexoticfeatures terms could also produce topological Kondo boxes29–33. of our system. Quantum criticality in a Ising chain has We also note that since a Kitaev superconducting wire also been observed in real materials42. can be engineered inultra-cold atomsthrough proximity The main result of this paper is the phase diagram of effect34 orviaaFloquettypeapproach35,thenacontrol- our model for ∆ (cid:54)= 0 in Fig. 1. We observe the survival 3 of the SPT phase, called 4MF, in the presence of finite II. LIMITING CASES interactions. ThisphaseischaracterizedbytwofreeMa- jorana fermions at each extremity of the ladder. Despite We discuss in this Section two limiting cases of our their proximity, the absence of an appropriate pairing system, the Kitaev model in Sec.IIA and the Hubbard or hopping term between the two ladders prevents a di- modelinSec.IIB.Then, wepresentapedagogicalmean- rectcouplingbetweentheMajoranafermionsateachex- field argument for the existence of the DCI in Sec. IIC. tremity. Weconsequentlyobserveafour-folddegeneracy of the ground state of the system with open boundary conditions. Each of these ground states has a different A. Kitaev’s wire combinationoffermionicparities. Atverylargecoupling and weak chemical potential, two similar phases appear. 1. Topological phases and Majorana fermions Both of them are Mott-Ising phases related to the Mott phases of the Hubbard model. For positive g (MI-AF), In the absence of interactions U = 0, our system re- the corresponding low-energy model is an Antiferromag- ducestotwouncoupledKitaevwires2. Thismodelisone netic Ising model, which presents orbital currents and ofthesimplestmodelspresentingatopologicalphase. In a spontaneous breaking of the time reversal symmetry. this Section, we consider a single chain and recall some For g negative (MI-F), the low-energy model is a Fer- of the main results for the sake of comparison and to romagnetic Ising model in a transverse field, which also introduce notations that will be used in the rest of the breaks time reversal symmetry and exhibits currents be- paper. tweenthetwowires. Atlargechemicalpotential,apolar- The Kitaev Hamiltonian for the fermionic species c, izedtrivialphaseopens,correspondingtodepletedorfull with open boundaries conditions (OBC), is written as: wires. At finite positive coupling, an intermediate phase opens between 4MF and the polarized phase. This is L L−1 (cid:88) (cid:88) the only gapless phase in this diagram and has a central H {c}=−µ c†c + −t(c†c +c† c ) K j j j, j+1 j+1 j charge c = 1. This phase is an extension of the critical j=1 i=1 point at g =0, whose critical model is two Ising models. +∆c†c† +∆∗c c . (4) WewilldenoteitDoubleCriticalIsing(DCI). Itisnoth- j j+1 j+1 j ing but a Luttinger Liquid (LL) of a complex mixture of The gauge invariance of the fermions allows for a can- fermions on the two wires. cellation of the phase of ∆: let ∆ = |∆|eiφ, the appli- Finally, it should be noted that other coupling forms cation cj → eiφ2cj transforms ∆ in |∆|. Consequently, in the two-chain model of spinless fermions give distinct in the rest of this paper, we consider ∆ real and posi- phase diagrams43–46. tive. Nonetheless, itcomesatthepriceofthelossofthis The outline of the paper is as follows. In Section II, gauge freedom, equivalent to the conservation of charge. wewillpresentsomewell-knownresultsonbothKitaev’s WeintroducetheMajoranaoperatorsαj bysplittingthe wire and the Hubbard model, two solvable limiting cases on-site fermions in their real and imaginary parts: of our system, and introduce the notations that will be 1 c = (α +iα ). used in the rest of the paper. We also present the re- j 2 2j+1 2j sults of some standard mean-field computations. It pro- In order for the fermions to respect their standard alge- vides some insight on the opening of the DCI phase at bra, the Majorana fermions verify the well-known Clif- a critical interaction strength. The next Sections will be ford algebra: devoted to the properties of the phase diagram and the DCI phase. In Section III, we detail an entanglement {α ,α }=2δ . j k j,k probe, the bipartite charge fluctuations47, that we use Placing ourselves at t = ∆ to keep expressions simple, forcharacterizingtheDCI phase. Afterremindingsome and focus on µ≤0, we obtain: previously obtained results on bipartite fluctuations in Luttinger Liquids48, we extend them to a more general iµ(cid:88) (cid:88) H =− α α −it α α (5) caseandcomputeitsbehaviorforthecriticalIsingmodel. K 2 2j+1 2j 2j+1 2j+2 Using exact and numerical computations and Conformal j j Field Theory (CFT) arguments, we find that the critical Atthepointµ=−2t,itreducestoasimplechainoffree Ising phase is characterized by a negative subdominant Majorana fermions of size 2L, H = −it(cid:80)2L α α , K j=1 j j+1 logarithmiccontributiontothebipartitefluctuations. In where L is the original number of fermionic sites. This SectionIV,weusebosonization49 toderivethephasedi- model is conformally invariant with a central charge c= agram close to half-filling, and explicit the properties of 150. 2 the 4FM and the two MI phases. In Section V, we de- We pose δµ = µ + 2t. The Kitaev’s model can be velop effective models to describe the emergence of the rewritten in terms of the α operators. DCI phase, and present numerical confirmations of its iδµ(cid:88) (cid:88) existence. WefinallyconcludeinSec.VI.Appendicesare H =− α α −it α α . (6) K 2 2j+1 2j j j+1 devoted to the details of computations. j j 4 Thus δµ favors the pairing of neighboring Majorana interactiontermsareinvariantunderthissymmetry(and fermions every two sites and Majorana dimerization in fermionic parity), the topological phase is preserved. the ground state. The sign of δµ differentiates two pairing,translatedbyonesite,correspondingtothestandard and topological superconducting phases. When δµ < 0, 3. Link with quantum Ising model thechemicalpotentialdominatesandcouplesMajoranas on the same site. The fermionic excitations in this case Let σx,y,z be the Pauli matrices. Kitaev wire can di- correspond to the physical fermions. In contrast to that, rectly be mapped on a spin model through the canonical when δµ>0, the hopping dominates and couples Majo- Jordan-Wigner transform36: ranasonneighbouringsites.ThereisaMajoranadangling uncoupled at each extremity of the wire (α2 and α2L+1 σz =c†c − 1 σx = 1(c +c†)exp(iπ(cid:88)j−1c†c ). for the limiting case µ=0). The phase transition occurs j j j 2 j 2 j j l l at µ=±2t. l=1 Inthe topologicalphase, theground stateis two-foldde- The obtained Hamiltonian can be written as: generate due to these two zero modes. No local operator can distinguish the two degenerate ground states of L(cid:88)−1(∆−t) (∆+t) (cid:88)L H = σxσx − σyσy − µσz. thetopologicalphase,whiletheyhaveoppositefermionic I 2 j j+1 2 j j+1 j parity. The corresponding operator is defined as: j=1 j=1 (9) L L L At ∆ = t, we recover the well-known Quantum Ising (cid:88) (cid:89) (cid:89) P =exp(iπ c†c )= (1−2c†c )= iα α . Model in a transverse field. Properties do not change j j j j 2j 2j+1 j=1 j=1 j=1 when we leave this special point, as long as ∆,t(cid:54)=0. (7) The two ground states of the topological phase have op- positefermionicparityduetothefollowingcommutation B. The Hubbard model rules: When ∆ = 0, our model reduces to the celebrated {αA,P}=0 [P,H]=0. (8) j and well-studied fermionic Hubbard model through the mapping: These results can be also obtained by solving the Hamil- tonianwithperiodicboundaryconditions(PBC),usinga c →c c →c . 1 ↑ 2 ↓ Bogoliubovtransform. DetailscanbefoundinAppendix A. The U(1) symmetry is restored and the Bethe ansatz method is applicable to solve the model at arbitrary chemical potential. In this Section, we present some 2. Symmetries known results on this stapled model52. Its exact phase diagram is displayed in Figure 2. The Kitaev’s model does not present any continuous The Hubbard model also has a SU(2) spin symmetry. symmetry. Indeed, the unusual pairing term breaks the Introducing ∆ (cid:54)= 0 breaks this spin rotation symmetry conservation of the charge (number of fermions) down down to the SO(2)∼U(1) rotations around the y-axis: to the conservation of fermionic parity. From Equation (cid:18) (cid:19) (cid:18) (cid:19)(cid:18) (cid:19) c cos(φ) sin(φ) c (A1), one can easily see that, in the case of PBC, the 1 → 1 (10) c −sin(φ) cos(φ) c ground state is odd (P = −1) in the topological phase 2 2 and even (P = 1) in the two trivial phases51. With which leave the model invariant. The associated con- OBC,thetwogroundstateshavedifferentfermionicpar- served charge is ities. With two non-interacting wires, denoting P and 1 P2 the fermionic parities in each wire, similar arguments J =i(cid:88)(cid:16)c† c −c† c (cid:17), (11) provethatthegroundstateisodd-oddinthetopological y j,1 j,2 j,2 j,1 j phase(P =P =−1)andeven-eveninthetrivialphase 1 2 with PBC (P1 =P2 =1). Similarly with OBC, the four corresponding to the total spin in the y direction. Non- ground states correspond to (P1 =±1, P2 =±1). zerovaluesofJy indicateabreakingoftime-reversalsym- There are two other discrete symmetries. The first one metry and the presence of a current between the two is the particle hole symmetry c →(−1)jc† occurring at wires. A phase with a unique ground state must conse- j j µ=0. This implies the symmetry of the phase diagram quently have (cid:104)J (cid:105)=0, and the symmetry is actually not y for µ → −µ, and half-filling at µ = 0. The second and broken in any of our phases, in agreement with Mermin- mostimportantoneisthetime-reversalanti-unitarysym- Wagner theorem for continuous symmetries. Nonethe- metry, where K is complex conjugation. Spinless real- less, as shown in Sec. VB, the competition between space fermions are left invariant by T with our choice of superconductivity and Mott physics in our model will gauge. Asmentionedintheintroduction,aslongaslocal produce orbital currents. 5 solutionsbyminimizingthetotalenergyaftersolvingthe following two consistency equations: 6 Charge M-H π g/t 4 ρ − 1 =−(cid:90) dk −(µ−g(ρ2/1− 21))−2tcosk . ngth 2 LL LL 1/2 2 0 4π(cid:113)(µ−g(ρ2/1− 12)+2tcosk)2+4∆2sin2k e n str 0 o Spin M-H Thereisnosimpleanalyticalexpressiontothesesolutions acti 2 but a numerical study reveals the appearance of a set of nter asymmetricalsolutionsatfiniteinteractionstrength. We I 4 Polarized Polarized findawholeparameterspacewherethereexistsanasym- metricalsolutionwhoseenergyislowerthanthesymmet- 6 rical solution. It is an indicator of the opening of a new 4 3 2 1 0 1 2 3 4 Chemical potential µ/t phase, androughlycorrespondstothelimitsofthe DCI phase. Nonetheless, as expected, the mean-field argument does not correctly describe its properties, obtained FIG. 2. Exact phase diagram of the Hubbard model at zero with numerical simulations. While the mean-field com- magnetic field, obtained from Bethe Ansatz52. Charge M-H putation predicts a finite difference in densities in each corresponds to the celebrated Mott-Heisenberg phase, open- wireeveninthethermodynamiclimit,numericalsimula- ing at arbitrarily low g > 0. The charge mode is gapped tions assert that the difference in electronic populations and the electronic density is fixed at half-filling, while the between the two wires is only around 2 fermions, what- spinmodeisfree(anditseffectivemodelisaSU(2)-invariant ever the number of sites we consider. Moreover, while Heisenberg model). Spin M-H is its equivalent for g < 0, the numerical simulations predict a gapless phase, here inversing the role of the two sectors. Both of them are char- thephaseisnecessarilygapped. However,themean-field acterized by a central charge c = 1. Luttinger liquid (LL0 phases) corresponds to two phases of free diluted electrons approach has the advantage to simply explain the spon- (with a central charge c = 2). Polarized phases are trivial taneous breaking of symmetry between the two wires we phases with totally empty or full wires. observeinnumericalsimulations: insteadofhavingasin- gle ground state, we obtain a doubly degenerate ground state with fermionic parity (even, odd) and (odd, even). C. Mean-field approach to the DCI phase The complete study of the DCI phase and its properties can be found in Section V. While in one dimension, mean-field computations are usually not reliable, they can give us some insights on the physical properties of a model. In the absence of Coulomb interaction g = 0, the transition between the III. BIPARTITE CHARGE FLUCTUATIONS topological and trivial phase is simply given by µ=±2t AND CENTRAL CHARGE and independent of ∆. Due to the conservation of the fermionic parity in each wire, the only partitioning of The DCI phase that appears at the transition be- the interaction that one can introduce without explicitly tweenthepolarizedandthe4MFphases(whichareboth breaking any symmetry is: gapped) is a gapless phase with a central charge c = 1. 1 1 We propose that this phase is an extension of two crit- g(ρ (j)− )n +g(ρ (j)− )n , 1 2 j,2 2 2 j,1 ical Quantum Ising points. In terms of central charge, there is no difference between a critical c = 1 bosonic where ρ1/2 is the fermionic density in each wire. In this field and two critical c = 1 Majorana fields. Indeed, it section, we will focus on ersatz where the density is con- 2 is trivial to show that one can divide the bosonic field stant in each wire. Assuming a symmetry between the intotwoMajoranasor, onthecontrary, combinethetwo two wires, we obtain a simple equation for the transition Majoranafieldstoformabosonicfield. Followingonthe lines: works in Refs 47, 48, and 53 , we present in this Section µ∓2t a criterion on the bipartite charge fluctuations to differ- g = (12) ± ρ − 1 entiate between these two cases, after some reminder on ± 2 the properties of critical phases in one dimension. Bi- where ρ is the density at the transition point µ = ±2t partite fluctuations can be measured in ultra-cold atoms ± forthenon-interactingHamiltoniancomputedinEq. A3. through the quantum gas microscope54,55, and in real materials from density-density correlation functions47 or Butonecanassumethebreakingoftheexpectedsym- though capacitance measurements56. In spin analogues, metry between the two wires c ↔ c and allow for dif- thebipartitefluctuationscanbemeasuredthroughspin- 1 2 ferent densities. One has now to compare the potential spin correlation functions. 6 A. Central charge and bipartite charge fluctuations α is the short distance cut-off. In a similar one channel model of gapped fermions, we will expect F (l) to scale A Inonedimension,conformalinvarianceofcriticalmod- as al + b + O(1), where a and b will be non-universal elshasledtoanumberofprogress,intermsofexpressing numbers. properties of physical effects in solvable models. One of Aproofcanbegiven48 : letφbeacriticalrealbosonic thefundamentalobjectofaconformallyinvarianttheory field. We recall from bosonization (more details will be is its central charge c. It is related, for example, to the giveninSectionIV)thatthechargedensityofafermionic scaling of the ground state energy for finite systems or fieldinonedimensioncanbeexpressedas−∂xφ. Onecan π the scaling of the entanglement entropy. We will use the then rewrite FA in function of the φ field: lattertocomputethiscentralcharge,inordertodiscrim- (cid:90)(cid:90) inate between possible critical models, as the entropy is π2F = dxdy(cid:104)∂ φ(x)∂ φ(y)(cid:105) ≈(cid:10)(φ(l)−φ(0))2(cid:11) A x x c c easily computed in DMRG. [0,l]2 (17) We remind the reader the definition of the entangle- To evaluate these correlators, one can either consider di- ment entropy for a non-degenerate ground state. Let |φ(cid:105) rectly the OPE of the primary fields ∂ φ or the OPE be the ground state wave function and ρ = |φ(cid:105)(cid:104)φ| the z of the vertex operators eiφ, or just the correlator of free associated density matrix. Let separate our system into bosons in a conformal theory (we note z = τ +ix and two connex parts A and B, let A be of length l. We de- ω = τ(cid:48) +iy, with τ and τ(cid:48) imaginary time considered finethereduceddensitymatrixρ ofAbyρ =Tr (ρ), A A B equal here): where Tr is the partial trace on B. We define the en- B tanglement entropy S of A by: A 1 1 ∂ φ(z)∂ φ(ω)=− S =−Tr(ρ log(ρ )). (13) z ω 4πε(z−ω)2 A A A 1 For a periodic system , the entropy scales as57: (cid:104)φ(x)φ(0)(cid:105)≈− log(|x|). 2πε (cid:18) (cid:19) c L lπ SA(l)= 3log π sin(L) +O(1). (14) Getting rid of unphysical terms due to the natural divergences of the theory, we find a term proportional to In this expression, l represents the size of the sub-region logl. To identify the coefficient ε, only free parameter A while tracing out the region B (the total system has of the critical theory, one can look at physical observ- the length L). Gapped phases are not conformally in- ablessuchasthecompressibilitytoobtainEquation1648. variant and, as a result, have a ”central charge” c = 0, while gapless phases and critical points have usually non This result can be extended to a wire with several trivial central charge. Of relevance in this paper are the channels and a central charge c = m, m ∈ N, corre- central charges of Quantum Ising or a free wire of Majo- sponding to m gapless bosonic modes and possibly other ranafermions,c= 1,andofaLLorafreescalarbosonic gapped modes described by Sine-Gordon models. For 2 mode, c = 1. When there exists several independent bipartite charge fluctuations that are quadratic in the modes, this central charge is additive. Finally, this cen- bosonic fields, it is shown in Appendix B that the log- tral charge can be challenging to compute in finite size arithmic term in bipartite charge fluctuations is always system, if one has no a priori knowledge of the critical positive. theoryandinparticularonthemappingbetweencritical theoryandoriginalmodel. Indeed, itcanbenumerically challengingtoextracttheexactlogarithmiccontribution B. Bipartite charge fluctuation for a c= 1 model 2 withouttakingintoaccounttheprecisefinite-sizecontri- butions. Wearenowinterestedincomputingthechargefluctu- Another entanglement observable we will be interested ations in the case of a c= 1 model. We refer the reader in is the scaling of the bipartite charge fluctuations of A. to Ref. 58 for a review on c2onformal field theory. Let QA be the charge operator of the subsystem A, then Before entering into the details of the bipartite fluctu- F =Tr(Q2ρ )−Tr(Q ρ )2 =(cid:10)Q2(cid:11)−(cid:104)Q (cid:105)2. (15) ations, there is a first point to address. In the Kitaev A A A A A A A model (and its counterpart Quantum Ising), the total We present a summary of some results detailed in 47. charge (total spin along z) is no longer a good quantum Differences of scaling in FA allow for differentiating be- number and a conserved quantity. In particular, the tween gapped and gapless phases. In particular, in the groundstatehasnolongeraproperelectronnumberbut case of a one-channel free fermion model with U(1) sym- is a superposition of states with a different number of metry (a c=1 model), the bipartite charge fluctuations particles. This has strong consequences on the bipartite scale logarithmically with the size of A and can be com- charge fluctuations: there are no longer any symmetry puted from the underlining conformal theory. With K between F (l) and F (L−l). The quickest way to show A A the LL parameter of the model, one obtains: this is to consider the two limiting cases F (0) and A K l FA(L). While FA(0) is still 0, FA(L) is no longer equal FA(l)= π2 logα +O(1). (16) to 0 as the charge is no longer fixed on the whole wire. 7 As a consequence, the standard trick from CFT to take Posing v = 2tα and m = αδµ, the Hamiltonian of the into account finite size effects, replacing the length l of system is now given by: thesegmentweconsiderby Lsin(lπ)isonlyvalidforthe π L (cid:90) iv direct conformal contributions. Consequently, we start H = dx (γ ∂ γ −γ ∂ γ )−imγ γ . (21) by considering an infinite wire. Moreover, we can expect K 2 L x L R x R L R a non-zero linear contribution to F , and not just a A Atm=0,onecanidentifythisHamiltonianwithitscon- logarithmic one. This change breaks the symmetry that formalactioncounterpart. Introducingz =τ+ix,γ and L existed between the bipartite charge fluctuations and γ correspondstotheholomorphicandanti-holomorphic R theentanglemententropyforstandardLuttingerLiquids. part of the conformal field. The action is given by: (cid:90) We will prove in this Section that the sign of the S =ε dzdz¯γ ∂ γ +γ ∂ γ . (22) logarithmic contribution to the charge fluctuations L z¯ L R z R is negative. The change in the behavior of these fluctuations comes both from the different underlying ε is the critical energy scale we will check with corre- critical theory, but also from the difference on how lation functions. to express the fermionic density in terms of the pri- The fermionic density operator (i.e the σz field for mary fields. Indeed, the difference will subsist in the the spins) can be written in terms of primary fields as case of a c = 2 × 1 theory. We propose two differ- iγL(z)γR(z¯), up to constant terms that will disappear entproofsforthecom2putationofthechargefluctuations. because we are interested in the connected correlators. One can consequently rewrite: The first one is simply based on exact computation on (cid:90)(cid:90) the regularized lattice model. Details can be found in F (l)=− dxdy(cid:104)γ (ix)γ (−ix)γ (iy)γ (−iy)(cid:105) A L R L R c AppendixC.Onecanrecovertheexactexpressionofthe [0,l]2 linear coefficient for all ∆/t, using Fejér theorem: (23) Using the Operator Product Expansion (OPE)58,59 |∆] F (l)= l+O(logl). (18) γ (z)γ (w) = 1 1 for the Majorana field directly A 2|∆|+2t L L 2πεz−w yields the result: Toobtainthesub-dominantlogarithmcoefficient,amore (cid:90)(cid:90) 1 1 involved computation is required. For ∆ = t, one can F (l)=− dxdy A 4π2ε2|x−y|2 proceed to the complete computation and obtain: [0,l]2 l 1 γ +2log(2) 1 F (l)= − log(l)− euler +O(1). (19) ≈− logl+αl+β. (24) A 4 2π2 2π2 2π2ε2 The logarithmic contribution for the bipartite charge The minus sign in front of the logarithmic contribution fluctuations are this time negative. Numerical compu- can be understood from the (i)2 = −1 pre-factor in Eq. tations of the relevant integrals confirm that the results (23) stemming from the definition of the (electron) den- stand for all ∆/t(cid:54)=0. sity operator as iγL(z)γL(w). As in the computation of Wecanalsorecoverthisresultdirectlyfromtheunder- Section IIIA, we directly got rid of the unphysical di- lying conformal theory. The critical conformal theory of vergences of the theory, that are due to the absence of QuantumIsingcanbeexpressedasatheoryofafreereal ultraviolet cut-off. α and β are a priori non-universal (Majorana) fermions ψ, as shown in Section II. We only constants that arise from the integration, and are linked consider the point t = ∆ for simplicity, but the analysis to the cut-off of our theory. The results of our two com- stands at all ∆(cid:54)=0. putations coincide for ε = 1. This can be confirmed by We can reformulate the Hamiltonian in the following the computation of the correlation function correspond- way. Let γj,b = α2j and γj,A = α2j+1 the two different ing to iγL(ix)γL(iy) in the original Bogoliubov particles species of Majorana in the wire. language. One can show that the coefficient of the lead- ingtermcorrespondingtotheOPEexpansionoftheCFT H =it(cid:88)γ (γ −γ )−iδµ(cid:88)γ γ . (20) corresponds indeed to ε = 1 and, moreover, that it does K j,B j+1,A j,A 2 j,A j,B not depend on ∆/t. In Kitaev’s wire, the logarithmic j j contributions to the bipartite charge fluctuations are ac- √ We first go to continuous limit γ → 2αγ (x), tually independent from the ratio ∆/t as long as ∆(cid:54)=0. j,A/B A/B where α is the lattice spacing. Then we introduce two Figure3providesacomparisonbetweentheexactresults chiral fermions: and some numerical computations using DMRG with Matrix Product States (MPS) from the ALPS project γ −γ γ = B√ A code60,61. While the agreement is nearly perfect for the R 2 linear contribution, the logarithmic contribution is more γ +γ γ = A√ B. complextocatchandslowtoconverge. Inparticular,for L 2 ∆/t≤0.2, finite-size effects are far from negligible. 8 Finally,Figure4presentstheevolutionofthetwocoef- ficientsaswecrossthetransitionlineforoneKitaevwire 0.40 at ∆=t. The introduction of the mass term m in Eq.21 formally cuts the logarithm at long range. The logarith- 0.35 mic contribution then goes to zero smoothly away from A quantum criticality. The change in the linear coefficient Fof 0.30 also marks the transition. The exact expression can be nt 0.25 found in Appendix C, obtained from the lattice model e ci using Fejér theorem. It is constant for |µ|≤2t. effi0.20 o c r 0.15 a e IV. CLOSE TO HALF-FILLING Lin0.10 Theoretical computation 0.05 DMRG In this Section, we are interested in the physics of the Integrals on infinite wire complete model of the two wires close to the line µ=0, 0.00 0.0 0.5 1.0 1.5 2.0 2.5 3.0 corresponding to half-filing. We present an analytical Pairing strength ∆/t approach to the problem, assisted by numerical verifica- (a)Linearcoefficientofthechargefluctuations. tions using both Exact Diagonalization (ED) and Den- sity Matrix Renormalization Group (DMRG). Compu- tations to check convergence of the different parameters 0.00 have been pushed up to 150 fermions per wire (a total of Theoretical computation 300 fermions). DMRG FA 0.01 Integrals on infinite wire of t en 0.02 A. Bosonization at half-filling ci effi o 0.03 c We will proceed with a standard Abelian bosonization c scheme49, considering both ∆ and g as perturbations. mi h 0.04 The fundamental property behind bosonization is the rit a correspondence between the excitations of a fermionic g system above the Fermi surface and the excitations of a Lo 0.05 well-chosen real bosonic field. 0.06 First we will assume that we are at half-filling and 0.0 0.5 1.0 1.5 2.0 2.5 3.0 Pairing strength ∆/t consequently far from the bottom of the energy band in each wire. We can then linearize the spectrum around (b)Logarithmiccoefficientofthechargefluctuations. the Fermi energy and separate each fermion field into its left- and right-moving part: FIG. 3. (Color online) Linear and logarithmic contributions to the bipartite charge fluctuations as a function of ∆/t for c =eikFjc +e−ikFjc . (25) j,σ 1,j,σ −1,j,σ a single Kitaev wire at the transition point, from theoreti- calcomputationsandDMRGsimulations60,61 witha90-sites Here, σ is the index of wire (σ = 1,2) and r = 1(−1) wire. Linear contributions are easily extracted from the sim- representstheright-(left-)movingelectrons. Goingtothe ulations and match very well the theoretical values. On the continuum limit and linearization of the hopping term other hand, logarithmic contributions are more complex to gives: extract and convergence is slower. We observe yet a good agreement for ∆/t>0.3. Below this value, finite-size effects (cid:88)(cid:90) (cid:16) (cid:17) cannotbeneglected. Fithavebeenrealizedusingtheinfinite −ivF dx c†1,σ(x)∂xc1,σ(x)−c†−1,σ(x)∂xc−1,σ(x) wire form of the fluctuations, as, at long range, there is a σ precision loss due to the dominant linear term. (26) v = 2tαsin(k ) is the Fermi velocity, α is a short dis- F F tancecut-offweintroducetogofromalatticetheorytoa continuous theory (Conventionally, α corresponds to the ify the following commutation rules: distance between two sites of the lattice) and k = π/2 F theFermimomentum. Wethenintroducetworealscalar π [φ(x),θ(x(cid:48))]=i sgn(x(cid:48)−x) and {U ,U†}=2δ . bosonicfieldsφσ andtheirconjugatefieldθσ correspond- 2 r r(cid:48) r,r(cid:48) ing to the excitations of each wire, and four Majorana (27) fermions U that will serve as Klein factors. They ver- The well-known mapping between the fermions and the r,σ 9 Term Dimension Bare value (cid:113) 0.26 0.00 vF,± - vF 1± πvgF K - K (0)= 1 0.24 FA √ ± √ ± (cid:113)1±πvgF FA0.22 0.01of cos( 2θ−)√cos( 2θ+) 12(K+−1+K−−1) ∆(1)(0)= 4∆πα r coefficient of 000...121608 00..0032mic coefficient ccccoooossss((((2222√√√2222φφθθ−+−+)))) 2222KKKK+−−−−+11 gg∆∆−+(+(−((2200))))((00==))==22−ππgg0022 a h e t Lin0.14 0.04gari TABLEI.Dimensionsofthedifferenttermsofthebosonized Linear coeff (theory) Lo model, and bare values in the RG flow. 0.12 Linear coeff (DMRG) Log coeff (DMRG) 0.10 0.05 3.0 2.5 2.0 1.5 1.0 theory, we do not consider the renormalization of the Chemical potential µ/t Fermi velocities. K and K are the Luttinger parame- + − ters. g andg appearduetotheCoulombcouplingbe- + − tweenthetwowires. BotharepresentinHubbardmodel FIG. 4. (Color online) Linear and logarithmic contributions andareresponsiblefortheMott-Heisenbergphases. The to the bipartite charge fluctuations as a function of µ/t for a singleKitaevwireacrossthetransitionpoint,at∆=t,from pairing ∆ plays now the role of a coupling ∆(1) between theoretical computations and DMRG simulations with a 90- the two charge and spin sectors, that cannot be a priori site wire. Linear contributions are easily extracted from the separated. ∆(2) and ∆(2) are not initially present in the + − simulations and match very well the theoretical values. This bareHamiltonian, but aregeneratedunderRGby ∆(1) . coefficientbecomesconstantatµ=−2t. Logarithmiccontri- In a diagrammatic language, they correspond to second butionsincreaseaswereachthetransition,andaremaximum order contributions in ∆. at this point. We define the renormalization length as: α(l) = αel. Bosonizationconsequentlyallowsustorecoverthefollow- bosonic fields is given, in the thermodynamic limit, by: ingrenormalizationflowequations,includingallrelevant orders: U cr,j,σ = √2rπ,σαe−i(rφj,σ−θj,σ). (28) dK± =−2π2g±2 K2 + 2π2(∆(±2))2 + π2(∆(1))2 dl v2 ± v2 4v v F,± F,± F,+ F,− After some computations, using the convention U†U† = dg L R ± =(2−2K )g −i, the Hamiltonian at g =0 is: dl ± ± H =(cid:88)(cid:90) dxv2Fπ((∂xθσ)2+(∂xφσ)2)+2∆sπinα(kF)cos(2θσ) d∆dl(1) =(2− 12(K+−1+K−−1))∆(1) σ d∆(2) 2π2(∆(1))2 (29) ± =(2−2K−1)∆(2)+ . (31) with vF = 2tasin(kF). We then introduce the charge dl ± ± vF,∓ φ+ = φ1√+φ2 and spin ( relative charge mode in our case) Dimensional analysis of these equations discriminates 2 modesφ− = φ1√−φ2. Thedensity-densityinteractionterm threedifferentphasesathalf-filling,tobeanalyzedbelow. renormalizes the2Fermi velocities and the Luttinger pa- To qualitatively compare the effects of gε and ∆(1), we rameters of both bands, and induces two cosine terms. compare the bare value of the latter, an a priori strongly The final Hamiltonian, including terms generated in the relevantcoupling, andtheeffectivemassmg obtainedby Renormalization Group (RG) process, is given by: Bethe-Ansatz in the Hubbard model (in other words the Mott gap in the charge sector), at low coupling52 H = (cid:88)(cid:90) dxv2Fπ,ε(Kε(∂xθε)2+Kε−1(∂xφε)2) mg =−2+ |g| +2(cid:90)∞dωJ1(ω)e−|g4|tω ε=± t 2t ω g √ ∆(2) √ 0 + ε cos(2 2φ )− ε cos(2 2θ ) (cid:114) α2 ε α2 ε ≈ 4 |g|e−2|πg|t for |g|(cid:28)t. ∆(1) √ √ π t + cos( 2θ )cos( 2θ ). (30) α2 + − J is the Bessel function of the first kind. 1 √ Scaling dimensions and bare values of the different pa- • If ∆ (cid:38) |m |, ∆(1) dominates the g cos(2 2φ ) g ± ± rameters can be found in Table I. At this order of RG terms and goes to strong coupling. Both θ modes ± 10 become massive and are locked to the minima of B. Characterization of the 4 Majorana phase the ∆(1) term. By continuity with the topological phaseatg =0,weexpectthisstrongcouplingfixed Based on adiabaticity, we expect the topological prop- point to correspond to the SPT phase presenting erties of the 4MF phase to be well-described by the four Majoranas, or 4MF phase. g = 0, |µ| < 2t case, i.e. two uncoupled Kitaev wires in their topologically non-trivial phases. Hence, four • If |m | (cid:38) ∆ and g > 0, g is renormalized to g + zero-energy Majorana end states should be present and large coupling before ∆(1) reaches significant val- remain uncoupled, corresponding to a fourfold degen- ues. g is irrelevant and renormalized to 0. φ − √ + erate ground state. We present in this section a few is consequently locked to 0 [π/ 2], corresponding analytical and numerical arguments that support this to a Mott ordered phase, and ∆(1) vanishes at this claim. Following the Z classification by Kitaev and 8 fixed point. ∆(2) is still relevant and acquires a Fidkowski22, the phase is, in this case, simply a sym- − non-zero value in the initial steps of the renormal- metry protected topological phase (SPT). Indeed, these ization. It consequently gaps the spin sector and Majoranas are not protected against any local interac- both modes (φ+,θ−) are eventually locked. This tion: an arbitrarily low pairing term between the two fixed point describes the MI−AF phase. wires i|∆ |(c† c† −c c ) would directly couple the ⊥ j,1 j,2 j,2 j,1 free Majoranas, destroying the phase. • If |m | (cid:38) ∆ and g < 0, the reasoning is the same g A first approach consists in considering the pertur- as for g >0 but for an inversion of the charge and bative effect of g on the extremity of two Kitaev wires the spin sector. (φ ,θ ) are locked, describing the − + in the topological phase. We assume t = ∆ to get a MI−F phase. simpler picture. Using the notation of Section IIA, we recall that the free Majorana fermions are γB and γA As long as we stay close to half-filling, one can use the 1,σ L,σ (the additional σ is the wire index). The interaction same bosonization scheme to determine the effects of a term Eq. (2) can be rewritten as: −g (cid:80)γA γB γA γB . chemicalpotential. Indeed,oneonlyneedstoaddaterm: 4 j,1 j,1 j,2 j,2 j √ Only terms at least of order (g)L will directly couple µ 2 t − ∂ φ , the free Majoranas. In the thermodynamic limit, we π x + consequently expect the survival of the 4 Majoranas phase. intheHamiltonian. Theeffectofthistermistwo-fold: it √ reducestheeffectivedimensionofcos(2 2φ )andrenor- + From a numerical point of view, we have studied sev- malizes the Fermi velocity62. When the renormalized eral markers for the topological phase. The first is ob- Fermi speed approaches 0, it indicates that we are too viously the change in degeneracy going OBC to PBC. far from half-filling and that the spectrum is no longer While in the latter the ground state present no degen- linear, leading to the breakdown of the bosonization ap- eracies, it is four-times degenerate in the former. The proximation. We will summarize the effect of the chemi- parity of the ground states follows the same rules as in cal potential on each of the previously obtained phases: the non interacting case: we go from an odd-odd ground • If ∆(cid:38)|m |, none of the fixed operators includes a state (P1 = P2 = −1) with PBC to a ground state in g each parity sector (P = ±1, P = ±1) with OBC. As φ term, meaning that no transition occurs before 1 2 + fermionicparityineachwireisagoodquantumnumber, wereachthebottomofthebandandthebosoniza- it is quite easy to observe this degeneracy with both our tion procedure breaks down. ED or DMRG simulations. Typical behavior for the en- • If|m |(cid:38)∆andg >0,theumklaaptermcontrolled ergy of the first few levels on the line µ=0 is presented g by g starts oscillating. In the Hubbard model oc- in Figure 5. + curring at ∆ = 0, g is renormalized to zero at a A second good marker for topology is the de- + finite ratio µ, corresponding to a vanishing charge generacy of the entanglement spectrum63 in the g gap and a commensurate-incommensurate transi- periodic ladder.These eigen-energies are obtained by tion to a gapless Luttinger phase62. At finite but re-interpreting the eigenvalues of the reduced density small ∆, µ weakens g+ by reducing its dimension matrix ρA = TrA(|Ψ(cid:105)(cid:104)Ψ|) as Gibbs exponential forms. until ∆(1) dominates the RG process and flows to Here A is simply the half of the wire. A comprehensive strong coupling. This leads to a resurgence of the and clear point of view of the properties of the entangle- 4 Majorana phase at finite µ. ment spectrum for topological fermionic phases can be found in Ref. 64. The presence of Majorana boundary • If |m | (cid:38) ∆ and g < 0, as with the first phase, a statestranslatesintoafour-folddegeneracyintheentan- g transitiondoesnotoccurinthebosonizationrange. glement spectrum, as cutting the system is analogous to creating new boundaries. As in the previous numerical We describe in the next two parts the properties of probe, this degeneracy can be observed in both ED these three phases. or MPS simulations. Typical behavior for the entan-