Classical spin-liquid on the maximally frustrated honeycomb lattice
aa r X i v : . [ c ond - m a t . s t r- e l ] O c t Classical spin-liquid on the maximally frustrated honeycomb lattice
J. Rehn, Arnab Sen, Kedar Damle, and R. Moessner Max-Planck-Institut f¨ur Physik komplexer Systeme, 01187 Dresden, Germany Department of Theoretical Physics, Indian Association for the Cultivation of Science, Kolkata 700032, India Department of Theoretical Physics, Tata Institute of Fundamental Research, Mumbai 400 005, India (Dated: October 2, 2018)We show that the honeycomb Heisenberg antiferromagnet with J / J = J , where J / / are first-, second- and third-neighbour couplings respectively, forms a classical spin liquid withpinch-point singularities in the structure factor at the Brillouin zone corners. Upon dilution withnon-magnetic ions, fractionalised degrees of freedom carrying 1 / PACS numbers:
Motivation. —The honeycomb lattice has – somewhatbelatedly – become one of the prime hunting grounds forspin liquids (SL) in d = 2 [1], in addition to the kagomeand the J − J square lattice Heisenberg models, whichhave been the focus of much attention over decades, con-tinuing until today. In both these latter cases [2–11],confidence in the existence of a quantum SL state for S = 1 / { N a, Li } IrO , provided aHeisenberg term is added [23–25].In fact, detailed studies of these materials suggest thatfurther nearest neighbor terms play an important role inexplaining spiral ordering at low temperatures [26], andone of the models studied in some detail is the J − J − J Heisenberg model, which had already been subjectto considerable earlier attention [27–30]. In determiningthe Hamiltonian appropriate to these materials, it hasturned out to be instructive to consider their response todisorder [31].Here, we identify and study in detail an unusual, hith-erto overlooked, classical SL state on the honeycomblattice, associated with the (known) degeneracy point J / J = J of the Heisenberg model on the honey-comb lattice. It exhibits remarkable new features. Thesearise from the fact that the dual lattice, as well as theunderlying Bravais lattice, is the tripartite triangular lat-tice. They include pinch points in the structure factor atthe zone corner wavevector Q (which distinguishes be- FIG. 1: (color online). Projection of the octahedron into thehexagon and the J − J − J model on the honeycomb lattice.The J interactions are differentiated with colors. tween the three sublattices), as well as novel disorder ef-fects whereby, upon dilution, fractionalised moments car-rying one third of the microscopic spin moment appear.These fractionalized moments interact via a frustrated,sublattice-dependent, long range interaction in the limitof low temperature, T .This model is further remarkable as it can be thoughtof the first realisation of a SL in d = 2 of edge-sharing simplices, which here take the form of octahedra. In ad-dition, its XY version does exhibit nematic order by dis-order, which turns out to be straightforwardly detectablein neutron scattering through the appearance of peaks inthe structure factor.The remainder of this paper is organised as follows.We first introduce the model and derive and describe itsSL, for which we formulate a novel low-energy descrip-tion. We then study its behaviour under dilution. Allour analytical predictions are supported by Monte Carlo(MC) simulations of the microscopic Hamiltonian. Weclose with an outlook, in which we argue that this modelis quite natural, as (i) the degeneracy point correspondsto a reasonably natural set of parameter values; and (ii)we expect the SL to fan out as T is increased, at the -2 π - π π π q x -2 π - π π π q y π - π π π q x -2 π - π π π q y FIG. 2: (color online). Structure factor as obtained in MonteCarlo simulations of the pristine Heisenberg (left) and XY(right) systems. Both results correspond to N = 1800 spinsat T /J = 0 . expense of adjacent phases exhibiting lower entropies. Model. —The Hamiltonian for classical O ( n ) spins ~S i of unit length on sites i of the honeycomb lattice reads: H = J X h i,j i ~S i · ~S j + J X hh i,k ii ~S i · ~S k + J X hhh i,l iii ~S i · ~S l = J X α ( ~S α ) + const. , (1)where h i, j i , hh i, k ii and hhh i, l iii refer respectively tofirst, second and third nearest neighbour pairs, while thesecond line follows from fixing J / J = J .This form shows that each and any configuration whereeach hexagon, labelled by α , has vanishing total spin, ~S α = 0, is a ground state. Such a rewriting is oftenhelpful for geometrically frustrated lattices. It is most of-ten used for ‘corner-sharing’ structure of elementary sim-plices [19, 20], examples being pyrochlore (corner-sharingtetrahedra) or kagome (corner-sharing triangles) lattices.It immediately allows to estimate the dimensionality ofthe ground state manifold, F . This proceeds by subtract-ing the number of constraints, K , imposed by Eq. 1, fromthe total number of degrees of freedom, D , of the spinsystem.For a system of n -component spins with N such sim-plices, and each spin part of b simplices, D = q ( n − /b per simplex, where the number of spins in a simplex q = 3 , , K = n constraints, as eachcomponent of its total spin must vanish. Hence, F = q ( n − b − n. (2)To maximize F , and hence enhance the chance of findinga SL [19, 20], one thus should minimize b , or maximise n and q . Indeed, b is minimal for corner-sharing arrange-ments, and q = 4, n = 3 result in the well-establishedclassical SL on the pyrochlore lattice. Triangle-based lat-tices (kagome has q = 3) need higher, n ≥
4, componentspins for a similar SL to arise [32]. The J − J model on the square lattice with J = J / q = 4; it does not support F > n . Indeed, nosuch Heisenberg model with F > F = 1 for q = 6 and b = 3,which corresponds to the frustration point of the honey-comb lattice, Eq. 1! It can be thought of as edge-sharing octahedra (Fig. 1), and thus presents the first instance ofa possible SL on an edge-sharing lattice. It is also thefirst with b >
2, a fact with significant consequences aswe explore in detail below.Before we do this, we demonstrate that this model doesindeed exhibit a SL with algebraic correlations for T → n or self-consistent Gaussianapproximation [33]) and comparing the result to clas-sical MC simulations. These yield a T = 0 structurefactor presenting pinch points , the defining characteristicof such algebraic SLs [34]. Somewhat unusually, in thiscase the pinch points are located at the corners of theBrillouin zone.From our MC simulations for Heisenberg and XYspins, we plot structure factor and specific heat onFigs. 2, 3. The MC simulations employ a combinationof heat-bath and microcanonical moves as well as par-allel tempering moves. The structure factor from MCsimulation of Heisenberg spins agrees with the analyticalprediction for soft spins. By contrast, for n = 2, the cor-responding XY model, low temperature peaks developin addition to the pinch points. This is an instance ofnematic (collinear) order by disorder, as is readily veri-fied by constraint counting [19, 20]. We note in passingthat the appearance of these peaks provides an unusuallydirect signature of collinear ordering.This interpretation is confirmed by a low- T specificheat of c = 0 . k B per spin, reduced from the value of c = bnk B / q expected from equipartition in the absenceof order by disorder [15–19], as is found in the Heisenbergmagnet with c = 0 . k B .The existence of pinch points in the Heisenberg casecomes as somewhat of a surprise given the non-bipartitenature of the dual triangular lattice. In the correspond-ing corner-sharing models, the bipartiteness of the duallattice (square, honeycomb or diamond lattice) is a cru-cial ingredient for such pinch points [34]. Indeed, in workclose in spirit to the present one, on bosons on a honey-comb and the dual triangular lattice [35], one finds anIsing emergent gauge field implying the absence of pinchpoints. The way this issue resolves itself in the presentcase is quite interesting: First, the pinch points are lo-cated at the Brillouin zone corners, corresponding to athree-sublattice wave vector Q . Second, the low-energydescription is naturally expressed in terms of a vectorfield that captures slow modulations near wavevector- Q ,reminiscent of the two dimensional height field acting as c T/J c T/J
FIG. 3: (color online). Specific heat obtained in Monte Carlosimulations of the pristine Heisenberg ( n = 3, top) and XY( n = 2, bottom) systems, with c = k B < bn q k B indicatingnematic order by disorder for the XY case. an emergent U (1) gauge field [34]. In detail, this proceedsas follows.Consider an A-sublattice (B-sublattice) site ~r A ( ~r B )of the honeycomb lattice, which sits at the center ofan “up-pointing” (“down-pointing”) triangle comprisingdual lattice points ~R a , ~R b and ~R c belonging to the threesublattices of the tripartite dual triangular lattice. Onewrites the corresponding O ( n ) spins ~S ~r in terms of ~ζ ~R and ~τ ~R , two O ( n ) vector fields on the dual triangularlattice. ~S ~r A = X α = a,b,c ( ~τ ~R α + ~ζ ~R α ) , ~S ~r B = X α = a,b,c ( ~τ ~R α − ~ζ ~R α ) . In the self-consistent Gaussian approximation, the parti-tion function for the Hamiltonian Eq. 1 can be writtenas a product of ~ζ and ~τ partition functions, with actions S ζ = F ( { ~ζ } ) , S τ = F ( { ~τ } ) + F ( { ~τ } ) , (3)where F ( { ~v } ) = ρ X ~r ( ~v ~R a ( ~r ) + ~v ~R b ( ~r ) + ~v ~R c ( ~r ) ) (4) for { ~v } = { ~ζ } or { ~τ } , and F ( { ~τ } ) = βJ X ~R (6 ~τ ~R + 2 X ~R n ∈ ∂ ~R ~τ ~R n ) (5)here, ~R n ∈ ∂ ~R denotes the six dual triangular lattice sites ~R n that are nearest neighbours of the dual triangularlattice site ~R . The stiffness constant ρ is adjusted toyield h ~S ~r i = 1.This action implies that ~ζ encodes the T = 0 fluc-tuations of the classical SL, while ~τ captures thermalfluctuations. The T → ~τ fieldsdo not contribute. This action, as well as the expres-sions for the physical spins ~S , are both invariant under ~ζ ( ~R ) → ~ζ ( ~R ) + Re( ~χ exp(2 πi Q · ~R )) for any ~χ . Dilution Effects. —The ground states of SLs often areless revealing of their topological nature than their ex-citations. An elegant way to visualise the latter as ef-fectively a ground state property is to introduce disorderwhich then nucleates excitations. In SLs, this is perhapsmost easily done by replacing some of the magnetic ionswith non-magnetic ones. For classical SLs, this dilutionproblem has been studied in some detail both experimen-tally [36–39], and theoretically [40–43]. In particular, forthe cases of SCGO, the checkerboard and the pyrochlorelattices, it was found that fractional impurity momentscarrying one half of the moment of a free spin arise as acooperative phenomenon. These so-called orphan spinsoccur when all but one of the spins of a simplex are re-placed – so that the total spin of that simplex (see Eq. 1)can no longer possibly vanish.These orphans turn out to provide a number of signa-tures of the new structure of the honeycomb SL. Firstof all, they directly reflect the fact that we have b = 3 edge-sharing octahedra meeting in each site – the frac-tional impurity moment is not one half but one third ofthat of a free spin! This is displayed in Fig. 4 (top panel)where a calculation based on a hybrid hard-soft spin the-ory [42, 43] is compared with numerical results for thelocal susceptibility. This is, to our knowledge, the firstinstance of fractionalisation into three items in a classicalspin model.Interactions between these orphans are entropic in na-ture and take the form of an effective Heisenberg ex-change J eff . They are mediated by the bulk SL, andhence reflect the structure of the latter. In the clas-sical SLs known so far, these effective interactions canbe written in a form which is uniformly antiferromag-netic [44]. Here, this is not possible: We now find thatthese interactions are antiferromagnetic / ferromagneticfor orphans residing on the same / different sublatticeof the dual triangular lattice, respectively, with the an-tiferromagnetic interactions being twice as strong as theferromagnetic ones. This intricate structure in the ef-fective exchange couplings follows from our field theory, ( M d - M u ) / S β hS JS β = 1002001000B (h,T) ( χ d - χ u ) / S T/(JS ) h/J = 0 JS /(27T)JS /(3T) 〈 S z ( ) S z ( x ) 〉 / T x T β = 6412825651210242048 FIG. 4: (color online).
Top: ‘Impurity magnetization’, de-fined as the difference of total magnetization in the dilutedand undiluted systems, as observed in MC simulations of themodel. The solid curve corresponds to the theoretical pre-diction for a free spin S/ h , i.e., the Langevinfunction B S/ ( h, T ). The inset shows the ‘impurity suscepti-bility’ at zero external field, consistent with a Curie law forfractionalized spins S/ Bottom:
Testing the scaling pre-diction for the charge correlations on a finite lattice of linearsize L = 210 using the expression for correlations within thesoft spin approach. Crucially, correlations between sites onthe same sublattice have been multiplied by an extra scalingfactor of − / which relates these entropic interactions to a bulk prop-erty of the pristine spin liquid, namely the correlationsbetween the thermally excited net spins ~S (Eq. 1): βJ eff ≈ −h ~S ,~r · ~S ,~r ih ~S ,~r · ~S ,~r i , (6)For low T and large distances | ~r − ~r | ≫ a , where a is the lattice spacing, this gives a scaling form: βJ eff = η ( ~r , ~r ) F (( ~r − ~r ) √ T ) (7) T → = 12 π η ( ~r , ~r ) log( ~r − ~r ) , (8) where η = +1 ( η = − /
2) if the orphans are on thesame (different) sublattices of the dual triangular lattice.This is verified using the analytical large-n result for afinite lattice, Fig. 4. In the limit T → βJ eff exhibits along-ranged logarithmic form. Outlook. —Our model, notwithstanding its simplicity,displays a plethora of phenomena of current interest; theunusual emergent ~τ fields and the new fractionalized be-havior of 1 / . . S = 1 /
2. Quitegenerally, at finite T , the classical SL behavior will befavoured over competing phases on account of its largeentropy, and in particular fan out from the degeneracypoint. We hope that this work will incite further investi-gation on appropriate honeycomb materials. Acknowledgements:
The work of AS is partly supported through the Part-ner Group program between the Indian Association forthe Cultivation of Science (IACS) and the Max Planck In-stitute for the Physics of Complex Systems. JR acknowl-edges the hospitality of IACS during the initial stagesof this work. This work was in part supported by DFGvia SFB 1143. We thank G. Jackeli, A. L¨auchli and J.Richter for discussions. [1] P. W. Anderson, Mat. Res. Bulletin , 153 (1973).[2] P. Chandra and B. Doucot, Phys. Rev. B , 9335(R)(1988).[3] Elbio Dagotto and Adriana Moreo, Phys. Rev. Lett. ,2148 (1989).[4] F. Figueirido, A. Karlhede, S. Kivelson, S. Sondhi, M.Rocek, and D. S. Rokhsar, Phys. Rev. B , 4619 (1990).[5] Subir Sachdev, Phys. Rev. B , 12377 (1992).[6] M. E. Zhitomirsky and Kazuo Ueda, Phys. Rev. B ,9007 (1996).[7] R. R. P. Singh, W. Zheng, J. Oitmaa, O. P. Sushkov, andC. J. Hamer, Phys. Rev. Lett. , 017201 (2003).[8] Matthieu Mambrini, Andreas L¨auchli, Didier Poilblanc,and Fr´ed´eric Mila, Phys. Rev. B , 144422 (2006). [9] J. Richter and J. Schulenburg, Eur. Phys. J. B , 117(2010).[10] Simeng Yan, David A. Huse, Steven R. White, Science , 1173 (2011).[11] Laura Messio, Bernard Bernu, and Claire Lhuillier, Phys.Rev. Lett. , 207204 (2012).[12] J. Villain, R. Bidaux, J.-P. Carton and R. Conte, J.Physique , 1263 (1980).[13] P. Chandra, P. Coleman, and A. I. Larkin, Phys. Rev.Lett. , 88 (1990).[14] J. N. Reimers and A. J. Berlinsky, Phys. Rev. B , 9539(1993).[15] J. T. Chalker, P. C. W. Holdsworth, and E. F. Shender,Phys. Rev. Lett. , 855 (1992).[16] I. Ritchey, P. Chandra and P. Coleman, Phys. Rev. B ,15342(R) (1993).[17] M. E. Zhitomirsky, Phys. Rev. B , 094423 (2008).[18] Gia-Wei Chern and R. Moessner, Phys. Rev. Lett. ,077201 (2013).[19] R. Moessner and J. T. Chalker, Phys. Rev. B , 12049(1998).[20] R. Moessner and J. T. Chalker, Phys. Rev. Lett. , 2929(1998).[21] Z. Y. Meng, T. C. Lang, S. Wessel, F. F. Assaad and A.Muramatsu , Nature , 847 (2010).[22] A. Yu. Kitaev, Annals of Physics , 2 (2003).[23] Yogesh Singh, S. Manni, J. Reuther, T. Berlijn, R.Thomale, W. Ku, S. Trebst, and P. Gegenwart, Phys.Rev. Lett. , 127203 (2012).[24] I. I. Mazin, Harald O. Jeschke, Kateryna Foyevtsova,Roser Valent´ı, and D. I. Khomskii, Phys. Rev. Lett. ,197201 (2012).[25] G. Jackeli and G. Khaliullin, Phys. Rev. Lett. ,017205 (2009).[26] Johannes Reuther, Ronny Thomale, and Stephan Rachel,Phys. Rev. B , 100405(R) (2014).[27] J. B. Fouet, P. Sindzingre, C. Lhuillier, Eur. Phys. J. B , 241 (2001).[28] Johannes Reuther, Dmitry A. Abanin, and RonnyThomale, Phys. Rev. B , 014417 (2011).[29] A. F. Albuquerque, D. Schwandt, B. Het´enyi, S. Cap-poni, M. Mambrini, and A. M. L¨auchli, Phys. Rev. B , 024406 (2011).[30] R. F. Bishop, P. H. Y. Li, O. G¨otze, J. Richter, and C.E. Campbell, arxiv:1504.02275 (2015).[31] Eric C. Andrade and Matthias Vojta, Phys. Rev. B ,205112 (2014).[32] David A. Huse and Andrew D. Rutenberg, Phys. Rev. B , 7536(R) (1992).[33] D. A. Garanin and Benjamin Canals, Phys. Rev. B ,443 (1999).[34] C. L. Henley, Annu. Rev. Condens. Matter Phys. , 179(2010).[35] O. I. Motrunich, Phys. Rev. B , 115108 (2003).[36] X. Obradors, A. Labarta, A. Isalgu´e, J. Tejada, J. Ro-driguez, and M. Pernet, Sol. State Commun. , 189(1988).[37] P. Schiffer and I. Daruka, Phys. Rev. B , 13712 (1997).[38] P. Mendels, A. Keren, L. Limot, M. Mekata, G. Collin,and M. Horvati´c, Phys. Rev. Lett. , 3496 (2000).[39] A. D. LaForge, S. H. Pulido, R. J. Cava, B. C. Chan, andA. P. Ramirez, Phys. Rev. Lett. , 017203 (2013).[40] C. L. Henley, Can. J. Phys. , 1307 (2001).[41] R. Moessner and A. J. Berlinsky, Phys. Rev. Lett. ,3293 (1999).[42] Arnab Sen, Kedar Damle, and R. Moessner, Phys. Rev.Lett. , 127203 (2011).[43] Arnab Sen, Kedar Damle, and R. Moessner, Phys. Rev.B , 205134 (2012).[44] J. Rehn, Arnab Sen, Alexei Andreanov, Kedar Damle, R.Moessner, and A. Scardicchio, Phys. Rev. B92