Global Phase Diagram of the Extended Kitaev-Heisenberg Model on Honeycomb Lattice

0.29 MB·English

by  Jie Lou

Global Phase Diagram of the Extended Kitaev-Heisenberg Model on Honeycomb
Lattice
Jie Lou, Long Liang, Yue Yu, and Yan Chen
Department of Physics and State Key Laboratory of Surface Physics, Fudan University, Shanghai 200433, China and
Collaborative Innovation Center of Advanced Microstructures, Nanjing 210093, China
(Dated: January 29, 2015)
WestudytheextendedKitaev-Heisenberg(EKH)quantumspinmodelbyaddingbond-dependent
off-diagonal Heisenberg term into the original KH model, which was recently proposed to describe
the honeycomb Iridates. A rigorous mathematical mapping of spin operators reveals the intrinsic
5
symmetry of the model Hamiltonian. By employing an unbiased numerical entanglement renor-
1
0 malization method based on tensornetwork ansatz, we obtain theglobal phasediagram containing
2 eightdistinctquantumphases. Byusingthedualmappingofspinoperators, eachoftheindividual
magnetic phase in the global phase diagram can be clearly understood. At last, we show that a
n
valencesolidstateemergesasthegroundstateinthequadro-criticalregionwheremultiplemagnetic
a
phases compete most intensively.
J

8 PACSnumbers: 75.10.Jm,75.10.Nr,75.40.Mg,75.40.Cx
2

] The family of layered honeycomb iridium oxides coupling constants. By using this dual transformation,
l
e Na2IrO3 and Li2IrO3 have attracted great interest both eachoftheindividualmagneticphaseintheglobalphase
-
r experimentallyandtheoretically[1–17]duetoitspossible diagram can be clearly understood. For instance, some
t
s relevance of Kitaev spin liquid [18]. Earlier theoretical seemingly complicated magnetic orders can be mapped
t. studies proposedthat the exchangefromdirectd-orbital into simple antiferromagnetic (AFM) order or ferromag-
a
overlapofIrionsmayresultinaKitaev-Heisenberg(KH) netic(FM)order. Furthermore,boththefinitesizeeffect
m
model [19]. Extensive study of such model has shown a andthe results for infinite system arediscussed. At last,
-
d variety of fascinating phenomena, including an extended weshowthatavalencesolidstateemergesastheground
n spin liquid phase and quantum phase transitions into state in the quadro-critical region where multiple mag-
o several well-understood magnetic ground states. Exper- netic phases compete most intensively.
c imentally, neutron scattering data revealed the presence Model Hamiltonian.— Our model is the EKH model
[
ofthezigzagphaseinNa2IrO3[1–3]. However,itappears defined on a honeycomb lattice, for which the Hamilto-

1 to be hard stabilize within the HK model; one may con- nian is
v
sider additionalexchange paths [13] or further neighbor-
90 ing hoppings [12]. Recently, a generic bond-dependent H = X[JSi·Sj +KSiγSjγ +Γ(SiαSjβ +SiβSjα)], (1)
9 off-diagonal exchange term is included in the extended hi,ji
6 Kiteav-Heisenberg (EKH) model [16]. By using a com-
0 bination of classical techniques and exact diagonaliza- wherehi,jidenotesthenearest-neighbor,andtermswith
1. tion, 120◦ and incommensurate spiral order, as well as strength J, K and Γ represent Heisenberg, Kitaev and
bond dependent off-diagonal exchange interactions, re-
0 extendedregionsofzigzagandstripyorderareobtained.
5 This EKH model can be regarded as the minimal model spectively. We use γ to distinguish one spin direction
1 on each bond for the Kitaev exchange, as shown in
todescribetheessentialphysicsofhoneycombiridiumox-
v: ides. However,the effect of bond-dependent off-diagonal Fig. 1. We have parameterized the energy scales so that
i exchangetermisunclearandtheobtainedphasediagram J2+K2+Γ2 = R2, J = Rcosθsinφ,K = Rcosθcosφ,
X
Γ = Rsinθ. φ controls the ratio between Kitaev and
is not clearly understood.
r Heisenberg terms.
a
In this letter, we perform a systematic study of the It is well known that at the Kitaev limit J = Γ =
EKH quantum spin model on honeycomb lattice by us- 0,K 6= 0, the ground state of system is a topological
ing the multiscale entanglement renormalization ansatz spin liquid which has no long range order. Away from
(MERA)[20]. We obtainedaglobalphasediagramcon- thelimit,severallongrangemagneticordersmayemerge
taining eightdistinct quantumphases. Arigorousmath- and compete with each other. In addition to the AFM
ematical mapping of spin operators reveals the intrinsic and FM order, the system also harbors stripy (ST) and
symmetry of the model Hamiltonian. In particular, by zigzag(ZZ)phases,inregioncenteredattwopointsK =
rotationofspinoperatorsobeyingcertainrules,theaddi- ±2J. It can be shown that exactly at these two points,
tional bond-dependent off-diagonalHeisenberg exchange the original KH hamiltonian can be mapped into a pure
term can be mapped to Heisenberg or Kitaev interac- Heisenberg model with the coupling constant J′ = −J.
tions and vice versa. So the original EKH model can be The trick is to divide the system into four sublattices,
mapped to the model of same form but with modified eachflip spin components in accordanceto certainrules,
2
-180 -90 0 90 180
90
[ZZ I] [ZZ II]
A B
E 120[FM] AF[FM]
FM[120] [ZZ IV]
[ZZ III] F
SL D
SL C 0
ZZ ST SL
AF[120]
FM[AF] [ST I]
G H
[ST II] 120[AF]
-90
FIG. 2: Global phase diagram from MERA calculations of
24-sitessystem. Longitudeφandlatitudeθoftheworldmap
parameterize the couplings of J, K and Γ terms. In the up-
per panel, phases are described from the perspective of the
originalspinmodel,whereasthelowerpanelshowsthephase
map when the rotated spin model is considered. Regions of
same type of magnetic order are shown in same color, la-
beled in the bottom. Diamonds notates special points where
theoriginal model can be simplified tosimple Heisenberg in-
teractions. Red squares mark parameter sets for which we
show spin configurations of ground states obtained from ER
FIG.1: (coloronline). MappingtheoriginalEKHmodeltoits
calculations. Phaseboundariesmarkedbyredwidelinescor-
counterpart through spin rotations. The upper panel stands
responds to first order phase transition.
for the original spin model, where spin directions γ(α,β) for
Kitaev and bond-dependent (Γ) exchanges are denoted by
different colors. By performing the spin rotation marked by
numbers 0 to 5, following the rule listed in the table, we can that its ground state is classical FM state with energy
maptheoriginal spinmodeltoitscounterpart. Differentcol- per-siteexactlyequalto−0.375. InFig.3,panel(A),we
ors are painted for each site to denote rules of sign change, showdetailsofthegroundstate,fromperspectivesbefore
followingwhichZZ/SPphasescanbemappedintoAFM/FM. and after rotations of spins, for comparison. In terms of
Shades in the background in the upper panel represent the
original representation, the spin configuration of ground
MERA structure we adopted in numerical calculation.
state is rather complicated, and long range magnetic or-
der can not be easily distinguished. On the other hand,
in the mapped representation, the spin configuration of
asdenotedbydifferentcolorsinFig.1. Inthatsense,the
ground state is a simple ferromagnet with spins align-
ground states of KH model can be well understood in
ing in same direction. Hence, the intrinsic symmetry
terms of the competitions among simple magnetic long-
of the model Hamiltonian may help us to better under-
range orders.
stand the seemingly complicated magnetic orders of the
Mapping ofSpin Operators.—Inordertobetterunder-
system induced by the additional off-diagonal exchange
standtheeffectofΓ-term,weappliedrotationofthespin
operatorsSγ →Sα (See the table inFig.1)so thatthe Γ term. There exists several such parameter sets where
the original model can be simplified to pure Heisenberg
term can be mapped back to Kitaev exchange term. By
models (shown in Fig. 2). These special cases offer some
doingso,theoriginalhamiltonianismappedtothemodel
additional hints for solving EKH model, but it becomes
of same form but with modified coupling constants,
rather complicated when various magnetic orders com-
H = X J′Sei·Sej +K′SfiγSfjγ +Γ′(SfiαSfjβ +SfiβSfjα)(,2) pnaetteess.with each other, especially when the Γ term domi-
<ij>
Globalphasediagram.—WeuseanunbiasedMERAal-
wherestrengthofeachtermisrelatedto originalparam- gorithmtofullyexplorethephasediagramofEKHmodel
etersasJ′ =−Γ,K′ =(K−Γ−J),andΓ′ =−J. Notice on honeycomb lattice. In the upper panel of Fig.1, we
that the Heisenberg and bond-dependent exchange cou- show the ER tensor network we construct for the hon-
plings are swapped through the mapping. eycomb lattice. This method is capable of capturing the
The mapping can be extremely helpful to understand Kitaevspinliquidandotherphaseswherequantumfluc-
the characteristics of ground state. For instance, we tuations play an important role. The system sizes we
take the limit J = 0,Γ = K = 1. Under rotations of calculated range from 24 to 96.
spin operators, the mapped Hamiltonian corresponds to The global phase diagram we numerically obtained is
a simple Heisenberg model with FM interaction H = shown in Fig. 2, which is presented in format of a world
−Phi,jiSiSj. Indeed, numerical calculation confirms map as functions of two parameters θ and φ. The anti-
3
ferromegnetKitaev point (spin liquid) is located exactly
at the map center, and the original Kitaev Heisenberg
model with Γ = 0 corresponds to the equator line. In
Fig. 2, we present two phase maps and each of them be-
longs to a different spin representation. The very same
A1 A2 B1
state show different magnetic long range orders without
(upper panel) or with (lower panel) spin rotations. Here
we denote these two different orders with labeling “X”
and “Y¯”. By performing the spin rotations, we notice
that Heisenberg and Γ exchanges are swapped, namely,
J =−Γ and vice versa.
C1 C2 D1
eAFM/FM and 120◦ phases—Due to the fact that the
EKHmodelismappedtoitself,weexpectthephasemap
tobeselfconsistentbyspinrotations. Suchmappingcor-
respondenceisclearlyobservedforAFMandFMorders,
which occupy a large portionof the phase map. It is im-
portant to notice that one block of certain order (AFM E1 E2 E3
for example) in the original map can be further sepa-
rated into sections of different orders (120 and FM) in
the phase map for rotated spins, vice vegrsa. g
This feature seems to be contradictory, but careful
analysis reveals that these two regions are indeed of dif-
F1 F2 F3
ferent magnetic orders. As a matter of fact, a first order
phase transitionis clearlyshownwhen the phase bound-
FIG. 3: Spin configurations and nearest neighbors’ correla-
ary is approached and eventually crossed (see Figure in tionshSi·Sjiobtainedfromnumericalcalculationfor24-sites
supplementary material). The physical reason can be system. On each site, spin hSii is plotted as combination of
understood as the following. For the FM order, spins a circle (proportional to hSizi) and a vector (proportional to
hSiiareparallelsforallsites. Naturally,gweexpectequal hSixyi). Spin correlations are notated by bonds connecting
e spin pairs, whose width are proportional to the strength of
spincomponentsforeverysite,namely Sfix =Sfjx. Onthe hSi · Sji. Red/green color stands for FM/AF correlation,
other hand, the very phase is of the AFM order in term
respectively. Panel A ∼ F show results achieved in differ-
of original spins, which require |Six| = |Sjx|. Notice that ent phases, their corresponding parameter sets can be read
rotations of spins (Sx →−Sy, for example) are involved from the marks on phase maps. The number of each nota-
i fi tionstandsfortheperspectivefromwhichthegroundstateis
whenoriginalhamiltonianismappedtothenewone. As
viewed: “1” notates the original spin picture, “2” stands for
aresult,the onlypossibleoutcomeforthese tworequire-
the rotated spin picture, and “3” means the result of zigzag
ments to be met at the same time is that all spins order mapping (for the rotated spin model), described in Fig. 1.
along (1,1,1) direction, namely Sx = Sy = Sz. In an- Results for larger systems can befound in Sup. material.
i i i
other word, the bulk of FM/AF phase is of Ising type.
This is to be expected, asgintroduction of Γ term breaks
SU(2) symmetry of original Heisenberg magnets. AFM north of the equator, and a U(1) AFM south. A
Similarsituationoccursforthe120◦/AFM phase. The firstorderphasetransitionshowsupatthephasebound-
120◦phaseisnamedafterthefactgthatallanglesbetween ary. They are different physicalorderswith distinct bro-
ngext-nearest neighbor spins S are 2π/3 in this order. an kensymmetry. SuchphenomenaisuniversalforallAFM
illustratedinFig.3. Suchphaeseappearsouthofeuqator, andFMorders,inbothoriginalandrotatedspinmodels.
which represent a positive nearest neighbor Heisenberg The finite size effect is less remarkable as shown in our
couplingsJ =−Γ,unfavoringtheFM phase. Therefore, numerical calculation. All magnetic orders remain the
spinsinorieginalpictureordersinthgeplaneperpendicular same, except a slightly shifts of phase boundaries [21].
to (1,1,1), as a compromise. (Notice, the system is still Zigzag/stripy phases andfrustration—Whileswitching
in the AFM phase with all nearest neighbors spins an- longitude and latitude is helpful to understand the rela-
tiparallel.) ItcanbeshownthatwhenceSx+Sy+Sz =0 tion between two phase maps, especially for AFM/FM
i i i
(true for any vector in plane perpendicular to (1,1,1)), orders,itdoesnotworkwellwhenzigzag/tripyphasesare
after rotations listed in tables, the 120◦ order emergies. considered. Asamatteroffact,thesephasesarefarmore
In this case, the U(1) symmetry, insgtead of SU(2) sym- complicatedthanNeel-likeordersdiscussedabove. From
metry, is broken in the AFM order. numerical calculation, we find that spins are much more
Due to effect of Γ term, the orginal SU(2) symmetric prone to tilting awaythan strictly order alongone direc-
NeelAFMorderbreaksintotwodifferenttypes: anIsing tion. Meanwhile, various candidates of magnetic orders
4
compete for the ground state, and energy gaps between 0.25
them are relatively small. The reason for such different
behaviors can be found in the mapping when rotations 0.2
of spins are done. As mentioned before, such discrep-
ancy may not affect AFM and FM orders much, (as we 0.15 NN_M
ozibgszeargvepdhainseoaunrdnsutrmipeyricpahlacsealsciuglnaitfiiocnan),tlby.utAsdoa cmhaatntgeer <SS>ij -J/K 2 234NNNNNN___MMM
of fact, the new patter^n of Kit^aev and Γ interactions is 0.1 SC2 LLL===234468 1.6 5NNNN__OM
frustrated in terms of zigzag/stripy order. C/1 2NN_O
0.05 S 3NN_O
For the original Kitaev Heisenberg model, a four- 1.2 4NN_O
sublattice mapping between the ST/ZZ order and the 5NN_O
0.7 0.8 0.9 1
0
AFM/FM order is well known and discussed [19]. Nat- 0.7 0.8 0.9 1 1.1 1.2
-J/K
urally, we expect same type of operation can be applied
for the rotatedEKH model. Unfortunately it is not pos-
FIG. 4: Evolution of spin correlations (SC) hSi·Sji when
sible to obtain such a perfect mapping in this case. One
system undergoes a transition from magnetic long range or-
can show that, when we flip spin operators one by one
ders to the VBS phase. K and Γ are fixed to be 1, and J
along an hexagon, according to arrangement of Kitaev is tunedthrough thephase transition. Different colors repre-
interaction γ on each bond, we can not obtain self- sentsvariousdistancesbetweentwospinsiandj,fromnearest
i,j
consistency after a full circle. As a result, we have to neighbor (NN) to 5th nearest neighbor (5NN), the maximal
sacrifice certain bonds in order to map the rotated spin distance in a periodic 24-site lattice. At −J/K = 0.95, a
g
model to its ZZ/ST counterpart. One way to make the phasetransitionbetweenFM andtheVBSorderoccurs,and
compromise igs shofwn in the lower panel of Fig. 1. We correlations forrotated spins(SCnotated bysuffixM) drops
abruptly. When J/K goes over 1, the system enter the FM
keep the sign-flipping rule along chains following certain
order, and SCs for original spins (suffix O) do not decay by
direction,andreversespinoperatorswhenjumpbetween
increasing theseparation. The inset shows theratio between
chains. Clearlythereexists3degenerateways(directions) average spin correlation intra and inter a plaquette hexagon
forsuchsignchangerules. Onlyone isplottedinthe fig- SC1/SC2. Results for various system sizes are plotted.
ure, and shades indicate the chosen direction of ZZ/ST
g f
chains.
Such a compromise sacrifices certain Kitaev interac- capability,wearenotabletoprobetheregionthoroughly.
tions intra-chain. Numerical study shows that such a Valence bond solid. —Asexpected,multiplemagnetic
mapping is energy-favorable for most zigzag region in long range orders compete most intensively close to two
therotatedspinmodel. ExamplesofsuchZZ phasesare points −J = K = Γ = 1 and J = −K = −Γ = 1.
shown in Fig. 3. In the rotated spin modgel, we expect These two special points corresponds to verging points
the ZZ phase to appear nearthe point K =−2J,Γ=0, of multiple phases, as shown in Fig. 2. It is not clear
whicghistranslatedtoJ =0,K =−1,Γ=e 1. Noteiceethat what is to be expected from analytical analysis, and our
thespin-spincorrelationshowsthatintra-chains’correla- numerical calculation shows that small regions surround
tions are significantly larger than inter-chains’, which is these multi-critical points harbor the VBS phase. A pe-
consistent with the mapping we made. It is necessary to riodic strong-weak shifting pattern for nearest neighbor
mentionthatduetothisfrustration,thesestatesbecome spin correlation can be observed in Fig. 3, panel D2.
lessstablewhenwetuneparametersawayfromthezigzag More precisely speaking, correlations inside a plaquette
point. spin bend away and new magnetic order appears. (hexagonunitcell)islargerthanintra-plaquetteones. As
^
One notice that for J < 0 and J > 0, the zigzag phase a result, we expect long range spin correlations to drop
forrotatedspinmodelalsoemergesasdifferentmagnetic fast when the distance between two separated spins is
order in originalspin picture. The reason is very similar increased. In contrast, for all magnetic orders discussed
to the AF case explained above. Due to perturbation of above, correlations decay much slower when spins are
Γ=−J term, system magnetize is different quadrant in further separated. In Fig. 4, we show such results from
e
these two cases. Details can be found in [21]. numericalcalculationsandatransitionbetweenthemag-
The mapping we choose is not necessarily the unique netically ordered phase (120/FM) and the VBS phase
g
one, however, especially in larger systems. In small sys- can be distinguished.
tem with 24 and 36 sites, numerical calculations shows To further clarify this, we define spin correlations in-
that the former mentioned states are the ground states tra/inter plaquette hexagon (P) as SC1 and SC2, where
spanning the zigzag region for the rotated spin model. SC1 = PPi=PjSi ·Sj and SC2 = PPi6=PjSi · Sj. In
However for larger systems, various magnetic order can- subset of Fig. 4, we show that for various system sizes
didates contradict to such mapping also appears, with we calculated, the ratio between SC1 and SC2 increase
reasonably good energy [21]. With current calculation significantly when the VBS region is approached. It is
5
worth mentioning that the local spin expectation value [5] J.Reuther,R.Thomale,andS.Trebst,Phys.Rev.B84,
hS i vanishes in these criticalregions(as shownin Fig. 3 100406 (2011).
i
for24-sitessystem),asimilarfeatureobservedfortheKi- [6] F. Trousselet, G. Khaliullin, and P. Horsch, Phys. Rev.
B 84, 054409 (2011).
taev spin liquid. However, long range spin correlations
[7] R.Schaffer,S.Bhattacharjee, andY.B.Kim,Phys.Rev.
reveal that they are intrinsically different. We have to
B 86, 224417 (2012).
clarify that our calculation is limited to relatively small [8] C.C. Price and N.B. Perkins, Phys. Rev. Lett. 109,
system sizes, as a result, finite size effect boundary con- 187201 (2012); C.PriceandN.B.Perkins,Phys.Rev.B
ditions may play important roles. We can not rule out 88, 024410 (2013).
the possibility of spin liquid in this critical region. [9] Y. Singh and P. Gegenwart, Phys. Rev. B 82, 064412
Conclusions. —We numerically studied the global (2010).
[10] Y. Singh,S. Manni, J. Reuther,T. Berlijn, R.Thomale,
phase diagram of the EKH model containing eight dis-
W. Ku, S. Trebst, and P. Gegenwart, Phys. Rev. Lett.
tinct quantum phases. A rigorous mathematical map-
108, 127203 (2012).
ping of spin operators reveals the intrinsic symmetry of [11] S.Bhattacharjee, S.-S.Lee, Y.B.Kim,NewJ.Phys.14,
the model Hamiltonian. By using the dual mapping of 073015 (2012).
spin operators, each of the individual magnetic phase [12] I. Kimchi and Y.Z. You, Phys. Rev. B 84, 180407(R)
in the phase diagram can be clearly understood. Fi- (2011).
nally, we showed that a valence solid state emerges as [13] J.Chaloupka,G.Jackeli, G.Khaliullin,Phys.Rev.Lett.
110, 097204 (2013).
the ground state in the critical region where multiple
[14] V.M. Katukuri, S. Nishimoto, V. Yushankhai, A. Stoy-
magnetic phases compete most intensively.
anova, H. Kandpal, S. Choi, R. Coldea, I. Rousochatza-
Acknowledgments.— This work was supported by
kis,L.Hozoi,J.vandenBrink,NewJ.Phys.16,013056
the State Key Programs of China (Grant No.
(2014).
2012CB921604), the National Natural Science Founda- [15] H. Gretarsson, J.P. Clancy, X. Liu, J.P. Hill, E. Bozin,
tion of China (Grant Nos. 11274069,and 11474064). Y. Singh, S. Manni, P. Gegenwart, J. Kim, A.H. Said,
D. Casa, T. Gog, M.H. Upton, H.-S. Kim, J. Yu, V.M.
Katukuri, L. Hozoi, J. van den Brink, and Y.-J. Kim,
Phys. Rev.Lett. 110, 076402 (2013).
[16] J.G.Rau,E.Lee,H.Y.Kee,Phys.Rev.Lett.112,077204
[1] X. Liu, T. Berlijn, W.-G. Yin, W. Ku, A. Tsvelik, Y.-J. (2014).
Kim, H. Gretarsson, Y. Singh, P. Gegenwart, and J.P. [17] N.B.Perkins,Y.SizyukandP.W¨olfle,Phys.Rev.B89,
Hill, Phys. Rev.B 83, 220403 (2011). 035143 (2014).
[2] F. Ye, S. Chi, H. Cao, B.C. Chakoumakos, J.A. [18] A. Kitaev, Ann.Phys.321, 2 (2006).
Fernandez-Baca, R. Custelcean, T.F. Qi, O.B. Korneta, [19] J. Chaloupka, G. Jackeli, and G. Khaliullin, Phys. Rev.
and G. Cao, Phys. Rev.B 85, 180403 (2012). Lett. 105, 027204 (2010).
[3] S.K.Choi,R.Coldea,A.N.Kolmogorov,T.Lancaster,I.I. [20] G.Vidal,Phys.Rev.Lett.99,220405(2007);G.Evenbly
Mazin, S.J. Blundell, P.G.Radaelli, Y.Singh,P.Gegen- andG.Vidal,ibid102,180406(2009);G.EvenblyandG.
wart, K.R. Choi, S.-W. Cheong, P.J. Baker, C. Stock, Vidal,ibid104,187203(2010);G.EvenblyandG.Vidal,
and J. Taylor, Phys.Rev. Lett.108, 127204 (2012). Phys. Rev.B 79, 144108 (2009).
[4] H.-C. Jiang, Z.-C. Gu, X.-L. Qi, and S. Trebst, Phys. [21] See the Supplemental Material for the detailed descrip-
Rev.B 83, 245104 (2011). tion of phase diagram.