Clues and criteria for designing Kitaev spin liquid revealed by thermal and spin excitations of honeycomb iridates Na$_2$IrO$_3$

1.9 MB·

by  Youhei Yamaji

Clues and criteria for designing Kitaev spin liquid revealed by
thermal and spin excitations of honeycomb iridates Na IrO
2 3
Youhei Yamaji1, Takafumi Suzuki2, Takuto Yamada2, Sei-ichiro Suga2, Naoki Kawashima3, and Masatoshi Imada4
1Quantum-Phase Electronics Center, The University of Tokyo, Tokyo, 113-8656, Japan
2Graduate School of Engineering, University of Hyogo, Hyogo, Himeji, 670-2280, Japan
3Institute for Solid State Physics, The University of Tokyo, Chiba, 277-8581, Japan and
4Department of Applied Physics, The University of Tokyo, Tokyo, 113-8656, Japan
(Dated: January 22, 2016)
Contrary to the original expectation, Na IrO is not a Kitaev’s quantum spin liquid (QSL) but
6 2 3
shows a zig-zag-type antiferromagnetic order in experiments. Here we propose experimental clues
1
andcriteriatomeasurehowamaterialinhandisclosetotheKitaev’sQSLstate. Forthispurpose,
0
we systematically study thermal and spin excitations of a generalized Kitaev-Heisenberg model
2
studied by Chaloupka et al. in Phys. Rev. Lett. 110, 097204 (2013) and an effective ab initio

n Hamiltonian for Na IrO proposed by Yamaji et al. in Phys. Rev. Lett. 113, 107201 (2014),
2 3
a by employing a numerical diagonalization method. We reveal that closeness to the Kitaev’s QSL
J
is characterized by the following properties, besides trivial criteria such as reduction of magnetic

1 orderedmomentsandN´eeltemperatures: (1) Twopeaksinthetemperaturedependenceofspecific
2 heat at T and T caused by the fractionalization of spin to two types of Majorana fermions. (2)
ℓ h
In between the double peak, prominent plateau or shoulder pinned at Rln2 in the temperature
] 2
dependenceofentropy,whereRisthegasconstant. (3)Failureofthelinearspinwaveapproximation
l
e at the low-lying excitations of dynamical structure factors. (4) Small ratio T /T close to or less
ℓ h
-
than 0.03. According to the proposed criteria, Na IrO is categorized to a compound close to the
r 2 3
t Kitaev’sQSL,andisproventobeapromisingcandidatefortherealizationoftheQSLiftherelevant
s
. material parameters can further be tuned by making thin film of Na2IrO3 on various substrates or
t
a applying axial pressure perpendicular to the honeycomb networks of iridium ions. Applications of
m thesecharacterization to (Na1−xLix)2IrO3 and other related materials are also discussed.
-
d PACSnumbers: 75.10.Jm,75.10.Kt,75.40.Gb,75.40.-s
n
o
c I. INTRODUCTION rana fermions propagating in the background of static
[ Z gauge fields (0 or π flux) that are also written by the
2
localized Majorana fermions4. Reflecting the fractional-
1 Enormous efforts to realize quantum spin liquids
ization of quantum spins into the Majorna fermions, the
v (QSLs) have been made since the pioneering work of
2 Anderson and Fazekas1,2. Geometrically frustrated in- spin correlations for further neighbor sites become ex-
1 teractions, for instance antiferromagnetic Heisenberg in- actlyzerointhe Kitaev’sQSLstates,while thoseforthe
5 nearest neighbor survive.
teractions on a variance of triangular lattice, have been
5
studied as one promising way to realize the QSL state. In addition to the ground state properties, low-energy
0
Recently, the Kitaev model on a honeycomb structure3 excitations4 anddynamics5 ofbothgapped/gaplessQSL
.
1 has attracted attention, because the ground state is ex- stateshavealsobeeninvestigatedbyanalyticalmethods.
0
actly provento be ina QSL phase. Inthe Kitaevmodel, The ground state is in the 0 flux sector and the lowest
6
two spins on the nearest neighbor sites i, and j inter- excitationbyaspinflipfromtheKitaev’sQSLstateisex-
1
v: aKctSbzyStzh,ewIhsienrge-tthypeeIsiinntgeraancitsiootnr,oKpyxSaixxiSsjxd,eKpeynSdiysSojyn,atnhde panresitsienderbayntadMdainjgoraan‘laocfearlimzeiodn’.π-Sfliunxceptahireaecxccoitmedpaπn-yfliunxg
i z i j
X three differentbonding directionsthat composethe hon- pairdoesnotpropagateintheKitaevmodel,nondisper-
eycomb structure. These anisotropic interactions cause sive mode appears in the low-lying excitation, where the
r
a a strong frustrated effect, distinct from the typical geo- excitationenergyisnonzero. Suchagappedexcitationis
metrical frustration. confirmed in the dynamical spin structure factor (DSF)
The QSL state in the Kitaev model contains two dis- that can be observed in inelastic-neutron-scattering ex-
tinct phases,namely QSLs with gapless excitations from periments. In the gapped QSL phase, a δ-function peak
thegroundstateandthosewithonlygapfulexcitations3. emergesatthelowestexcitationenergyoftheDSF,while
The gapless QSL appears, when the system is located a sharp peak with a long tail generating gapless excita-
around the symmetric point where the bond-depending tions derivedfrom the incoherentpart is observedin the
interactions have the same magnitude, K = K = K . gapless QSL phase.
x y z
If one coupling constant becomes much larger than the Quiterecently,thermalpropertiesfortheKitaevmodel
otherones,forexample|K |≫|K |=|K |,thegapped on the honeycomb structure have been investigated by
z x y
QSL state is stabilized. Interestingly, the both QSL the quantum Monte Carlo calculations6. The Kitaev
states, regardless of whether the excitation is gapped model always exhibits a two-peak structure in the spe-
or gapless, can be described by noninteracting Majo- cificheat,associatedwiththefractionalizationofasingle
2
quantum spin into two types of Majorana fermions; one HamiltonianofNa IrO ,weemployasimple generalized
2 3
is the itinerant (dispersive) Majorana fermion and the Kitaev-HeisenbergHamiltonian proposed by Chaloupka,
other is the localized (dispersionless) Majorana fermion. Jackeli, and Khaliullin, in Ref. 9 and an ab initio effec-
The two peaks representthe two crossovertemperatures tive Hamiltonian proposed in Ref. 14. In the previous
associated with the growth of short-range spin correla- works14,16, we discussed the accuracy of these effective
tions (or thermal excitations of the itinerant Majorana models and concludedthat the ab initio model proposed
fermions) around the high temperature peak and freez- in Ref. 14 explains not only thermal properties but also
ing of flux (or the thermal excitation of the localized dynamics of this compound.
Majorana fermions) around the low-temperature peak, Bycomparingthetemperaturedependencesofspecific
respectively. The above characteristic features and pic- heatC andtheequal-timespincorrelationsaswellasthe
tures for the Kitaev’s QSL states are expected to be ro- dispersionof the spin excitation basedon the linear spin
bust against small perturbations such as magnetic fields wave approximation and the dynamical spin structure
andHeisenberg-typeinteractions. However,detailsofthe factors S(Q,ω), we identify three distinct characteristic
stability have not been well understood yet for real ma- regions in the phase diagram of the effective Hamiltoni-
terials. ans: In addition to the spin liquid phase, the magneti-
InthesearchforrealizationoftheKitaev’sQSLstates, cally ordered phase is classified to two distinct regions.
Na IrO has attracted attention as one candidate mate- Within the magnetically ordered phases of the general-
2 3
rial. In Na IrO , Ir4+ ion can be expressed as a pseu- ized Kitaev-Heisenbergmodel, the system is classifiedto
2 3
dospin with the total angular momentum one-half7. We the category I, when the quantum spin system shows a
call this pseudospin just as ‘spin’ hereafter. IrO oc- single peak structure in the temperature dependence of
6
tahedrons in Na IrO are built into a planar structure C. Ifaquantumspinsystemshowsatwo-peakstructure
2 3
parallel to the ab plane and Ir4+ ions constitute a hon- in C despite its magnetic order, the system is classified
eycomb structure7,8. Furthermore, the neighboring IrO to the category II. When the ground state is the Ki-
6
octahedrons are connected by sharing two oxygenatoms taev’sQSL,thesystemisclassifiedtothethirdcategory,
on an edge. The two oxygen atoms make bridges con- namely,thecategoryIII.AssummarizedinTableI,C for
necting the neighboring Ir atoms and both Ir-O-Ir bond the systems in the category III have two-peak structure
angles are nearly equal to 90◦. Because of the Ir-O-Ir commonly to the category II.
bonds, the perturbative process generates three differ- As shown later, in the category I, the low-lying ex-
entanisotropicinteractionsbetweenIr4+ ions depending citations of the generalized Kitaev-Heisenberg model in
on the Ir-Ir bond directions. In addition, Ir4+ ions in- S(Q,ω), which is induced by flipping a spin, are well
teract via direct overlap of their orbitals, generating the describedbyusingtheconventionallinearspinwavethe-
Heisenberg-type interaction as well. Thus, both of the ory,i.e.,successfullyinterpretedasdispersionof(nearly)
Kitaev-typeandHeisenberg-typeinteractionsemergebe- free magnons. In contrast to the category I, a system
tween the neighboring Ir4+ ions9, leading to the Kitaev- categorized as the category II shows that the low-lying
Heisenberg model. excitationsinS(Q,ω)arenotcapturedbythelinearspin
Contrary to the initial predictions7,8, Na2IrO3 under- wavetheoryonaqualitativelevel. Thebreakdownofthe
goes a magnetic phase transition to a zigzag antiferro- linear spin wave theory is a common property in both
magnetic order at TN ∼ 15K10,11. In order to under- categories II and III, although the spin wave analysis
stand the zigzag ordering, several effective models have has been employed in comparing effective Hamiltonians
been proposed and examined so far9,10,12–15. Some of withexperimentalresultsofA IrO3(A=Na,Li)17,18 and
2
them have succeeded in explaining the thermodynamic another Kitaev’s QSL candidate α-RuCl 19. The above
3
quantities such as the specific heat and/or the magnetic different categories are not necessarily separated by the
susceptibility9,12,14. phase transition: Categories I and II are separated from
Several ab initio derivations of effective spin Hamilto- III by a quantum phase transition, while the states in
nians for Na IrO have shown that low-energy physics the categories I and II can be connected smoothly. All
2 3
of Na IrO is roughly described by dominant Kitaev’s ofhighlygeneralizedKitaevmodelstreatedinthispaper
2 3
Ising-type exchange interactions, K = K & K ≃ −30 can be represented by one of these three empirical cate-
x y z
meV14, while other much smaller interactions including gories consistently in all the physical quantities studied.
theHeisenbergexchangeeventuallyholdthekeyfordriv- Thus, the spin excitation spectra would provide us with
ing the zigzag order. Therefore, one intriguing challenge useful supplementary data for the classification.
is to figure out a guideline for materials design that en- Wethenexaminethenatureofthegroundstateofthe
ables the Kitaev’s QSL, against the small interactions abinitioHamiltonianofNa IrO inlightofourproposed
2 3
through control of them by starting with Na2IrO3. categorization. We show that the ab initio Hamiltonian
In this paper, we first report the characterization of ofNa IrO belongstothecategoryII.Thissupportsthat
2 3
Na IrO based on our numerical studies. Focusing on a better chance ofmaterialdesignto realize the Kitaev’s
2 3
ab initio effective Hamiltonians for Na IrO , we evalu- QSL may exist through a realistic tuning of the mate-
2 3
ate how close Na IrO is located from the Kitaev’s QSL rial parameters of Na IrO . In this context, we propose
2 3 2 3
phase by introducing several criteria. For the effective that the temperature dependence of the entropy as well
3
as the energy scales measuredby the ratioof the higher- 2. Ab initio Hamiltonian of Na IrO
2 3
and lower-temperature peaks in C give a quantitative
measureoftheclosenesstotheKitaev’sQSLsinthecat- LetusconsiderahighlygeneralizedformoftheKitaev-
egory II. We believe that this proposal offers a guideline Heisenberg model on the honeycomb structure for the
for the materials design. purpose of bridging to the realistic and ab initio Hamil-
tonian. The Hamiltonian is given as
II. MODEL AND METHOD Hˆ = S~ˆTJ (λ)S~ˆ , (5)
λ i Γ j
Γ=X,Y,ZX,Z2nd,C3rdhiX,ji∈Γ
A. Effective Hamiltonian
where the matrices of the exchange couplings for the
three different nearest-neighbor (X-, Y-, and Z-bonds),
In this paper, we examine magnetic and thermal exci-
thesecondneighbor(Z -bond),andthethirdneighbor
tations of the generalized Kitaev-Heisenberg model pro- 2nd
posedinRef.9andtheab initioHamiltonianofNa IrO (C3rd-bond) bonds are given as
2 3
proposed in Ref. 14.
K′ 0 0 0 I′′ I′
2 2
JX(λ)= 0 0 0 +λI2′′ J′′ I1′ , (6)
0 0 0 I′ I′ J′
1. Generalized Kitaev-Heisenberg model    2 1 
ThegeneralizedKitaev-Heisenbergmodelisoneofthe 0 0 0 J′′ I′′ I′
2 1
simplest models that describe both Kitaev’s QSL and JY(λ)= 0 K′ 0+λ I2′′ 0 I2′ , (7)
zigzag magnetic orders. The model is parameterized by 0 0 0 I′ I′ J′
   1 2 
twoexchangecouplings,namelytheKitaev-typecoupling
K = 2Asinϕ and the Heisenberg coupling J = Acosϕ,
0 0 0 J I I
as 1 2
JZ(λ)=0 0 0 +λ I1 J I2 , (8)
Hˆ = S~ˆTJ S~ˆ , (1) 0 0 K I2 I2 0
CJK i Γ j    
Γ=XX,Y,ZhiX,ji∈Γ
where A has the dimension of energy and S~ˆ is an SU(2) J(2nd) I1(2nd) I2(2nd)
spin operators (Sˆx,Sˆy,Sˆz) at the i-th site,iand the ma- JZ2nd(λ)=λ I1(2nd) J(2nd) I2(2nd) , (9)
i i i  I(2nd) I(2nd) K(2nd) 
trices of the exchange couplings for the three nearest-  2 2 
neighbor bonds, X-, Y-, and Z-bond (see Fig. 1(a)), are
defined as J(3rd) 0 0
K+J 0 0 JC3rd(λ)=λ 0 J(3rd) 0 , (10)
JX = 0 J 0 , (2)  0 0 J(3rd) 
0 0 J
  respectively. The details of the bond are illustrated in
Fig. 1 (a). Here a parameter λ is introduced to inter-
J 0 0 polate the ab initio Hamiltonian for Na2IrO3 at λ = 1
JY = 0 K+J 0 , (3) and the Kitaev model at λ = 0. The ab initio estimates
0 0 J of the exchange couplings are summarized in Table II.
  This model at λ = 1 well explains not only the thermal
and properties, such as the specific heat for 5K < T < 40K
and the susceptibility for 5K < T < 400K, but also the
J 0 0
low-lying magnetic excitations14,16.
JZ = 0 J 0 . (4)
WhilethegroundstateoftheinterpolatedHamiltonian
0 0 K+J
  is in the gapless QSL phase at λ = 0, the zigzag-type
The ground states of the generalized Kitaev- antiferromagneticorder is stabilized at λ=114. Then, a
HeisenbergmodelrangefromtheKitaev’sQSLstotrivial quantumphasetransitionfromthetopologicalQSLstate
magneticallyorderedstatesdependingonthecontrolpa- to the magnetic orderedstate has to occur at least once,
rameter ϕ for fixed A taken positive. For ϕ ∼ 90◦ and between λ=0 and 1.
ϕ ∼ 270◦, the Kitaev’s QSLs appear. When the system
size N = 24 is considered, the stripy, N´eel, zigzag, and
ferromagnetic orders appear for −76,1◦ . ϕ . −33.8◦, B. Specific heat
−33.8◦ . ϕ . 87.7◦, 92.2◦ . ϕ . 161.8◦, and 161.8◦ .
ϕ . 251.8◦, respectively, determined in Ref.9 by exact The specific heat of the spin Hamiltonians, Hˆ and
CJK
diagonalization with N =24. Hˆ , is calculated by using exact energy spectra up to
λ
4
TABLEI: Threecategoriescharacterizedbythermalandmagneticexcitations;peakstructuresinthetemperaturedependences
of the specific heat C and nature of quasiparticles (QP). The nature of QP can be discussed from the comparison between
dynamical spin structurefactors S(Q,ω) and linear spin wave theory (SW).
C LRO QP S(Q,ω) vs. SW
I. single peak magnetic free magnon consistent
II. two peaks magnetic correlated magnon discrepant
III. two peaks no Majorana inconsistent
FIG. 1: (Color online) (a) Honeycomb structure model. Red, blue, and green lines denote Z-, Y-, and X- bond, respectively.
Dotted(dashed)linesrepresentthesecond(third)neighborbonds. Forthesecondneighborbonds,onlythebondperpendicular
toZ-bond,namely,Z -bondisshownbecausetheamplitudeforthesecondneighborinteractionsintheabinitioHamiltonian
2nd
is quite small and can be ignored14. Here, C denotes the set of the third neighbor bonds. (b-e) N=12, 16, 24, and 32-site
3rd
clusters, respectively.
theLanczosstepswithaHamiltonianHˆ,theTPQstates
TABLE II: Exchange couplings of the ab initio effective
at lower temperatures are constructed as follows: Start-
Hamiltonian for Na IrO derived in Ref. 14.
2 3 ing with an initial vector |Φ i = |φ i, the k-th step
0 +∞
JZ (meV) K J I1 I2 Lanczos vector |Φki (k ≥1) is constructed as
-30.7 4.4 -0.4 1.1
JX,Y (meV) K′ J′ J′′ I1′ I2′ I2′′ Hˆ |Φk−1i
|Φ i= . (12)
-23.9 2.0 3.2 1.8 -8.4 -3.1 k
JZ2nd (meV) K(2nd) J(2nd) I1(2nd) I2(2nd) qhΦk−1|Hˆ2|Φk−1i
-1.2 -0.8 1.0 -1.4
The above k-th step Lanczos vector is a TPQ state at a
JC3rd (meV) J(3rd) finitetemperatureT. Thecorrespondinginversetemper-
1.7
ature β = (k T)−1 is determined through the following
B
formula20,
N = 16 sites, and is estimated by employing thermal 2k k
B
pure quantum states20,21 for the 24- and 32-site clusters β = Λ−hΦ |Hˆ |Φ i +O(1/N), (13)
with the periodic boundary condition. The finite size k k
clusters used in the following are illustrated in Fig.1(b)- where k is the Boltzmannconstantand Λ is a constant
B
(e). larger than maxima of hHˆi. In other word, a TPQ state
Herewebrieflysummarizetheconstructionofthermal at T is given as,
pure quantum (TPQ) states following Ref. 20. A TPQ
stateatinfinitetemperaturesissimplygivenbyarandom |φTi=|Φki. (14)
vector,
The specific heatandentropyofHˆ arethenestimated
2N−1 byusingTPQstates|φ i. Thethermodynamicsandsta-
T
|φ+∞i= ci|ii, (11) tisticalmechanicstellusseveralprescriptionstocalculate
Xi=0 the specific heat and the entropy. Here, we calculate the
specific heat C by using the derivativeof internalenergy
where |ii is represented by the real-space S = 1/2 basis
with respect to the temperature as
and specified by a binary representation of decimal and
{ci} is a set of random complex numbers with the nor- dhφ |Hˆ |φ i
malizationcondition 2N−1|c |2 =1. Then,byutilizing C = T T , (15)
i=0 i dT
P
5
which is empirically known to be intruded by less statis- In this paper, we focus on the sum of diagonal ele-
tical errors in comparison with results obtained through ments, namely S(Q,ω)= Sµµ(Q,ω). In the gen-
thermal fluctuations of Hˆ. In the present paper, the en- µ=Xx,y,z
tropy S is estimated by integrating C/T from high tem- eralized Kitaev-Heisenberg model, the symmetry of the
peratures as model Hamiltonian ensures that the off-diagonal com-
ponent of the DSF becomes exactly zero. However, in
+∞ C general, Sµν(Q,ω) may have non-zero off-diagonal ele-
S =Nk ln2− dT′ , (16)
B Z T′ ments, if the off-diagonal elements of Jˆµν are non-zero.
T Γp
The contribution from such off-diagonal element is pro-
where C ∝ T−2 is assumed in the above integral for the portional to the Fourier transform of the correspond-
high temperature asymptotic behavior of C. Here, we ingtime-displacedspincorrelation,suchas Sˆx(t)Sˆy(0)
note that, for the specific heat and entropy defined in i j
D E
Eq.(15) and Eq.(16), respectively, of the lattice models, and Sˆx(t)Sˆz(0) . Sincetheamplitudeofthespincorre-
i j
it is convenient to use Nk , instead of the gas constant D E
B lationisscaledbytheamplitudeofthematrixelementof
R used in experiments. Jˆµν, the diagonal elements are dominant. The diagonal
Γp
elementsoftheab initioHamiltonianisindeeddominant
over the off-diagonal elements, and S(Q,ω) is expected
C. Equal-time spin correlation
to contain the main contribution of the spin excitations.
In comparisonwith the peak structures of the specific
heat, we examine temperature dependence of the equal-
III. THERMAL AND SPIN EXCITATIONS
time spin correlations. For short-rangespin correlations,
we calculate expectation values of spin operators SˆµSˆµ
i j A. Results of generalized Kitaev-Heisenberg model
atafinitetemperatureT withthethermalpurequantum
states20,21 |φ i as
T
To gain insights into the nature and proximity of the
SˆµSˆµ ≡hφ |SˆµSˆµ|φ i, (17) Kitaev’s QSL, we compare the specific heat, the linear
D i jET T i j T spin wave dispersion, and the dynamical spin structure
factor (DSF), which allows a classification of the spin
forthenearest-neighborpairshi,ji. Long-rangespincor- Hamiltonians Hˆ to three distinct categories as sum-
CJK
relations are characterized by the peak value in the mo-
marized in Table I. For several choices of the parameter
mentumdependenceoftheequal-timespinstructurefac-
ϕ for the generalized Kitaev-Heisenberg model, we clas-
torST(q)definedbyFouriertransformationof SˆiµSˆjµ , sify the nature of the spin and thermal excitations. The
as D ET resultant categorizationis summarized in Table III.
The specific heat C for the generalized Kitaev-
N−1
S (q)= 1 SˆµSˆµ cos(q·R ), (18) Heisenberg model is shown in Fig. 2. Here, we show the
T N µ=Xx,y,z Xℓ=0 D 0 ℓET ℓ rQeSsuLlstsaftorϕN==901◦2,,w16e,s2e4e,naondst3r2o.ngExNcedpetptehnedKenitcaeevfo’sr
at q =Q where R is the position vector of the ℓ-th site the second-largest and largest system sizes N = 24 and
ℓ
and Q is the momentum at the maximum. 32. For the trivially ordered states at ϕ = 0◦ and 180◦,
the temperaturedependences ofC showsaSchottky-like
singlepeak within the errorbars. Inthe thermodynamic
limit, the peak may evolve into the anomaly (divergence
D. Dynamical spin structure factors
or peak) expected by the growth of spin correlations ac-
companied by the transition to the long-range order. In
Todiscussmagneticexcitationsbyaspinflip,wefocus contrast, for the Kitaev’s QSL states at ϕ = 90◦ ( and
on the dynamical spin structure factor (DSF) at zero equivalently at ϕ = 270◦), there are two peaks in the
temperature. The DSF is defined as
temperature dependences of C, which is a hallmark of
thermal fractionalization proposed in Ref. 6, as is intro-
1 1
Sµν(Q,ω)≡− lim Imhφ |Sˆµ† Sˆν|φ i, duced in Sec. I. The low-temperature peak of C may be
π ǫ→+0 0 Q ω+E0+iǫ−H Q 0 associated with the contribution from the thermal flux
(19)excitations (or the thermal excitations of the localized
Majorana fermions), while the high-temperature peak
whereφ isthe groundstate ofHwith the groundstate
0 may represent the excitation of the itinerant Majorana
energyE0. ThespinoperatorSˆQµ istheFouriertransform fermions. The category II represented by ϕ = 100◦ and
of Sˆµ, where µ,ν stand for x,y or z component. After 240◦ is the same as the category I as to the presence
i
calculatingφ andE bytheLanczosmethod,Sµν(Q,ω) of the magnetic order, while the two-peak structure of
0 0
is obtained by the continued fraction expansion22,23. C exists similarly to the category III. The ordered state
6
FIG.2: (Coloronline)TemperaturedependencesofthespecificheatC ofthegeneralized Kitaev-HeisenbergmodelHˆ . The
CJK
results for N =12 and N = 16 obtained by fully diagonalizing Hˆ (denoted as “full ED”) are illustrated with thin broken
CJK
andthinsolid(blue)curves,respectively. ForN =24andN =32(thickbrokenandthicksolid(red)curves),thethermalpure
quantum(TPQ) states20 areemployed. Thepossible errors of TPQ duetothetruncationof theHilbert space areshown in C
byshaded (gray) belts, which isestimated by usingthestandard deviation of theresults obtained from 4to 36 initial random
wave functions at the high temperature limit, kBT/A→ +∞. From the top leftmost to top rightmost panels, C/N is shown
for ϕ=0◦, 90◦, 100◦, and 120◦ in this order. The same quantities are shown for ϕ=180◦, 210◦, 240◦, and 300◦ in this order
from the bottom leftmost to bottom rightmost.
QSL, we will propose later that the entropy at temper-
atures between the two peaks in C serves as a hallmark
of the closeness to the Kitaev’s QSL. In Fig. 3, we show
theschematicillustrationforthecategorizationobtained
from the temperature dependence of C and the presence
of the magnetic order.
In Fig. 4, we compare the temperature dependence of
thelong-rangepartofthespincorrelationrepresentedby
the peak in S (Q) and the short-range part represented
T
by SˆµSˆµ . While, in the category I, the short-range
i j
D ET
spin correlations SˆµSˆµ and the long-range spin cor-
i j
D ET
relations S (Q) at the ordering wave vector Q grow si-
T
multaneously around the temperature where C has the
singlepeak, SˆµSˆµ andS (Q)growindependentlyin
i j T
D ET
the category II represented by ϕ = 100◦ and 240◦. The
FIG.3: (Coloronline)Categorizationofgroundstatesofgen-
eralizedKitaev-Heisenbergmodel. Phaseboundariesdepicted growthofthe spin correlationchangesfromthe category
by the solid lines are drawn by using the results in Ref. 9. I to the categoryII aroundϕ=300◦. In the categoryII,
SoliddotsrepresenttheparametersshowninFigs. 2,4,and5. SˆµSˆµ growsas temperature falls to T , which corre-
i j h
Dashed lines separate whether or not the two-peak structure D ET
spondstothehigh-temperaturepeakinC,andsaturates
inthetemperaturedependenceofthespecificheatisobserved
withinthemagneticorderedphase,namelythecrossoverbor- below Th. On the other hand, the long-rangespin corre-
der between the categories I (blue (lightly shaded) area) and lations represented by ST(Q) grow significantly around
II(red(darkshaded)area). Notethatthedottedlinesdonot the temperatures Tℓ where C has the low-temperature
represent thephase boundary. peak, in the category II. In the category III, the short-
range spin correlation grows at T , while the long-range
h
part does not show appreciable temperature dependence
evenatT incontrasttothecategoryII.Thelowtemper-
ℓ
at ϕ = 300◦ is located almost on the border between aturepeakinC inthecategoryIIIarisesfromanentirely
the category I and category II. Although the two-peak different mechanism from that in the category II, as we
structure itself is not a unique feature of the Kitaev’s already mentioned. The difference between the category
7
FIG. 4: (Color online) Temperature dependence of short-range spin correlation SˆxSˆx for X-bond and long-range spin
D i jE
T
correlation ST(Q) in comparison with the specific heat C. The results for N = 24 by employing the thermal pure quantum
states20 areshown. Thepossible errorsofTPQduetothetruncationoftheHilbertspaceareshowninDSˆixSˆjxE ,ST(Q),and
T
C byshaded(gray)belts,whichisestimatedbyusingthestandarddeviationoftheresultsobtainedfrom4to36initialrandom
wavefunctionsatthehightemperaturelimit, kBT/A→+∞. ForST(Q),following pointsinFig.5(i)areselected asthewave
vector Q; the Γ, Y, Γ∗, and X points are selected since ST(Q) for each value of ϕ has maxima at the momentum. These
points represent the Bragg points of the ferromagnetic (FM), zigzag, N´eel, and stripy orders, respectively. The momentum Q
isconsistentwiththephasediagram inFig.3. Fromthetopleftmost totoprightmostpanels,theresultsareshownforϕ=0◦,
90◦, 100◦, and 120◦ in this order. The same quantities are shown for ϕ = 180◦, 210◦, 240◦, and 300◦ in this order from the
bottom leftmost to bottom rightmost.
II and the category III is evident in the temperature de- magnetic/ferromagnetic Heisenberg model. Therefore,
pendence of the peak of S (Q) in comparison with that the low-lying excitation of the DSF is expected to be
T
of SˆµSˆµ , as shown in Fig.4. well explained by the spin wave mode. At ϕ = 180◦,
i j
D ET all poles in the DSF are perfectly located on the spin
Basedontheaboveresults,wecategorizethequantum wave mode. In the antiferromagnetic case at ϕ = 0◦,
phases obtained for the generalized Kitaev-Heisenberg
the low-lying excitation agrees with the linear spin wave
model. The summary is shown in Table III and Fig. 3.
mode by introducing a renormalizationfactor a. a is es-
timated from the best fitting of spin wave dispersions
for the poles of the low-lying excitations in the DSF;
TABLE III: Categorization of ground states for several
S(Q,ω ) ∼ a × ω (Q), where S(Q,ω ) de-
choicesofϕofthegeneralizedKitaev-Heisenbergmodel. Def- lowest LSW lowest
inition of thecategories I – III is shown in Table I. notes the poles of the lowest excitation in the DSF and
ω (Q) is the linear spin wave mode. At ϕ = 0◦, we
ϕ quantumphase category LSW
◦ obtain a∼1.3, which is the upper limit because the po-
0 N´eel I.
90◦ Kitaev’s QSL III. sitions ofpoles areaffectedby the systemsize. (The size
100◦ zigzag II. dependence becomes strong especially at the symmetric
120◦ zigzag I. wave-number points. ) As well studied in Refs. [24,25],
180◦ ferromagnetic I. the exact low-lying excitations in the Heisenberg models
◦
210 ferromagnetic I. aredescribedby the linearspinwavemodewithanO(1)
240◦ ferromagnetic II. renormalization factor a. The renormalization factor is
270◦ Kitaev’s QSL III. a = π in the S=1/2 spin chain case24. This value can
300◦ stripy I/II. 2
be an indicator for the renormalization of the quantum
fluctuation. As shown below, we observe the two-peak
By comparing the low-lying excitations of the dynam- structure in the temperature dependence of the specific
ical spin structure factors (DSFs) with the linear spin heat C, when a'1.5.
wave approximation, we further confirm that the cate- In contrast, magnetic excitations described by the
gorization is robust. In Fig.5, we show both results of poles in the DSF are completely different from coher-
the DSF S(Q,ω) for N =24 and of the linear spin wave ent magnetic excitations in the Kitaev’s QSL phase at
calculations. ϕ = 90◦(270◦). This is due to diverging quantum fluc-
We start from the results for the simplest case. At tuations in the Kitaev’s QSL phase and massive degen-
ϕ = 0◦ and 180◦, the model becomes the antiferro- eracy of classical spin orders spoils the linear spin wave
8
FIG. 5: (Color online) Dynamical spin structure factors of generalized Kitaev-Heisenberg models for 24-site cluster. Area of
eachcircleisproportionaltotheintensityinlogarithmicscale. Abscissarepresentslabeledpointsin(i). Fromthetopleftmost
to top rightmost panels, S(Q,ω) is shown for (a) ϕ =0◦, (b) 90◦, (c) 100◦, and (d) 120◦ in this order. S(Q,ω) is shown for
◦ ◦ ◦ ◦
(e) ϕ = 180 , (f) 210 , (g) 240 , and (h) 300 in the order from the bottom leftmost to bottom rightmost. Here, the results
of the linear spin wave theory are shown in solid curves. In order to compare the spin wave dispersions with the low-lying
excitations, each spin wave result in the top panels is multiplied by a constant. The multiplication constants for (a), (c) and
(d) are 1.3, 1.65, and 1.45, respectively. The spin wave results in the lower panels are drawn without such tuning constant.
For the results in thezigzag and stripy phase, (c), (d), and (h), we also plot the spin wave dispersions obtained from the case
where theordered state is rotated by 2π/3.
analysis. Below, we explain the low-lying excitations of modeexcepttheM/Γ(Γ∗)pointcanexplainthelow-lying
theDSFs,whenthesystemapproachestheKitaev’sQSL excitations of the DSF without the renormalization fac-
phasefromthedeepinsideofthemagneticorderedphase. tor a discussed above.
First, we focus on the positive Kitaev-coupling case Here,we detailthe discrepancybetweenthe spinwave
for 0≤ ϕ≤180◦. For 0 ≤ϕ≤180◦ in the magnetic or- modeandthe low-lyingexcitationoftheDSFforK <0,
deredphases,the low-lyingexcitationofthe DSF canbe which is another clue to categorizethe magnetic ordered
explained by the correlated/renormalized magnon exci- phasesintothecategoriesIandIIinadditiontothetem-
tation. The fitting parametera is alwayslargerthan the perature dependence of specific heat. When the system
unity; the low-lying excitation becomes hard in compar- approaches the Kitaev’s QSL phase around ϕ = 270◦,
ison with the linear spin wave mode. When the system the categorization of the categories I and II can be dis-
approaches the Kitaev’s QSL phase around 90◦, a dras- cussedfromthediscrepancybetweenthespinwavemode
tically increases. Though not shown in the figure, we and the low-lying excitation of the DSF at the M point.
obtaina∼1.85atϕ=86◦. At ϕ=100◦ in the category At ϕ = 210◦ in the category I, the low-lying excitation
II, the ground state is the zigzag ordered state. From is well explained by the spin wave mode. However, the
Fig. 2,weconfirmthe two-peakstructureinthetemper- discrepancybetween the low-lyingexcitation ofthe DSF
ature dependence of C. The renormalization factor a at and the spin wave mode develops clearly at the M point
ϕ=100◦ isestimatedasa∼1.65andislargerthanthat atϕ=240◦ in the categoryII. The low-lyingexcitations
of the S=1/2 chaincase. When the system goes into the of the DSF are located in the lower energy region than
deep inside ofthe magnetic phase,the factor a decreases that of the spin wave mode. (The discrepancy between
and crosses a ∼ 1.5, where the two-peak structure is al- thebothisalsoclearattheYpoint. However,theexcita-
most smeared out in C. At ϕ = 120◦ in the category tion ofthe DSF atthe Y point is identical to that at the
I, a is about 1.45. While the ground state is still in the M point. This is due to to the finite size effect; the sym-
zigzag ordered phase, the temperature dependence of C metry breaking is prohibited in the finite size systems.)
shows the usual single peak. We regard ϕ=210◦ is located in I near in the crossover
Next, we see the results for the negative Kitaev cou- region of the categories I and II, where both characters
pling case for 180◦ ≤ ϕ ≤ 360◦. In contrast to the posi- are mixed. Based on these observations, here we catego-
tive Kitaev coupling case, the low-lyingexcitationin the rize ϕ = 210◦ and ϕ = 240◦ as the category I and cat-
DSF can be well explained by free magnon picture for egory II, respectively. In the stripy phase for ϕ > 270◦,
K < 0 in the magnetic ordered phases: The spin wave thewavevector,wherethespinwavemodedeviatesfrom
9
the low-lyingexcitationsofthe DSF, movestothe Γ and thermodynamic limit, defined as
Γ∗ points. At ϕ = 300◦, the low-lying excitation of the
DSF shows slight softening at the Γ and Γ∗ points in m(Q)≡ lim ST(Q), (20)
T→+0
comparison with the spin wave mode. In addition, the p
temperaturedependenceofC showsasinglepeakwitha where we take the T =0 limit.
prominent shoulder, which is in between the single peak The quantum phase transitions are expected to cause
structureofthe categoryIandthe two-peakstructureof divergenceordiscontinuity(orasharppeakinfinite-size
the categoryII. Consequently, the system at ϕ=300◦ is systems) in the second derivatives of the ground state
concluded to be located in the crossover region between energy E with respect to λ,
the category I and the category II.
d2E/N
α≡− . (21)
dλ2
B. Results for ab initio Hamiltonian of Na IrO
2 3
In Fig.6, the λ-dependences of m(Q) and α for N = 24
site cluster are given.
Fromthegrowthofm(Q)attheYpointcorresponding
to the zigzag order observed in the experiments and a
peak in α atλ∼0.59,we conclude that the zigzag order
Y * appears for λ&0.6. For λ/0.4, we expect the Kitaev’s
0.3 M
QSL ground states. Since the phase transition around
X zigzag (Z) λ ∼ 0.6 appears to be continuous with the reduction
of m toward the transition point λ ∼ 0.6, the distance
0.25
fromthe Kitaev’s QSL phase may be measuredfromthe
ordered moment.
0.2 The presence of the phase boundary to the zigzag or-
deredphaseisalsoconfirmedfromtheDSFresultsshown
in Fig. 7. At λ = 0, we observe a characteristic non-
0.15 dispersivemodeatω/A∼0.3reflectingtheKitaev’sQSL
0.4 ground state5. For λ < 0.6, some poles with the weak
0.3 intensity appear below ω/A∼0.3 and the peak with the
0.2 largest intensity is not located in the lowest excitation
0.1 mode. Thesepropertiesofthelow-energyexcitationsare
0 contrast to those in the magnetic ordered phase, where
0 0.2 0.4 0.6 0.8 1
the lowest excitation usually shows the largest intensity.
We also observe the lack of well developed peaks in the
equal-time spin structure factors shown in Fig. 6. Thus,
FIG. 6: (Color online) Square root of the equal-time spin
structure factors, m(Q) and the second derivative of ground we conclude that the ground state at λ< 0.6 is still the
state energy α as functions of λ calculated with theN =24- Kitaev’s QSL state.
site cluster in Fig.1(d). For 0.6 < λ ≤ 1.0, the lowest excitation with the
largest intensity appears at the M or Y point and the
peak at the Y point develops as λ increases. In this re-
In the light of the categorizationexamined in the gen- gion, the equal-time spin structure factors at the M and
eralized Kitaev-Heisenberg models, we examine the in- Y point well develops as shown in Fig. 6. Therefore, the
terpolated Hamiltonian between the Kitaev limit, Hˆ , ground state becomes magnetically ordered for λ > 0.6
λ=0
andtheabinitioHamiltonianofNa IrO ,Hˆ ,givenin and the presence of the phase boundary is expected for
2 3 λ=1
Eq.(5). First, we show that the peak value of the equal- λ∼0.6.
time spin structure factor for the zigzag order starts To confirm that the interpolated Hamiltonian Hˆ is
λ
growingatan onsetvalue λ . Here λ is around0.6. For categorized into the category II for λ > 0.6, we also ex-
c c
all 0≤λ≤1, the dynamical spin structure factor of Hˆ amine the temperature dependences of the specific heat
λ
is not captured by the linear spin wave theory16. Next, C. The results are shown in Fig. 8. First of all, for the
weshowthatthe temperaturedependence ofthe specific entire parameterrange,0≤λ≤1, two peaksare seenin
heatfor Hˆ alwayshas a two-peakstructure irrespective thetemperaturedependencesofC. Thelow-temperature
λ
of λ. peak is at T/|K| . 0.03, and the higher (high-T) one is
The quantum phases of the interpolated Hamiltonian at T/|K| ∼ 0.3, where |K|=30.7 meV (356 K) (see Ta-
areexaminedbyequal-timespinstructurefactorsandthe ble II). As already discussed in Ref. 5 and Ref. 6, the
secondderivativesofthe groundstateenergydefinedbe- originofthe high-T peakiswellexplainedbythe growth
low; herewe introduce the square rootof the normalized of magnetic correlations for the nearest neighbor pairs,
equal-time spin structure factor, which is extrapolated which is determined by the Kitaev couplings, K ∼ K′.
to the magnetic orderparameteratmomentum Qin the Thelow-T peakinthepureKitaevmodelcorrespondsto
10
1
ω/A
10-1
λ=0.0 λ=0.2 λ=0.4 λ=0.6 λ=0.8 λ=1.0
10-2 X Γ Y Γ* M X Γ Y Γ* M X Γ Y Γ* M X Γ Y Γ* M X Γ Y Γ* M X Γ Y Γ* M
Q Q Q Q Q Q
FIG. 7: (Color online) Dynamicalspin structurefactors of theinterpolated Hamiltonian Hˆ between theKitaev limit andthe
λ
ab initio spin Hamiltonian proposed in Ref. 14. From the leftmost to rightmost panels, S(Q,ω) is shown for λ = 0, λ = 0.2,
λ=0.4, λ=0.6, λ=0.8, and λ=1. Abscissa representsthe labeled points in Fig.5(i).
FIG.8: (Coloronline)TemperaturedependencesofthespecificheatC oftheinterpolatedHamiltonianHˆ betweentheKitaev
λ
limitandtheabinitiospinHamiltonianproposedinRef.14. TheresultsforN =12andN =16obtainedbyfullydiagonalizing
Hˆ (denoted as “full ED”) are illustrated with thin broken and thin solid (blue) curves,respectively. ForN =24 and N =32
λ
(thick broken and thick solid (red) curves), the thermal pure quantum (TPQ) states20 are employed. The possible errors of
TPQduetothetruncationoftheHilbertspaceareshowninC byshaded(gray)belts,whichisestimatedbyusingthestandard
deviation of the results obtained from 4 to 36 initial random wave functions.
thethermalfluctuationofthe localZ gaugefieldthatis quantitative measure to estimate the distance between
2
one of two Majorana fermions yielded via the fractional- the real material and the Kitaev’s QSL.
ization of an originalquantum spin. The non-monotonic
As is shown in Figs. 2 and 8, the specific heat C in
λ-dependencesinthelow-T peakcorrespondtothequan-
the category II shows two peaks similar to the cate-
tum phase transition around λ∼0.6 and possible emer-
gory III. However, the two-peak structure in the tem-
gence of the intermediate phase for 0.4.λ.0.6.
perature dependence of C itself is not a unique feature
in the vicinity of the Kitaev’s QSL. For example, ge-
ometrically frustrated and quasi-one-dimensional quan-
IV. DISTANCE FROM THE KITAEV SPIN tum spin systems also show the two-peak structure in
LIQUID PHASE the temperature dependence of C26. In that case, the
high- and low-temperature peaks depend on the dimen-
In the previous section, Sec. III, magnetically ordered sional anisotropy, where the low-temperature peak in C
materials categorized as the category II are expected to representsthe entropy release arising fromthe reallong-
beclosetotheKitaev’sspinliquidintheparameterspace range order. If the anisotropy increases, the release of
of the effective Hamiltonians. We further propose more the entropy at the low temperature peak may decrease.