Chiral topological insulating phases from three-dimensional nodal loop semimetals

0.25 MB·

by  Linhu Li

Chiral topologicalinsulating phases from three-dimensional nodal loopsemimetals
Linhu Li,1,2 Chuanhao Yin,3 Shu Chen,3,4,5 and Miguel A. N. Arau´jo1,2,6
1Beijing Computational Science Research Center, Beijing 100089, China
2CeFEMA,InstitutoSuperiorTe´cnico,UniversidadedeLisboa,Av. RoviscoPais,1049-001Lisboa,Portugal
3Beijing National Laboratory for Condensed Matter Physics,
Institute of Physics, Chinese Academy of Sciences, Beijing 100190, China
4School ofPhysical Sciences, UniversityofChineseAcademyof Sciences, Beijing, 100049, China
7 5Collaborative Innovation Center of Quantum Matter, Beijing, China
1 6Departamento de F´ısica, Universidade de E´vora, P-7000-671, E´vora, Portugal
0
Weidentifyatopological Z indexfor three dimensional chiral insulatorswith P∗T symmetry wheretwo
2
Hamiltoniantermsdefineanodalloop. SuchsystemsmaybelongintheAIIIorDIIIsymmetryclass. TheZ

n invariant isawindingnumber assignedtothenodal loopandhasacorrespondence tothegeometricrelation
a betweenthenodalloopandthezeroesofthegapterms.Diracconeedgestatesunderopenboundaryconditions
J
areincorrespondencewiththewindingnumbersassignedtothenodalloops. Weverifyourmethodwiththe

2 low-energyeffectiveHamiltonianofathree-dimensionalmaterialoftopologicalinsulatorsintheBi Te family.
2 3
2
PACSnumbers:03.65.Vf,71.20.-b,71.10.Fd
]
l
l
a
Introduction.-Topologicalinsulators(TIs)in threedimen- this case, although the NL itself does not carry topological
h
sions (3D) having time reversal symmetry can be character- charge, the addition of some chiral gap terms can make the
-
s izedbyZ numbersdefinedonsomediscretemomenta[1–5], system become topological and insulating, where the topo-
e 2
m which is equivalent to a quantized invariant expressed as an logicalinvariantisgivenbyawindingnumberdefinedalong
integral over the entire Brillouin Zone (BZ)[6]. There are theNL[29]. Aninterestingquestiontoaddressis,whateffect
.
t thentwotypesofTI’s,strongandweak,accordingtowhether cangaptermshaveona3DNLsemimetal?
a
m there is an odd or even numberof Dirac cone surface states, Inthispaper,westudya spin-1/2systemwith P∗T sym-
respectively[5]. BesidesthetimereversalTIs,thereisaclass metry,andshowthatanticommutingmassgaptermscandrive
-
d of chiral TIs which are described by a Z-type topological a3DNL semimetalintoa chiralTI,whichcanbecharacter-
n invariant[7]. Thegeometricalrepresentationofa topological izedbyaintegerwindingnumberdefinedalongeachNL.This
o
invariant in some vector spaces provides an intuitive way to windingnumberis determinedby the geometricrelation be-
c
analyzethetopologicalnatureofmanysystems[8–12]. tweentheNL andthezeroesof thegapterms. Althoughthe
[
Ontheotherhand,atransitionpointbetweentopologically gaptermsmaybeinitiallyconsideredsmall, ourresultsonly
1
depend on their zeroes, so that it holds valid for any finite
v differentinsulatingphasescanbeviewedasasemimetalphase
termsthatgapouttheNLs. Thesystem’ssurfacestatesmay
3 with nontrivial topology in its gap closing points[13, 14].
3 Topologicalsemimetals(TSMs)[15–21]haveaFermisurface hold an odd or even number of Dirac cones, and their exis-
1 tence has a correspondence to the NL winding number. In
(FS) with reduceddimension. While a 3D normalmetalhas
6 thissense,theNLcanservesasanindicatorofthetopological
a two-dimensional (2D) FS, a TSM has a one-dimensional
0
propertiesofa3Dinsulator.Inordertoshowtheutilityofour
. (1D)orzero-dimensional(0D)FSathalf-filling. 3Dsystems
1 theory,weapplyittoa3DmaterialofthetopologicalBi Te
with0DFSareknownastheWeyl[15]orDirac[16]semimet- 2 3
0 familyandgiveabriefdiscussion.
7 als. Inthesesystems,thetwobandstouchlinearlyatdiscrete
Minimal model.- We begin our discussion with a simplest
1 gapclosingpointsintheBZ,andholdtopologicallyprotected
twobandmodelforP∗T symmetry-protectedNLsemimetals:
: edge states under open boundary conditions (OPC), such as
v
Xi fsoemr iinmsetatanlcse[2,2th–e27F]ehrmavieaartctrsa.ctMedogrerorwecinegntalytt,e3nDtionn.odInalsluicnhe H0 =(m−k2x−ky2)σx+kzσz, (1)
r systems, the linear bandtouchingpointsformone or several
a withσthePaulimatricesactingonanorbitalspace. Forspin-
1DlinesintheBZ.
lesssystem,theP∗T symmetryissimplygivenbythecom-
One of the most interesting cases is when the nodal lines plexconjugationandaunitarymatrix,P,suchthattheHamil-
form closed, nodal loops (NLs). A NL can be classified in tonian satisfies PH∗(k)P−1 = H (k). In this case P = 1 and
0 0
either of two classes, according to whether it carries a Z the P∗T symmetryensuresthe absence of the second Pauli
2
monopolechargeornot. Theonewithoutamonopolecharge matrix.ThenodesofH yielda1Dsolution,aNLwithk =0
0 z
canshrinkintoapointanddisappear,andistopologicallytriv- andk2+k2 = m. Introducingaσ term,H = H +h σ ,not
x y y 0 y y
ialinthissense. NLsareprotectedbythecombinationofin- onlybreaksP∗T symmetrybuth alsoservesasaneffective
y
versionandtimereversalsymmetries,P∗T,forspinlesssys- massterm,whichcanbeeitherk-independentorrelatedwith
tems,whileadditionalsymmetriesarerequiredtoprotectNLs k. Intheformercase, itopensa stablegapintheBZ, which
in3Dsystemswithspin-orbitalcoupling[26,27].Ontheother drivesthesystemintoatrivialinsulator. Ifh isafunctionof
y
hand,NLsemimetalshavealsobeenstudiedin2D[28,29].In k,thenodesof H maybepairsofpoints,andthesystemisa
2
Weylsemimetal.
Weextendthismodelbyincludingthespindegreeoffree-
domandwritetheHamiltonianwithP∗T symmetryas
H = H +h σ s +h σ s +h σ s , (2)
0 yx y x yy y y yz y z k
z
ν=1
wherethePaulimatricess actinspinspace. P∗T symmetry
j
now reads Ps H∗(k)s P−1 = H (k), with P = 1 andsatisfies ν=−1
y 0 y 0
(P∗T)2 =−1,Hereweusethelabelsh torepresenttheterm
ij
of σ and s , with i = 0 or j = 0 for the identity matrix in
i j
the corresponding subspace. These five terms form an anti-
commuting set of Dirac matrices, but we note that there are FIG. 1: A sketch of the winding number of a NL. The blue circle
indicatestheNLgivenbyH = 0,theredlinesindicatethe1Dso-
alsootherequivalentchoices[30]. Thespectrumof(1)-(2)is 0
lutionofh = 0for modelsof thechiralclass, andeachlinegives
simplygivenby y
a singularity in the integrand of (4) in the NL plane. The arrows
aroundeachsingularityshowthedirectionofh nearthispoint,and
y
E± =±qh2x0+h2z0+h2yx+h2yy+h2yz. (3) thecorrespondingwindingnumberofthesepointsarelabeledinthe
sketch.ThewindingnumberoftheNL,ν ,isgivenbythesumma-
NL
TheeffectivemassgappingouttheNLisnow|hy|,wherehy = tionofwindingnumbersofeverysingularitywithintheNL,asthe
(hyx,hyy,hyz). Requiringhy = 0maygiveasolutionofpoints integralpathcanbesmoothly transformedintocirclesaroundeach
(0D),lines(1D)orsurfaces(2D),dependingonthenumberof singularity.
non-zerok-dependenttermsithas.
If |h | contains two non-zero k-dependent terms, h = 0
y y
model described by an anticommuting set of Dirac matrices
shallgiveoneorseveral1Dlines.Thusthesystemisgenerally
Γ=(σ s ,σ s ,σ s ,σ s ,σ s ),as
aninsulator,asthegapclosingconditionrequiresthecrossing x 0 z 0 y x y y y z
oftheNLandthese1Dlines,whichisaccidental.Suchafour-
H = h(k)·Γ, (5)
componentDirac Hamiltonian describes a chiral topological
insulator[7], as the modelsatisfies SH(k)S−1 = −H(k), with
with
the chiral operator S given by the absent fifth Dirac matrix.
Intheabsenceoftimereversalsymmetry,thesystembelongs h(k) = (h ,h ,h ,h ,h ),
x0 z0 yx yy yz
totheAIIIclassandcanbecharacterizedbyaZinvariant[3].
h = µ−t (cosk +cosk ),
Wedefineawindingnumberofh alongtheNL[29], x0 k x y
y
h = −t cosk ,
z0 ⊥ z
h dh −h dh
ν = 2 1 1 2, (4)
NL INL |hy|2 whichformtwoNLsinkx −ky planeforkz = π/2and−π/2
when0 < |µ/t | < 2. The positionandshapeof theNLs are
k
canbeshowntobeequivalenttoaBerryphaseoftheoccupied
onlyassociatedwiththeratioofµandt ,hencewecanchoose
k
Blochbandsathalffilling[31]. Here,h andh denotethetwo
1 2 t =t =1forthesakeofsimplicity.Inthefollowingweonly
k ⊥
k-dependenttermsofh . Ifweconsiderthek −k planethat
y x y considerthe case with positiveµ, with the centerofthe NLs
containstheNL,theintersectionofthe1Dlinesandtheplane
givenbyk = k = 0. For negativeµ, the centerof the NLs
x y
producesaseriesofsingularitiesintheplane,andthewinding
is at k = k = π, and a similar discussion applies. Thetwo
x y
number(4)oftheNLissimplythesummationofthewindings
NLs of k = ±π/2 give two independent winding numbers,
z
aroundthesingularitieswithintheNL,asshowninFig.1.This
ν± , respectively,andwedefinethetotalwindingnumberof
windingnumbermay take on anyinteger value, as it is only NL
thesystemas
associatedwiththenumberoflinesgoingthroughtheNL.We
alsonotethatinthepresenceoftimereversalsymmetry,this ν =ν+ +ν− . (6)
sum NL NL
model would fall into the CII class, which is described by a
Z2topologicalindexinstead[3]. Withoutlossofgenerality,herewechoosehyx =0topreserve
Finally,ifallthethreetermsofhyaredependedonk,hy =0 achiralsymmetrywiththeoperatorS =σysx. Theothertwo
shallgiveoneorseveral0Dpoints. Inthiscase,theintersec- gaptermsofh arefunctionsofk,andh = 0gives1Dlines
y y
tionofthesepointsandtheplaneoftheNLisalsoaccidental. in the BZ. We considerthe followingform of h thatbreaks
y
Fromthesymmetryclassificationpointofview,thepresence timereversalsymmetry:
of the fifth Diracmatrix breaksthe chiralsymmetry, andthe
modelfallsintothe Aclass, whichisnon-topologicalin3D. hyy = −tysinky (7)
Inotherwords,wecouldsmoothlymoveasingularityoutof h = µ −t sink −t sink , (8)
yz s x x z z
theNLwithoutclosingthegap.
Winding numbers and geometry of the loops for a lat- and the system falls into the AIII class. |h | = 0 givessome
y
tice model.- In order to reveal the topological propertiesde- 1Dlinesinthek −k planewithk = 0orπ,whichmayor
x z y
scribedbytheNLwindingnumber,wenextconsideralattice maynotbeenclosedbytheNLs. Wewouldalsoliketopoint
3
1 (a) 1 (b) 1 (c) trajectoryofhinthe1-5planeinΓspace,
πk/z 0 πk/z 0 πk/z 0 h(k )=(−cosk ,0,0,0,−sink ). (9)
x x x
−1 −1 −1
1 1 1
0 0 0
ky/π−1 −1 kx/π0 1 ky/π−1 −1 kx/π0 1 ky/π−1 −1 kx/π0 1 TheremainingHamiltonianterms,
(d) (e) (f)
1 1 1 h0 =(µ−cosk ,−cosk ,0,−sink ,µ −t sink ), (10)
y z y s z z
πk/z 0 πk/z 0 πk/z 0
canbeviewedasthevectorfromtheoriginofthevectorspace
−1 −1 −1
1 1 1
ky/π0−1 −1 kx/π0 1 ky/π0−1 −1 kx/π0 1 ky/π0−1 −1 kx/π0 1 tuolatrhceocmenptoenreonftsthoefehl0liptosethhe(k1x-)5.pTlahneepaarrealgleivleanndbyperpendic-
FIG.2: TheNLsandthe1Dlinesgivenby|hy| =0,representedby h0 = (µ−cosk ,0,0,0,µ −t sink ), (11)
blueandredlinesrespectively.Forthelatter,weonlyshowthelines k y s z z
withky = 0,astheonesforky = πlieoutsidetheNLs. (a)µ =0.2, h0⊥ = (0,−coskz,0,−sinky,0). (12)
µ = 0.8andt = 0; (b) µ = 1, µ = 0.8andt = 0; (c)µ = 1.8,
s z s z
µs =0.8andtz =0;(d)µ =1.6,µs =0.5andtz =0.2;(e)µ =1.6, The existence of edge states depends on whether the ellipse
µs =0.5andtz=0.4;(f)µ=1.6,µs =0.5andtz=0.6. h(k )enclosesthepointh0,andthisconditionreads
x k
|h0|= (µ−cosk )2+(µ −t sink )2 <1. (13)
outthatalthoughthe gaptermsneedbesmallforthesystem k q y s z z
topreserveaNLlikestructure,thetopologicalpropertiesare
ProvidedthatEq.(13)holds,theedgestateenergiesaregiven
not related to the exact value of these terms, but only to the
by
ratios between them. For the sake of simplicity, we choose
t =t =1andpositiveµ hereafter.
y x s E =|h0|=± cos2k +sin2k . (14)
Wefirstconsiderasimplecasewitht =0. Inthiscase,the ± ⊥ q z y
z
1Dsolutionofhy =0givesfourlinesperpendiculartokx−ky CandidateDiracconesat(k ,k ) = (0,±π)mustsatisfythe
plane, two with ky = 0 andtwo with ky = π. Forpositiveµ, inequality(13). Whent =0,Eyq.(z13)becom2es
z
thepairoflineswithk = πarealwaysoutsidetheloop. By
y
tuningµsandµ,theNLsmayenclose2,1or0lineswithky = (µ−1)2 <1−µ2s, (15)
0, asshown byFig.2(a)-(c). However,the windingsof these
linesinthek −k planehaveoppositevalues(asinFig.1),and forbothk =±π. Inotherwords,itgiveseithernoDiraccone,
x y z 2
theNLenclosingeither0or2lineswillresultinν =0. On or a pair of Dirac cones at (k ,k ) = (0,±π). For nonzero
NL y z 2
theotherhand,asthesystempreservesareflectionsymmetry t , there may exist 0, 1 or 2 Dirac cones depending on the
z
alongzdirection,eachlinewillbeenclosedbyeithertwoor parameters,asthecondition(13)becomes
noneoftheNLs,hencethetotalwindingnumber,ν ,inthis
sum
caseisalwayseven. (µs∓tz)2 <2µ−µ2, (16)
Inthe presenceofa nonzerotz, thereflectionsymmetryis fork =±π,respectively.InFig.2wedisplayphasediagrams
broken, and the lines of |h | = 0 will change shape with t z 2
y z showingthenumberofDiracconeedgestates. Theseresults
andeventuallyformaclosedring,asshowninFig.2(d)-(f).In
arealsoinconsistencewiththewindingnumbersoftheNLs,
thiscase, anenclosedlinewillcrossoneoftheNLsatsome
asshowninthefigure.
point, resulting in a topological phase transition. After this
In orderto visualize the edge states, nextwe choose open
transition, thesystem hasanoddwindingnumber,ν = 1,
sum boundaryconditioninxdirectionandrewritetheHamiltonian
asonlyoneoftheNLsenclosesasingularity.
asatight-bindingbetweenplanes
edgestatesandphasediagram.-Thetopologicalproperties
ofa3Dtopologicalinsulatorcanberepresentedbythenum- H = cˆ†U cˆ + 1cˆ† U cˆ +h.c., (17)
ber of Dirac cones in the edge states under OBC. Next, we X n d n X2 n+1 od n
n n
apply the method in Ref.[32] to study the edge states in our
model. The existenceof edgestates and their eigenenergies, withcˆn isavectorofannihilationoperatorscˆn,σ,s onplanen,
underOBC,areassociatedwiththebulktopologyofthesys- cˆn =(cˆn,+,+,cˆn,+,−,cˆn,−,+,cˆn,−,−)n,and
tem, which can be seen by the trajectory of h(k) in the 5-
componentvectorspaceformedbytheDiracmatricesΓ.Here Ud = (µ−cosky)σxs0−coskzσzs0
wechooseasurfaceplaneperpendiculartothe xdirectionby −sinkyσysy+(µs−tzsinkz)σysz, (18)
fixing ky and kz, and study the corresponding edge states as Uod = −σxs0+iσysz. (19)
an example. Edge states in the othertwo directionscan also
be studied in this way, and give similar results are obtained. We numerically diagonalize this Hamiltonian and show the
The Hamiltonian terms associated with k give an elliptical four closest doubly degenete energy bands above and below
x
4
2 2 ν+ =ν - =0 with k± = kx ±iky, ǫ0(k) = C + D1kz2 +D2|k+|2 and M(k) =
NL NL M− B k2 − B |k |2. Using Dirac matrices, this Hamiltonian
1 z 2 +
canbewrittenas
µ 1 ν+ N L =ν N- L =1 ν+ N L =ν N- L =0 tz 1 ν +N L =1, ν -N L =0
H(k)= M(k)σ s +A k σ s +A k σ s +A k σ s , (21)
z 0 1 z x z 2 x x x 2 y x y
ν+ =ν - =1 ν+ =ν - =0 where we left out the identity matrix as it only changes the
NL NL NL NL
0 0 shape of the energy bands, not the topology of the system.
0 0.5 1 1.5 2 0 0.5 1 1.5 2
µs µs Thechiraloperatoris givenby theabsentfifth Dirac matrix,
S =σ s .
y 0
FIG.3:Phasediagramwithtz=0(left)andµ=1(right).Thenum- P∗T symmetryishereimplementedbyP=σz. However,
berofDiracconesinedgestatesequalstothetotalwindingnumber
we note that particle-hole symmetry also exists in this case,
ofthemodel,ν = ν+ +ν− . Intheleftpanel,theyellowregion
sum NL NL whichreadsCH∗(−k)C−1 =−H(k),withC =σ s ,andsatis-
withν+ = ν− = 0includesthreedifferentsituations: i)twolines y y
NL NL fiesCC∗ =1.ThusthemodelEq.(21)fallsintotheDIIIclass,
areenclosedbytheNLswhenµ <1andµ<1;ii)allthelinesare
s
outoftheNLswhenµ < 1andµ > 1; andiii)notanylineexsits whichisalsocharacterizedbyaZ invariantin3D[3]. Similar
s
whenµs>1. toourpreviousdiscussion,wewriteH = H1+H2,with
H (k) = M(k)σ s +A k σ s , (22)
1 z 0 1 z x z
H (k) = A k σ s +A k σ s . (23)
2 2 x x x 2 y x y
ThenH hasaNLink −k plane,andthenodesofH produce
1 x y 2
a single line enclosed bythe loop. This givesa NL winding
numberν = 1,whichcorrespondstothetopologicalprop-
NL
ertiesoftheΓpoint.
Summary.-Insummary,wehavestudiedHamiltonianswith
P∗T symmetrywheretwotermsdefineaNLwhichisgapped
out by the other terms. In the presence of chiral symmetry,
these gap terms can drive the system into a chiral TI, which
canbedescribedbyawindingnumberdefinedalongtheNL.
FIG.4:ThespectraofthetwodoubletbandsnearesttoE=0,(a)-(c) This winding number is associated with the geometric rela-
foropenboundary conditionalong x, and(d)-(f)forperiodbound- tion between the NL and the zeroes of the gap terms. We
arycondition.Thethreecolumnsareforthreetopologicallydifferent
investigated a lattice model in detail, which has two NLs in
phases,whilethelefttwoaretopologicallynontrivial. Theparame-
the BZ, each of them with a winding number of 1 or 0 due
tersareµ=1and(a),(d)µ =t =0;(b),(e)µ =t =1;and(c),(f)
s z s z to the gap terms. This winding number corresponds to the
µ =2,t =0.
s z
emergenceofaDiracconeforthesurfacestatesunderOBC.
Finally,we appliedourmethodtothe 3Dtopologicalinsula-
zeroenergyinFig.4,withOBCalongxfrom(a)to(c),andpe- torsoftheBi2Te3family,andshowedtheconnectionbetween
riodicboundaryconditionfrom(d)to(f)ascomparison.Pan- theirtopologicalnatureandtheNLwindingnumber.
els(a)and(b)areforthetopologicallynontrivialphaseswith Acknowledgments.- Financial support from FCT through
ν+ = ν− = 1andν+ = 1, ν− = 0, andthespectrashow grant UID/CTM/04540/2013is acknowledged. S. C. is sup-
NL NL NL NL
edge states with two or one Dirac cone, respectively. Panel portedbyNSFCunderGrantsNo. 11425419,No. 11374354
(c)isforthetopologicallytrivialphasewithν+ = ν− = 0, andNo. 11174360.
NL NL
where there is no edge state connecting the conduction and
valencebands.
A real material example.-Finally, we apply our method to
3DtopologicalinsulatorsoftheBi Te family[33].Thesema-
2 3 [1] M.Z.HasanandC.L.Kane,Rev.Mod.Phys.82,3045(2010).
terials possess time reversalsymmetry and are characterized
[2] X.L.QiandS.C.Zhang,Rev.Mod.Phys.83,1057(2011).
bya Z2 number. Nevertheless,theirtopologicalnatureisde- [3] A.P.Schnyder,S.Ryu,A.Furusaki,andA.W.W.Ludwig,Phys.
terminedbythephysicsnearthetime-reversal-invariantpoint Rev.B78,195125(2008);
Γ(0,0,0)intheBrillouinzone,aroundwhichthelow-energy [4] J.E.MooreandL.Balents,Phys.Rev.B75,121306(2007).
effective Hamiltonian also satisfies a chiral symmetry. This [5] L.Fu,C.L.Kane,andE.J.Mele,Phys.Rev.Lett.98,106803
(2007).
Hamiltonianisgivenby
[6] Z. Wang, X.-L.Qi, and S.-C.Zhang, New Journal of Physics
12,065007(2010).
M(k) A k 0 A k
1 z 2 − [7] P.Hosur,S.Ryu,andA.Vishwanath,Phys.Rev.B81,045120
H(k)=ǫ (k)I + A1kz −M(k) A2k− 0 , (20) (2010).
0 4×4  A20k+ A20k+ −MA(1kk)z −−MA1(kkz)  [8] P1.95D4e5lp2la(2c0e,11D)..Ullmo,andG.Montambaux, Phys.Rev.B84,
5
[9] G.Zhang,Z.Song,Phys.Rev.Lett.115,177204(2015). 235126(2011).
[10] L.LiandS.Chen,Phys.Rev.B92,085118(2015). [23] K. Mullen, B. Uchoa, and D. T. Glatzhofer, Phys. Rev. Lett.
[11] L.Li,C.Yang,S.Chen,Europhys.Lett.112,10004(2015). 115,026403(2015).
[12] L.Li,C.Yang,S.Chen,Eur.Phys.J.B(2016)89:195. [24] R.Nandkishore,Phys.Rev.B93,020506(2016).
[13] B.-J. Yang, N. Nagaosa, Nature Communications, 5, 4898 [25] Z.YanandZ.Wang,Phys.Rev.Lett.117,087402(2016).
(2014). [26] C. Fang, Y. Chen, H.-Y. Kee, and L. Fu, Phys. Rev. B 92,
[14] S.Ganeshan,S.DasSarma,Phys.Rev.B91,125438(2015). 081201(2015).
[15] X.Wan,A.M.Turner,A.Vishwanath,andS.Y.Savrasov,Phys. [27] C.Fang, H.Weng, X. Dai, andZ.Fang, ChinesePhys.B 25,
Rev.B83,205101(2011). 117106(2016).
[16] S.M.Young, S.Zaheer, J.C.Y.Teo, C.L.Kane, E.J.Mele, [28] J.L.Lu,W.Luo,X.Y.Li,S.Q.Yang,J.X.Cao,X.G.Gong,
andA.M.Rappe,Phys.Rev.Lett.108,140405(2012). andH.J.Xiang,arXiv:1603.04596.
[17] P.Delplace,J.Li,andD.Carpentier,Europhys.Lett.97,67004 [29] L.Li,M.A.N.Arau´jo,Phys.Rev.B94,165117(2016).
(2012). [30] Arfken, G. Mathematical Methods for Physicists, 3rd ed. Or-
[18] J-.M.Hou,Phys.Rev.Lett.111,130403(2013). lando,FL:AcademicPress,pp.211-217,(1985).
[19] T.MorimotoandA.FurusakiPhys.Rev.B89,235127(2014). [31] M.V.Berry,Proc.R.Soc.London,Ser.A392,45(1984).
[20] S. M. Young and C. L. Kane, Phys. Rev. Lett. 115, 126803 [32] R. S. K. Mong and V. Shivamoggi, Phys. Rev. B 83, 125109
(2015). (2011).
[21] S. Ganeshan and S. Das Sarma, Phys. Rev. B 91, 125438 [33] H.Zhang,C.-X.Liu,X.-L.Qi,X.Dai,Z.FangandS.-C.Zhang
(2015). Nat.Phys.5,438(2009).
[22] A.A.Burkov,M.D.Hook,andLeonBalents,Phys.Rev.B84,