Ground State of the Spin-1/2 Kagome Lattice Heisenberg Antiferromagnet Rajiv R. P. Singh Department of Physics, University of California, Davis, CA 95616, USA David A. Huse Department of Physics, Princeton University, Princeton, NJ 08544, USA (Dated: February 1, 2008) Using series expansions around the dimer limit, we find that the ground state of the spin-1/2 Heisenberg Antiferromagnet on the Kagome Lattice appears to be a Valence Bond Crystal (VBC) with a 36 site unit cell, and ground state energy per site E = −0.433±0.001J. It consists of 8 a honeycomb lattice of ‘perfect hexagons’. The energy difference between the ground state and 0 other ordered states with the maximum number of ‘perfect hexagons’, such as a stripe-ordered 0 state, is of order0.001J. The expansion is also donefor the36 site system with periodic boundary 2 conditions;itsenergypersiteis0.005±0.001J lowerthantheinfinitesystem,consistentwithExact n Diagonalization results. Every unit cell of the VBC has two singlet states whose degeneracy is not a liftedto6thorderintheexpansion. Weestimatethisenergydifferencetobelessthan0.001J. The J dimerizationorderparameterisfoundtoberobust. Twoleadingordersofperturbationtheorygive 8 lowest triplet excitations to bedispersionless and confined to the‘perfect hexagons’. ] PACSnumbers: 75.10.Jm l e - r The spin-1/2 antiferromagnetic Kagome-Lattice t s Heisenberg Model (KLHM) with Hamiltonian, . E E E E t H H a m =J XSi Sj, (1) H · E E hi,ji - d n is a highly frustrated quantum spin model [1]. Its prop- E E E o erties have been studied by a wide variety of numerical H P H c and analytical techniques [2, 3, 4, 5, 6, 7, 8, 9]. Yet, [ E E the precise nature of the ground state remains a subject 3 of debate. Proposals have included a number of Valence v Bond Crystals (VBC) [10, 11, 12, 13] as well as spin- E E E E H H 2 liquid states with algebraic correlations [14, 15]. Recent 9 experimental work on the material ZnCu3(OH)6Cl2has E E 8 attracted further interest to this model [16, 17, 18, 19], 0 although this material is likely to also have significant . E E E 7 Dzyloshinski-Moria anisotropy. [20] H P H 0 Here,weshowthatthegroundstateofKLHMappears 7 to be a Valence Bond Crystal with a 36 site unit cell. It E E 0 consists of a Honeycomb lattice of perfect hexagons as : v initially proposed by Marston and Zeng [10], discussed i X in more detail by Nikolic and Senthil [12], and shown in Fig 1. In a dimer covering, all triangles that have a sin- FIG. 1: (Color online) Ground state ordering pattern of r a gletvalencebondarelocallyinagroundstate. Ascanbe low-energy (“strong”) bonds (blue) for the Kagome Lattice readilyshown,anydimercoveringleavesone-fourthofall Heisenberg Model. The perfect hexagons are denoted as H, the empty triangles by E and the pinwheels as P. The two trianglesintheKagomelatticeempty. Allquantumfluc- dimer coverings of the pinwheels that remain degenerate to tuations in the ground state originate from these empty highordersof perturbationtheoryare denotedbythicksolid triangles,sinceitisonlytherethatthesingletdimercov- (blue) and dotted (magneta) bonds. ering is not locally a ground state of the Hamiltonian. We develop series expansions around an arbitrary dimer covering of the infinite lattice using a Linked Cluster method [21] and compare the energies of vari- all bonds are equivalent in the Hamiltonian. Follow- ous dimer coverings. To carry out the expansions, all ing recent development of the Numerical Linked Cluster (“strong”) bonds that make up the dimer covering are scheme[22], we grouptogetherallweakbonds belonging given an interaction strength J and all other (“weak”) to each triangle. This significantly simplifies the calcu- bonds aregivena strengthλJ. Expansionsarethen car- lations: only 5 graphs contribute to the ground state ried out in powers of λ and extrapolated to λ=1 where energy to 5th order in λ (see Fig 2). The resulting series 2 expansionforthe groundstateenergyshowssurprisingly mentintothisVBCalsoleavesa“pinwheel”ofdimersin rapid convergence even at λ=1. everyunitcell. Therearetwodegeneratepinwheeldimer coverings,and this degeneracy is not lifted to quite high orders of perturbation theory, implying that in a lattice of N sites there are at least 2N/36 states whose energy difference, per site, is much less than 0.001 J. We have carried the expansion to 5th order for the (a) (b) ground state energy. For the infinite lattice, the contri- bution at 5th order is quite small (about 0.0005 J per site) andslightly increasesthe energydifference between thehoneycombandstripeVBC[12]thatwasestablished at 4th order. (d) In addition to carrying out the expansions for the in- (c) finite system, one can also study the ground state prop- erties of finite systems using the series expansion tech- niques. The 36 site cluster with periodic boundary con- ditions (PBC) studied by Leung and Elser [4] and Le- cheminant et al. [6], does accommodate the unit cell of theHoneycombVBCbutnotthatofthestripeVBC.For the honeycomb VBC dimer covering on this cluster, the (e) energyagreeswiththeinfinitesystemthrough3rdorder. However, at 4th order there are additional graphs that enter due to the periodic boundary conditions, permit- ting 4 dimers to resonate along a path through 4 empty FIG. 2: Topologies of the graphs that contribute to the trianglesthatwindsoncearoundtheperiodicboundaries. ground state energy to 5th order of perturbation theory. From comparing the series through 5th order and using Graphs (a), (b), (c), (d), and (e) begin contributing in or- simple Pade estimates, we find the ground state energy ders 2nd,3rd, 4th, 4th, and 5th, respectively. forthe36siteclustertobelowerthanthatoftheinfinite system by 0.005 0.001 J per site, consistent with the To 2nd order in our expansion for the ground state ± ExactDiagonalizationresultofanenergyof 0.438377J energy, all dimer coverings have the same energy. At − per site [4]. 0th order, the singlet dimers each have energy 0.75 J. − The λ = 1 energies for the infinite-lattice Honeycomb Thesedimerscanlowertheirenergyfurther by“resonat- VBC,the stripeVBC ofNikolicandSenthil[12],andfor ing”, but only on the empty triangles. At 2nd order, the finite 36 site cluster with PBC and honeycombVBC these pairwise dimer-dimer interactionslowerthe energy are listed in Table I, summed through 5th order. Note by 0.28125 J per empty triangle. The degeneracy of all the apparentlyrapidconvergenceofthe series,especially the dimer configurations is finally lifted at third order, for the infinite lattice VBC states. where there is binding of 3 empty triangles into “per- fect hexagons”; the 3 dimers around a perfect hexagon (see Fig. 2(b)) can then resonate together. This pro- TABLEI: Ground state energy persite in unitsof J for var- cess lowers the energy of each perfect hexagon at 3rd ious dimer states in perturbation theory. At each order we order by 0.0703125J. At 4th order, resonances between sum the contributionsthrough that order for λ=1. twoempty triangleswhenthey areconnectedbya single Order Honeycomb VBC StripeVBC 36-site PBC dimerbondasinFig. 2(c)produceanadditionalbinding energy of amount 0.01953125J. 0 -0.375 -0.375 -0.375 Thus at 4th order in our expansion we find that the 1 -0.375 -0.375 -0.375 ground state comes from the dimer covering with the 2 -0.421875 -0.421875 -0.421875 maximum number of perfect hexagons, with neighbor- 3 -0.42578125 -0.42578125 -0.42578125 ing hexagons arranged such that their empty triangles 4 -0.431559245 -0.43101671 -0.43400065 share a dimer bond as much as possible. It follows that 5 -0.432088216 -0.43153212 -0.43624539 once the positions and dimer coverings of two neighbor- ing perfect hexagonsarepicked, the restuniquely fall on a Honeycomb lattice of perfect hexagons. Furthermore, For the 36 site cluster, our considerations imply 48 the dimerizations of all hexagons are simple translations spin-singlet states with very low energies per site within of one another. This leads to the Valence Bond Crystal oforder0.001J ofthe groundstate,24correspondingto (VBC) arrangementshown in Fig. 1. For an infinite lat- the ground state degeneracy of the thermodynamic sys- tice, the groundstate ofthis arrayofperfecthexagonsis temand2eachforthepin-wheelstates. Buttheother48 24-fold degenerate, with 12 distinct translations, each of states with two ‘perfect hexagons’ should also fall below which has 2 reflections. As shownin Fig. 1, the arrange- the lowest triplet state, whose excitation energy for the 3 36 site cluster translates to a per site value larger than sponding to the 18 dimers. In 0th order, these triplets 0.004 J. In the Exact Diagonalizationstudies about 180 are dispersionless and have excitation energy J. To singlet states are found below the lowest triplet. Since second order in perturbation theory, only 9 of the 18 thebindingenergyforaperfecthexagontranslatesintoa triplets,consistingofthesixbelongingtothetwoperfect persiteenergyof0.002J forthe36sitecluster,somesin- hexagons,andthethreethatcoupletheperfecthexagons gletstateswithonlyoneperfecthexagonpresumablycan via empty triangles, become dispersive and form a net- also have energies below the triplet gap in that cluster. workonwhichexcitationscanmove. Theother9triplets, We have also studied the local bond energies to 3rd consistingofthe6pinwheeldimersandthe3otherdimers order in perturbation theory. To this order only four that do not touch empty triangles, can only perform types of nearest neighbor bonds, representing only 27 of virtual hops at this order and thus remain dispersion- the 72 bonds per unit cell, have their expectation values less. The effective Hamiltonian for the nine dispersive changed from the 0th order values of 0.75 J for strong states can be expressed in terms of a 9 9 matrix. Let bonds and 0 for weak bonds. These −are: bond type A, z1 = expi~k ~r1 with ~r1 = 4√3yˆ, and×z2 = expi~k ~r2 the strong bonds inside the perfect hexagons; bond type with ~r2 = ·6xˆ 2√3yˆ, and−z1∗ and z2∗ be their comp·lex − − B, the weak bonds inside the perfect hexagons; bond conjugates, with the lattice oriented as in Fig. 1, with a type C, the weak-bonds in the empty triangles which do nearest-neighborspacing ofunity. Then, in secondorder not form part of the hexagon; and bond type D, the perturbation theory, the triplet excitation energies are strong bonds that join two empty triangles between per- the eigenvalues of the matrix fect hexagons. Their expectation values are listed in Ta- 5 1 1 1 pbleerIbIo.nTdh.eAtltohtaolugghrotuhnedsesrtiaetsefoerntehrgeytoistaalbenouertg−y0is.2c2onJ- ∆/J =(1− 8λ2)+(4λ+ 8λ2)M1+(32λ2)M2, (2) verging quickly to this value, the energies of individual where the matrices M1 and M2 as a function of z1 and bonds remain very different, as is expected for a valence z2 are given in Table III and IV. bond crystal. TABLE III: Matrix M1 TABLEII: ExpectationvalueshSi·Sjiforvariousbondtypes 0 -1 -1 -1 1 0 0 0 0 Bond 0th order 2nd order 3rd order -1 0 -1 1 0 -1 0 0 0 A -0.75 -0.5625 -0.515625 -1 -1 0 0 -1 1 0 0 0 B 0 -0.1875 -0.257812 -1 1 0 0 0 0 -1 0 1 C 0 -0.1875 -0.1875 1 0 -1 0 0 0 0 z1 −z1 D -0.75 -0.5625 -0.5625 0 -1 1 0 0 0 z2 −z2 0 0 0 0 -1 0 z2∗ 0 -1 -1 Leung and Elser[4] had also calculated the energy- 0 0 0 0 z1∗ −z2∗ -1 0 -1 energy correlations in the ground state of the 36 site 0 0 0 1 −z1∗ 0 -1 -1 0 cluster. They had noted that the correlations at largest distances qualitatively matched the Valence Bond Crys- Notethatthisimpliesthattwotripletshavethelowest talpatternafteraveragingoverthe 48degeneratestates. energy andthey are eachconfined to one of the two per- However,onaquantitative levelthere was no correspon- fecthexagonsandthusremaindispersionless. Thetriplet dence with the VBC, as the dimerization pattern seen gapbecomes1 (1/2)λ (7/8)λ2,whichaddsupto 3/8 in the Exact Diagonalization was considerably weaker. at λ = 1 if we−truncat−e at this order. This can be−un- The modified values of the bond energy expectation val- derstoodasarisingfromtwofactors: Thehexagonislike ues, calculated by us, do not significantly change this a one-dimensional Alternating Heisenberg Chain and in conclusion. We believe the reason for the discrepancy is that case the gap is known to be 1 (1/2)λ (3/8)λ2 the following: In the 36 site cluster, energy-energy cor- [23]. There is a slight difference here−with resp−ect to the relations do not factorize even at the largest distances. Alternating Heisenberg Chain because of the additional Due to the Periodic Boundary Conditions, even bonds neighbors. The diagonal term is doubled but the sec- which are furthest away on the 36 site cluster can get ond neighbor hop is absent, so that the overallresult for correlated already in 2nd order of perturbation theory. the gap is the same. In addition, each strong bond has Thusitisnotappropriatetosimplycomparecorrelations its two ends both connected to the same end of another (S~i S~j)(S~k S~l) tothe products S~i S~j S~k S~l . Inthe dimer to which a triplet on that strong bond may make h · · i h · ih · i future, we hope to calculate and compare the energy- a virtual hop. This is like the triplets in the Shastry- energy correlations for the thermodynamic system and Sutherlandmodel[24,25]. Due tothis,andthefactthat the 36 site cluster within this dimer expansion. no such process happens in the ground state, there is We have also calculated two leading orders of pertur- a dispersionless reduction of the spin gap of magnitude 2 bation theory for the triplet excitation spectra. There (1/2)λ . At 2ndorderthese two types ofcontributions − are 18 elementary triplet excitations per unit cell, corre- simply add. Thus we find that although the series for 4 the groundstateenergyshowsapparentlystrongconver- cificheatandentropyoftheKLHMshouldhavestructure gence, the series for the spin gap does not, although the down to T/J < 0.001. It was found in the high temper- latter statement is about only the first two orders in the ature expansion study [5] that a naive extrapolation of expansion. Further study ofthe triplet excitationsis left the hightemperatureseriesdowntoT =0ledtoa finite for future work. ground state entropy. Furthermore, several studies have Our results have important implications for the finite suggestedmultiple peaksinthespecific heat[9,26]. Our temperature properties of the model. First of all, since work shows at least 2N/36 very low energy states of the thedifferenceinenergyoftheHoneycombVBCstateand ‘pinwheels’, for an N-site system, and many other very other dimer states (such as stripes) is less than 0.001J low lying states, giving further support to these results. per site,anyfinite temperaturetransitioninto aHoney- comb VBC phase should occur at a temperature of this order or lower. Secondly, given the large number (24) of degenerateorderingpatterns,the phasetransitionseems Acknowledgments likely to be first order. Third, since the attraction be- tween the empty triangles is only of order 0.02J, above T 0.02J, all dimer configurations should have com- This work was supported by the US National Science ≈ parable Boltzmann weight, giving rise to a dimer liquid Foundation,GrantsNo.DMR-0240918(R.S.)andDMR- regime at intermediate temperatures. Fourth, the spe- 0213706(D. H.) and PHY05-51164. TABLE IV: Matrix M2 0 0 0 1 -1 0 -1 −z2 1+z2 0 0 0 -1 0 1 1+z1 −z1 -1 0 0 0 0 1 -1 −z1 z1+z2 −z2 1 -1 0 0 1+z2∗ 1+z1∗ 1 0 -1 -1 0 1 1+z2 0 1+z1∗z2 0 −z2 z2 0 1 -1 1+z1 1+z1z2∗ 0 −z1 z1 0 -1 1+z1∗ −z1∗ 1 0 −z1∗ 0 0 0 −z2∗ −z1∗ z1∗+z2∗ 0 −z2∗ z1∗ 0 0 0 1+z2∗ -1 −z2∗ -1 z2∗ 0 0 0 0 [1] For a recent review see C. Lhuillier, (2003). arXiv:cond-mat/0502464. [13] R.BudnikandA.Auerbach,Phys.Rev.Lett.93,187205 [2] C. Zeng and V. Elser, Phys.Rev.B 42, 8436 (1990). (2004). [3] R.R.P.SinghandD.A.Huse,Phys.Rev.Lett.68,1766 [14] M.Hermele,T.Senthil,andM.P.A.Fisher, Phys.Rev. (1992). B 72, 104404 (2005). [4] P.W.LeungandV.Elser,Phys.Rev.B47,5459(1993). [15] Y. Ran, M. Hermele, P. A. Lee, and X. -G. Wen, Phys. [5] N.ElstnerandA.P.YoungPhys.Rev.B50,6871(1994). Rev. Lett.98, 117205 (2007). [6] P. Lecheminant, B. Bernu, C. Lhuillier, L. Pierre and P. [16] J. S. Helton, K. Matan, M. P. Shores, E. A. Nytko, B. Sindzingre, Phys.Rev.B 56, 2521 (1997). M. Bartlett, Y. Yoshida, Y. Takano, A. Suslov, Y. Qiu, [7] F. Mila, Phys. Rev.Lett. 81, 2356-2359 (1998). J.-H. Chung, D. G. Nocera, and Y. S. Lee, Phys. Rev. [8] C. Waldtmann, H. U. Everts, B. Bernu, C. Lhuillier, P. Lett. 98, 107204 (2007). Sindzingre, P. Lecheminant and L. Pierre, Eur. Phys. J. [17] O. Ofer, A. Keren, E. A. Nytko, M. P. Shores, B. B 2, 501 (1998). M. Bartlett, D. G. Nocera, C. Baines, and A. Amato, [9] G. Misguich and B. Bernu, Phys. Rev. B 71, 014417 cond-mat/0610540. (2005). [18] P. Mendels, F. Bert, M. A. de Vries, A. Olariu, A. Har- [10] J. B. Marston and C. Zeng, J. Appl. Phys. 69, 5962 rison, F. Duc, J. C. Trombe, J. Lord, A. Amato, and C. (1991). Baines, Phys.Rev.Lett. 98, 077204 (2007). [11] A.V.SyromyatnikovandS.V.Maleyev,Phys.Rev.B66, [19] T. Imai, E. A. Nytko, B. M. Bartlett, M. P. Shores, D. 132408 (2002). G. Nocera, cond-mat/0703141. [12] P. Nikolic and T. Senthil, Phys. Rev. B 68, 214415 [20] M.RigolandR.R.P.Singh,Phys.Rev.Lett.98,207204 5 (2007). [24] B. S. Shastry and B. Sutherland, Physica B 108, 1069 [21] J.Oitmaa,C.HamerandW-H.Zheng,SeriesExpansion (1981). Methods for strongly interacting lattice models, (Cam- [25] W. H. Zheng, J. Oitmaa and C. J. Hamer, Phys. Rev. bridge UniversityPress 2006). B65, 014408 (2002). [22] M.Rigol,T.Bryant,andR.R.P.Singh,Phys.Rev.Lett. [26] V. Elser, Phys. Rev.Lett. 62, 2405 (1989). 97, 187202 (2006). [23] A.B. Harris, Phys. Rev.B7, 3166 (1973).