Reentrance of disorder in the anisotropic shuriken Ising model
Abstract
For a material to order upon cooling is common sense. What is more seldom is for disorder to reappear at lower temperature, which is known as reentrant behavior. Such resurgence of disorder has been observed in a variety of systems, ranging from Rochelle salts to nematic phases in liquid crystals. Frustration is often a key ingredient for reentrance mechanisms. Here we shall study a frustrated model, namely the anisotropic shuriken lattice, which offers a natural setting to explore an extension of the notion of reentrance between magnetic disordered phases. By tuning the anisotropy of the lattice, we open a window in the phase diagram where magnetic disorder prevails down to zero temperature. In this region, the competition between multiple disordered ground states gives rise to a double crossover where both the low and hightemperature regimes are less correlated than the intervening classical spin liquid. This reentrance of disorder is characterized by an entropy plateau, a multistep Curie law crossover and a rather complex diffuse scattering in the static structure factor. Those results are confirmed by complementary numerical and analytical methods: Monte Carlo simulations, Husimitree calculations and an exact decorationiteration transformation.
pacs:
75.10.Hk,75.30.Kz,75.10.KtRecent progress in frustrated magnetism has delivered entire maps of longrange ordered and disordered phases, obtained for example via the variation of bond anisotropy Savary and Balents (2012); Lee et al. (2012); Hao et al. (2014); Chernyshev and Zhitomirsky (2014); Götze and Richter (2015); Oitmaa and Singh (2016); Essafi et al. (2016) or further nearestneighbour couplings Messio et al. (2012); He and Chen (2015); He et al. (2015); Bieri et al. (2015); Iqbal et al. (2015); McClarty et al. (2015); Henelius et al. (2016). Such phase diagrams have allowed to put a series of frustrated materials onto a global and connected map, that can be experimentally explored via physical or chemical pressure Zhou et al. (2012); Dun et al. (2014); Wiebe and Hallas (2015); Rau and Gingras (2015). On such phase diagrams, when two ordered phases meet, an enhancement of the classical groundstate degeneracy takes place Yan et al. (2013). This degeneracy can either be lifted by thermal fluctuations, giving rise to multiple phase transitions Jaubert et al. (2015); Robert et al. (2015), or may destroy any kind of order down to (theoretically) zero temperature. This is where spin liquids appear. But this picture is less clear at the frontier between ordered and (possibly multiple) disordered ground states. In particular how do disordered phases compete with each other at finite temperature ?
The frustrated shuriken lattice Nakano and Sakai (2013) – also known as squarekagome Siddharthan and Georges (2001); Richter et al. (2004); Derzhko and Richter (2006); Richter et al. (2009); Derzhko et al. (2013); Nakano et al. (2014, 2015); Ralko and Rousochatzakis (2015), squagome Tomczak and Richter (2003); Glaetzle et al. (2014), squakagome Rousochatzakis et al. (2013) or L4L8 Rousochatzakis et al. (2013) lattice – provides an interesting modelexample for such competition. Being made of cornersharing triangles, it is locally similar to the famous kagome lattice, but with the important difference that the shuriken lattice is composed of two inequivalent sublattices [see Fig. 1]. Such asymmetry offers a natural setup for lattice anisotropy. In the asymptotic limits of this anisotropy, a promising zerotemperature phase diagram has emerged for quantum spin, ranging from a bipartite longrange ordered phase to a highly degenerate ground state made of tetramer clusters of spins Rousochatzakis et al. (2013). However, while the quantum ground states Richter et al. (2009); Nakano and Sakai (2013); Rousochatzakis et al. (2013); Ralko and Rousochatzakis (2015) and the influence of a magnetic field Richter et al. (2004); Derzhko and Richter (2006); Richter et al. (2009); Derzhko et al. (2013, 2014); Nakano and Sakai (2013); Nakano et al. (2014, 2015) have been studied to some extent, little is known about the finitetemperature properties in zero field Siddharthan and Georges (2001); Tomczak and Richter (2003).
In this paper, our goal is to develop a comprehensive and precise understanding of the frustrated phase diagram of the Ising model on the anisotropic shuriken lattice, relying on a combination of numerical and analytical methods (Monte Carlo simulations, Husimi tree calculations and decorationiteration transformation). Using the lattice anisotropy as a tuning parameter, we find that this model supports two longrange ordered phases (ferromagnet (FM) and staggered ferromagnet (SFM)), two classical spin liquids (SL) characterized by complex static structure factors, and a zerotemperature paramagnet composed of two kinds of isolated (super)spins with strictly zero correlations between them. We shall refer to this latter phase as a binary paramagnet. Over an extended region of the phase diagram, there is a double crossover from the hightemperature paramagnet to the spin liquids and finally into the lowtemperature binary paramagnet. This double crossover gives rise to a nonmonotonic behavior of the correlation length, which can be seen as an analogue of reentrant behavior between disordered phases. As a byproduct, we notice an essentially perfect agreement between Husimitree analytics and Monte Carlo simulations in the disordered regimes.
The paper is divided as follows. The model is introduced in section I, followed by its phase diagram in section II. In section III, we analyze in details the doublecrossover region between disordered regimes. We conclude the paper by discussing possible experimental realizations of the shuriken lattice and summarizing our results in sections IV and V respectively. Most technical details are given in the appendices.
I Anisotropic shuriken model
The shuriken lattice is made of cornersharing triangles with 6 sites per unit cell [see Fig. 1]. As opposed to its kagome parent where all spins belong to hexagonal loops, the shuriken lattice forms two kinds of loops made of either 4 or 8 sites. As a consequence, of the spins in the system belong to the Asublattice, while the remaining of the spins form the Bsublattice. Let us respectively define and as the coupling constants between Asites on the square plaquettes, and between A and Bsites on the octagonal plaquettes. The Hamiltonian of the model can be written as:
(1) 
where we consider Ising spins with nearestneighbor coupling.
There is no frustration for ferromagnetic where the system undergoes a phase transition with spontaneous symmetry breaking for . We shall thus focus on antiferromagnetic , which will be our energy and temperature scale of reference. The thermodynamics will be discussed as a function of the coupling ratio [Derzhko and Richter, 2006,Rousochatzakis et al., 2013,Ralko and Rousochatzakis, 2015]
(2) 
with ferro and antiferromagnetic .
Ii Phase Diagram
The Hamiltonian of equation (1) is invariant under the transformation
(3) 
All quantities derived from the energy, and especially the specific heat and entropy , are thus the same for and . Their respective magnetic phases are related by reversing all spins of the Asublattices.
ii.1 Longrange order:
When the octagonal plaquettes are dominating (), the shuriken lattice becomes a decorated square lattice, with Asites sitting on the bonds between Bsites. Being bipartite, the decorated square lattice is not frustrated and orders via a phase transition of the 2D Ising Universality class ChunFeng et al. (2006) by spontaneous symmetry breaking. Nonuniversal quantities such as the transition temperature can be exactly computed by using the decorationiteration transformation Fisher (1959); ChunFeng et al. (2006); Strecka and Jascur (2015) [see appendix D.2]
(4) 
The lowtemperature ordered phases, displayed in Fig. 2.(b) and 2.(c), remain the ground states of the anisotropic shuriken model for and respectively. The persistence of the 2D Ising Universality class down to is not necessarily obvious, but is confirmed by finitesize scaling from Monte Carlo simulations [see appendix C].
These two ordered phases are respectively ferromagnetic (FM, ) and staggered ferromagnetic (SFM, ) [see Fig. 2.(b,c)]. The staggering of the latter comes from all spins on square plaquettes pointing in one direction, while the remaining ones point the other way. This leads to the rather uncommon consequence that fully antiferromagnetic couplings – both and are negative for – induce longrange ordered (staggered) ferromagnetism, reminiscent of Lieb ferrimagnetism Lieb (1989) as pointed out in Ref. [Rousochatzakis et al., 2013] for quantum spins. The existence of ferromagnetic states among the set of ground states of Ising antiferromagnets is not rare, with the triangular and kagome lattices being two famous examples. But such ferromagnetic states are usually part of a degenerate ensemble where no magnetic order prevails on average. Here the lattice anisotropy is able to induce ferromagnetic order in an antiferromagnetic model by lifting its groundstate degeneracy at (see below). This is interestingly quite the opposite of what happens in the spinice model Harris et al. (1997), where frustration prevents magnetic order in a ferromagnetic model by stabilizing a highly degenerate ground state.
ii.2 Binary paramagnet:
The central part of the phase diagram is dominated by the square plaquettes. The ground states are the same for all . A sample configuration of these ground states is given in Fig. 2.(d), where antiferromagnetically ordered squareplaquettes are separated from each other via spins on sublattice B. The antiferromagnetic squareplaquettes locally order in two different configurations equivalent to a superspin with Ising degree of freedom.
(5) 
where the site indices are given in Fig. 1. These superspins are the classical analogue of the tetramer objects observed in the spin model Rousochatzakis et al. (2013). At zero temperature, the frustration of the bonds perfectly decouples the superspins from the Bsites. The system can then be seen as two interpenetrating square lattices: one made of superspins, the other one of Bsites. We shall refer to this phase as a binary paramagnet (BPM).
The perfect absence of correlations beyond square plaquettes at allows for a simple determination of the thermodynamics. Let and be respectively the total number of unit cells and spins in the system, and be the statistical average of . There are square plaquettes and Bsites, giving rise to an extensive groundstate entropy
(6) 
which turns out to be half the entropy of an Ising paramagnet. As for the magnetic susceptibility , it diverges as . But the reduced susceptibility , which is nothing less than the normalized variance of the magnetization
(7)  
converges to a finite value in the BPM
(8) 
ii.3 Classical spin liquid:
There is a sharp increase of the groundstate degeneracy at , when the binary paramagnet and the (staggered) ferromagnet meet. As is common for isotropic trianglebased Ising antiferromagnets, 6 out of 8 possible configurations per triangle minimize the energy of the system. As opposed to the BPM one does not expect a cutoff of the correlations [see section III.3], making these phases cooperative paramagnets Villain (1979), also known as classical spin liquids.
Due to the high entropy of these cooperative paramagnets, the SL phases spread to the neighboring region of the phase diagram for and , continuously connected to the hightemperature paramagnet [see Fig. 2]. Hence, for , the anisotropic shuriken model stabilizes a cooperative paramagnet above a nondegenerate^{1}^{1}1besides the trivial timereversal symmetry longrange ordered phase. This is a general property of classical spin liquids when adiabatically tuned away from their highdegeneracy point, as observed for example in Heisenberg antiferromagnets on the kagome Elhajal et al. (2002) or pyrochlore Canals et al. (2008); Chern (2010); McClarty et al. (2014) lattices, and possibly in the material of ErSnO Yan et al. (2013). For on the other hand, multiple crossovers take place upon cooling which deserves a dedicated discussion in the following section III.
Iii Reentrance of disorder
iii.1 Double crossover
First of all, panels (a) and (c) of Fig. 3 confirm that the classical spin liquids and binary paramagnet persist down to zero temperature for and respectively, and that all models for have extensively degenerate ground states. For there is a double crossover indicated by the double peaks in the specific heat of Fig. 3.(b). These peaks are not due to phase transitions since they do not diverge with system size. The double crossover persists for . Upon cooling, the system first evolves from the standard paramagnet to a spin liquid before entering the binary paramagnet. The intervening spin liquid takes the form of an entropy plateau for [see Fig. 3.(b)], at the same value as the lowtemperature regime for [see Fig. 3.(a)]. All relevant thermodynamic quantities are summarized in Table 1.
While the mapping of equation (3) ensures the invariance of the energy, specific heat and entropy upon reversing to , it does not protect the magnetic susceptibility. The build up of correlations in classical spin liquids is known to give rise to a Curielaw crossover Jaubert et al. (2013) between two asymptotic regimes of the susceptibility, as observed in pyrochlore Isakov et al. (2004); Ryzhkin (2005); Conlon and Chalker (2010); Jaubert et al. (2013), triangular Isoda (2008) and kagome Isoda (2008); Li et al. (2010); Macdonald et al. (2011) systems. This is also what is observed here on the anisotropic shuriken lattice for [see Fig. 4]. But for intermediate models with , the double crossover makes the reduced susceptibility nonmonotonic. first evolves towards the values of the spin liquids SL (resp. SL) for (resp. ) before converging to in the binary paramagnet.
Beyond the present problem on the shuriken lattice, this multistep Curielaw crossover underlines the usefulness of the reduced susceptibility to spot intermediate regimes, and thus the proximity of different phases. From the point of view of renormalization group theory, the coordinates of the phase diagram are fixed points which deform the renormalization flows passing in the vicinity.
Monte Carlo  Husimi tree  exact  

n/a  
n/a  
n/a  
iii.2 Decorationiteration transformation
The phase diagram of the anisotropic shuriken model and, in particular, the double crossover observed for [see Fig. 2] can be further understood using an exact mapping to an effective model on the checkerboard lattice, a method known as decorationiteration transformation [see Ref. [Strecka and Jascur, 2015] for a review]. In short, by summing over the degrees of freedom of the Aspins, one can arrive at an effective Hamiltonian involving only the Bspins, which form a checkerboard lattice. The coupling constants of the effective Hamiltonian are functions of the temperature and for they vanish at both high and low temperatures, but are finite for an intermediate regime. This intermediate regime may be identified as the SL cooperative paramagnets of Fig. 2, whereas the lowtemperature region of vanishing effective interaction corresponds to the binary paramagnet (BPM). This mapping is able to predict a nonmonotonic behavior of the correlation length.
In this section we give a brief sketch of the derivation of the effective model, before turning to its results. Details of the effective model are given in Appendix D.
To begin, consider the partition function for the system, with the Hamiltonian given by Eq. (1)
(9) 
where is the inverse temperature and the sums are over all possible spin configurations. Since in the Hamiltonian of Eq. (1) the square plaquettes of the Asites are only connected to each other via their interaction with the intervening Bsites, it is possible to directly take the sum over configurations of Aspins in Eq. (9) for a fixed (but completely general) configuration of Bspins. Doing so, we arrive at
(10) 
where the product is over all the square plaquettes of the lattice and is a function of the four Bspins immediately neighbouring a given square plaquette. The Bspins form a checkerboard lattice, and Eq. (10) can be exactly rewritten in terms of an effective Hamiltonian on that lattice:
(11)  
(12) 
where is a sum over checkerboard plaquettes of Bspins. The effective Hamiltonian contains a constant term , a nearest neighbour interaction , a second nearest neighbour interaction , and a foursite ring interaction . All couplings are functions of temperature and are invariant under the transformation because the degrees of freedom of the Asites have been integrated out. Expressions for the dependence of the couplings on temperature are given in Appendix D.
The temperature dependence of the effective couplings can itself give rather a lot of information about the behavior of the shuriken model.
First we consider the case . In this regime of parameter space, all effective interactions vanish exponentially at low temperature . For intermediate temperatures the effective interactions in Eq. (12) become appreciable before vanishing once more at high temperatures. This is illustrated for the case in the upper panel of Fig. 5. Seeing the problem in terms of these effective couplings gives some intuition into the double crossover observed in simulations. As the temperature is decreased the effective couplings increase in absolute value and the system enters a short range correlated regime. However, as the temperature decreases further, the antiferromagnetic correlations on the square plaquettes of Aspins become close to perfect, and act to screen the effective interaction between Bspins. This is reflected in the exponential suppression of the couplings .
In the case , the effective interactions no longer vanish exponentially at low temperature, but instead vanish linearly
(13) 
The ratio of effective couplings to the temperature thus tends to a constant below , as shown in the lower panel of Fig. 5. Thus, the zero temperature limit of the shuriken model can be mapped to a finite temperature model on the checkerboard lattice for and to an infinite temperature model for .
The behavior of the spin correlations in the shuriken model can be captured by calculating the correlation length between Bspins in the checkerboard model. Since is small for all of the interactions , at all temperatures (see Fig. 5), this can be estimated using a perturbative expansion in . For two Bspins chosen such that the shortest path between them is along nearest neighbour bonds we obtain to leading order
(14)  
(15) 
where we choose units of length such that the linear size of a unit cell is equal to 1. Details of the calculation are given in Appendix D.
The correlation length between Bspins, calculated from Eq. (15), is shown for the cases and in Fig. 6. For the correlation length shows a nonmonotonic behavior, vanishing at both high and low temperature with a maximum at . On the other hand for , the correlation length enters a plateau for temperatures below and the system remains in a short range correlated regime down to . The extent of this plateau agrees with the lowtemperature plateau of the reduced susceptibility in Fig. 4
iii.3 Correlations and Structure factors
The nonmonotonic behavior of the correlation length estimated in the previous section III.2 can be measured by Monte Carlo simulations. Let us consider the microscopic correlations both in real () and Fourier () space. The function measures the correlation between a central spin and all spins at distance . Because of the nature of the binary paramagnet, one needs to make a distinction between central spins on the A and B sublattices. Let be the ensemble of sites at distance from a given spin on the sublattice. The correlation function is defined as
(16) 
where the absolute value accounts for the antiferromagnetic correlations. As for the static structure factor , it is defined as
(17) 
and are respectively plotted on the left and right of Fig. 7. Let us first consider what happens in absence of reentrant behavior. For [see panels (a,b)], the system is ferromagnetic at low temperature with over long length scales. Above the phase transition, the correlations are exponentially decaying.
When [see panels (c,d)], the correlations remain exponentially decaying down to zero temperature. The correlation length reaches a maximum in the spinliquid regime with . The quantitative superimposition of data for and is in agreement with the lowtemperature plateau of the correlation length in Fig. 6. The spin liquid remains essentially unchanged all the way up to , when defects are thermally excited. However even if the correlations are exponential, they should not be confused with paramagnetic ones, as illustrated by their strongly inhomogeneous structure factors [see Fig. 8 and Supplementary Materials].
Once one enters the doublecrossover region [see Fig. 7.(e,f) for ], the correlation function becomes nonmonotonic with temperature, as predicted from the analytics of Fig. 6. In the binary paramagnet, the Bsites are perfectly uncorrelated, while the Asites have a finite cutoff of the correlation that is the size of the square plaquettes (superspins). This is why takes the form of an array of dots of scattering, whose width is inversely proportional to the size of the superspins [see Fig. 8]. Please note that the dip of correlations for the nearestneighbors in Fig. 7.(e) is because half of the nearestneighbors of any Asite are on the uncorrelated B sublattice.
The intervening presence of the spin liquids between the two crossovers is conceptually reminiscent of reentrant behavior Hablützel (1939); Vaks et al. (1966); Cladis (1988). Not in the usual sense though, since reentrance is usually considered to be a feature of ordered phases surrounded by disordered ones. But the present scenario is a direct extension of the concept of reentrance applied to disordered regimes. This reentrance is quantitatively characterized at the macroscopic level by the doublepeak in the specific heat, the entropy plateau and the multistep Curielaw crossover of Fig. 3.(b), and microscopically by the nonmonotonic evolution of the correlations [see Figs. 6, 7 and 8]. As such, it provides an interesting mechanism to stabilize a gaslike phase “below” a spin liquid, where (a fraction of) the spins form fully correlated clusters which i) can then fluctuate independently of the other degreesoffreedom while ii) lowering the entropy of the gaslike phase below the one of the spin liquids.
Iv The shuriken lattice in experiments ?
Finally, we would like to briefly address the experimental situation. Unfortunately we are not aware of an experimental realization of the present model, but several directions are possible, each of them with their advantages and drawbacks.
The shuriken topology has been observed, albeit quite hidden, in the dysprosium aluminium garnet (DAG) Landau et al. (1971); Wolf et al. (1972) [see Ref. [Wolf, 2000] for a recent review]. The DAG material has attracted its share of attention in the 1970’s, but its microscopic Hamiltonian does not respect the geometry of the shuriken lattice – it is actually not frustrated – and is thus quite different from the model presented in equation (1). However it shows that the shuriken topology can exist in solid state physics.
Cold atoms might offer an alternative. Indeed, the necessary experimental setup for an optical shuriken lattice has been proposed in Ref. [Glaetzle et al., 2014]. The idea was developed in the context of spinice physics, i.e. assuming an emergent Coulomb gauge theory whose intrinsic Ising degrees of freedom are somewhat different from the present model. Nonetheless, optical lattices are promising, especially if one considers that the inclusion of “proper” Ising spins might be available thanks to artificial gauge fields Struck et al. (2013).
But the most promising possibility might be artificial frustrated lattices, where ferromagnetic nanoislands effectively behave like Ising degreesoffreedom. Since the early days of artificial spin ice Wang et al. (2006), many technological and fundamental advances have been made Nisoli et al. (2013). In particular, while the thermalization of the Isinglike nanoislands had been a longstanding issue, this problem is now on the way to be solved Kapaklis et al. (2012); Farhan et al. (2013a, b); Morgan et al. (2013); Marrows (2013); Anghinolfi et al. (2015); Arnalds et al. (2016). Furthermore, since the geometry of the nanoarray can be engineered lithographically, a rich diversity of lattices is available, and the shuriken geometry should not be an issue. Concerning the Ising nature of the degreesoffreedom, nanoislands have recently been grown with a magnetization axis perpendicular to the lattice Zhang et al. (2012); Chioar et al. (2014).
To compute their interaction Zhang et al. (2012); Chioar et al. (2014), let us define the Ising magnetic moment of two different nanoislands: and . The interaction between them is dipolar of the form
(18) 
where is the strength of the dipolar interaction and is the vector separating the two moments. The resulting coupling is thus antiferromagnetic and quickly decays with distance. Hence, at the nearestneighbour level, a physical distortion of the shuriken geometry – by elongating or shortening the distance between A and B sites – would precisely reproduce the anisotropy of equation (1) for . However, the influence of interactions beyond nearestneighbours has successively been found to be experimentally negligible Zhang et al. (2012) and relevant Chioar et al. (2014) on the kagome geometry. Thus the phase diagram of Fig. 2.(a) could possibly be observed at finite temperature, but would likely be influenced by longerrange interactions at relatively low temperature.
V Conclusion
The anisotropic shuriken lattice with classical Ising spins supports a variety of different phases as a function of the anisotropy parameter : two longrange ordered ones for (ferromagnet and staggered ferromagnet) and three disordered ones [see Fig. 2]. Among the latter ones, we make the distinction, at zero temperature, between two cooperative paramagnets SL for , and a phase that we name a binary paramagnet (BPM) for . The BPM is composed of locally ordered square plaquettes separated by completely uncorrelated single spins on the Bsublattice [see Fig. 2.(d)].
At finite temperature, the classical spin liquids SL spread beyond the singular points , giving rise to a double crossover from paramagnet to spin liquid to binary paramagnet, which can be considered as a reentrant behavior between disordered regimes. This competition is quantitatively defined by a doublepeak feature in the specific heat, an entropy plateau, a multistep Curielaw crossover and a nonmonotonic evolution of the spinspin correlation, illustrated by an inhomogeneous structure factor [see Figs. 3, 4,6, 7 and 8]. The reentrance can also be precisely defined by the resurgence of the couplings in the effective checkerboard model [see Fig. 5].
Beyond the physics of the shuriken lattice, the present work, and especially Fig. 3, confirms the Husimitree approach as a versatile analytical method to investigate disordered phases such as spin liquids. Regarding classical spin liquids, Fig. 4 illustrates the usefulness of the reduced susceptibility [Jaubert et al., 2013], whose temperature evolution quantitatively describes the successive crossovers between disordered regimes. Last but not least, we hope to bring to light an interesting facet of distorted frustrated magnets, where extended regions of magnetic disorders can be stabilized by anisotropy, such as on the Cairo Rousochatzakis et al. (2012); Rojas et al. (2012), kagome Li et al. (2010); Apel and Everts (2011) and pyrochlore Benton and Shannon (2015) lattices. Such connection is particularly promising since it expands the possibilities of experimental realizations, for example in Volborthite kagome Hiroi et al. (2001) or breathing pyrochlores Okamoto et al. (2013); Kimura et al. (2014).
Possible extensions of the present work can take different directions. Motivated by the counterintuitive emergence of valencebondcrystals made of resonating loops of size 6 [Ralko and Rousochatzakis, 2015], the combined influence of quantum dynamics, lattice anisotropy Rousochatzakis et al. (2013); Ralko and Rousochatzakis (2015) and entropy selection presented here should give rise to a plethora of new phases and reentrant phenomena. As an intermediary step, classical Heisenberg spins also present an extensive degeneracy at Richter et al. (2009); Rousochatzakis et al. (2013), where thermal orderbydisorder is expected to play an important role in a similar way as for the parent kagome lattice, especially when tuned by anisotropy . The addition of an external magnetic field Derzhko and Richter (2006); Nakano et al. (2015) would provide a direct tool to break the invariance by transformation of equation (3), making the phase diagram of Fig. 2.(a) asymmetric. Furthermore, the diversity of spin textures presented here offers a promising framework to be probed by itinerant electrons coupled to localized spins via doubleexchange.
Acknowledgements.
We are thankful to John Chalker, Arnaud Ralko, Nic Shannon and Mathieu Taillefumier for fruitful discussions and suggestions. This work was supported by the Theory of Quantum Matter Unit of the Okinawa Institute of Science and Technology Graduate University.Appendix A Methods
Classical Monte Carlo simulations have been performed based on the singlespinflip algorithm. Let a Monte Carlo step (MCs) be the standard Monte Carlo unit of time made of attempts to flip a spin chosen at random. Typical simulations in this paper consist of

MCs, including MCs for equilibration;

1 measurement every 50 MCs for (total of samples);

1 measurement every 10 MCs for (total of samples);

system sizes varying from to . Fig. 2 has been obtained for sites.
In order to improve the statistics, a large number of temperatures were simulated, and data were averaged over 4 neighboring temperatures.
To avoid any potential problems of ergodicity breaking in Monte Carlo simulations, we combined this numerical approach to analytical calculations on a Husimi tree Husimi (1950), a method that has already demonstrated success in frustrated magnets Chandra and Doucot (1994); Yoshida et al. (2002); Jaubert et al. (2013); Strečka and Ekiz (2015). In a nutshell, the Husimi tree is a recursive approach on a Bethe lattice where all vertices are replaced by a cluster of spins. The clusters are connected to each other via their external corners, without making any closed loops. This allows to correctly take into account the interactions within each cluster, where frustration can be encoded.
Because the shuriken lattice is made of cornersharing triangles, a natural choice would have been to consider triangles as building blocks of the Husimitree recursion. However a single triangle does not properly include the geometry of the anisotropy presented in Fig. 1. This is why, in the same way as for the 16vertex model Foini et al. (2013); Levis et al. (2013), we chose a larger building block made of four triangles forming a “shuriken” [see Fig. 1], which includes the anisotropy between A and Bsites.
On the other hand, it neglects correlations on the length scale of the octagonal plaquettes and beyond. As such, the Husimi tree remains a mean field approximation which can only be qualitative in the vicinity of a critical point below its upper critical dimension. Since the 2D Ising Universality class is obviously not mean field, the Husimi tree underestimates the transition temperatures for by a factor of . This is why the boundaries of the FM and SFM phases have been determined with Monte Carlo simulations [open circles in Fig. 2.(a)].
Appendix B Pauling entropy of the spin liquids
In the isotropic case (), and by symmetry for as well [see equation (3)], a simple Pauling argument is possible for the calculation of the entropy Pauling (1935). If is the number of Ising spins, then there are triangles in the system. Out of the possible configurations, the Pauling argument states that approximately are allowed in the ground state, giving a total number of ground states in the spin liquids SL
(19) 
giving an entropy
(20)  
Appendix C 2D Ising Universality class
For , the anisotropic shuriken model orders at low temperature via a spontaneous symmetry breaking [see Fig. 2]. We know it is a critical point of the 2D Ising Universality class for large [see section II.1]. In this appendix, our goal is to confirm numerically that it remains in the same Universality class as , by considering two different values of the coupling ratio: and . By symmetry of equation (3), the results directly apply to also.
In Fig. 9, we analyze the specific heat for four different system sizes . The transition temperature scales like to its thermodynamic limit found at
(21)  
(22) 
Based on these values of the transition temperature, we can define the reduced temperature . Following standard finite size scaling Landau and Binder (2009), we confirm in Figs. 10 and 11 that the nature of the phase transition is consistent with the 2D Ising Universality class with critical exponents
(23)  
Appendix D Details of the decorationiteration transformation
In this Appendix we give the details of the map to the effective model on the checkerboard lattice derived in Section III.2. We give the derivation in Section D.1 and then give details of the calculation of the correlation length in Section D.3.
d.1 Derivation of the effective model on the checkerboard lattice
Consider the partition function of the anisotropic shuriken model
(24) 
where and are respectively the Hamiltonian of the square plaquettes of A spins and the Hamiltonian coupling the intermediate Bspins to the square plaquettes. Summing over configurations of A spins, we obtain:
(25) 
where the product is over all the square plaquettes of the lattice and depends on the configuration of the four Bspins immediately neighbouring a given square plaquette.
There are sixteen possible arrangements of the four spins Bspins surrounding a square plaquette of which only four are inequivalent from the point of view of symmetry. These give rise to four possible values for :
(26)  
(27)  
(28)  
(29) 
From these we can assign “free energies” to each of the four possible inequivalent configurations of B spins around a square plaquette, i.e.
(30)  
(31)  
(32)  
(33) 
The Bspins form a checkerboard lattice as illustrated in Fig. 12. Using Eqs. (26)(33) we can rewrite Eq. (25) in terms of an effective Hamiltonian on the checkerboard lattice
(34) 
The sum is a sum over the elementary units of the checkerboard lattice. The function is a function only of the four Bspins around a checkerboard unit and returns one of the four defined in Eqs. (30)(33) as appropriate to the configuration of those four spins.
We can rewrite explicitly in terms of interactions between the spins on the checkerboard lattice. The resultant effective Hamiltonian for the spins on the checkerboard lattice contains a constant term , a nearest neighbour interaction , a second nearest neighbour interaction , and a foursite ring interaction .
(35) 
All couplings are functions of temperature .
The relationship between the temperature dependent couplings appearing in Eq. (35) and the free energies defined in Eqs. (30)(33) is:
(36)  
(37)  
(38)  
We have thus succeeded in mapping the original model on the shuriken lattice, onto an effective model on the checkerboard lattice [Eq. (35)].
d.2 Transition temperature of the decorated square lattice
In the limit , one obtains the decorated square lattice. Applying to Eqs. (26)(33) and then injecting the results into Eqs. (37)(LABEL:eq:Jring), one obtains
(40)  
(41) 
The term does not cancel, but it only appears as a prefactor in the partition function of Eq. (34) and thus does not influence the critical point.
Our effective model thereby becomes a square lattice with a temperature dependent nearestneighbour coupling . It is exactly soluble and the transition temperature is obtained by injecting Eq. (40) into Onsager’s solution of the Ising square lattice Onsager (1944)
(42)  
which gives the result of Eq. (4)
(43)  
d.3 Correlation length
We observed in Section III.2 that for the couplings of the effective model are small compared to the temperature, for all values of temperature.
An expansion of the partition function of the effective model in powers of is thus justified. Where this expansion is assymptotically exact in both high and low temperature regimes.
Here we show how to use this expansion to calculate the correlation function for a pair of Bspins. For simplicity and concreteness we will do the calculation for a pair separated by a path such as that in Fig. 13, where the shortest route between them traverses only bonds and contains such bonds. However, there is no difficulty in making the calculation for other cases.
We have
where is the total number of spin configurations of the checkerboard model.
The leading nonzero term in Eq. (LABEL:eq:HTE) comes from the part of the sum in the numerator, and corresponds to covering the shortest path between and with interactions. There are ways of ordering the product of terms, which cancels the occurring in the denominator. We thus obtain
In our choice of units of length, made such that the linear size of a unit cell equals 1, the distance between the spins is
(46) 
We therefore have a correlation length
(47) 
References
 Savary and Balents (2012) L. Savary and L. Balents, Phys. Rev. Lett. 108, 037202 (2012).
 Lee et al. (2012) S. B. Lee, S. Onoda, and L. Balents, Phys. Rev. B 86, 104412 (2012).
 Hao et al. (2014) Z. Hao, A. G. R. Day, and M. J. P. Gingras, Phys. Rev. B 90, 214430 (2014).
 Chernyshev and Zhitomirsky (2014) A. L. Chernyshev and M. E. Zhitomirsky, Phys. Rev. Lett. 113, 237202 (2014).
 Götze and Richter (2015) O. Götze and J. Richter, Phys. Rev. B 91, 104402 (2015).
 Oitmaa and Singh (2016) J. Oitmaa and R. R. P. Singh, Phys. Rev. B 93, 014424 (2016).
 Essafi et al. (2016) K. Essafi, O. Benton, and L. D. C. Jaubert, Nature Communications 7, 10297 (2016).
 Messio et al. (2012) L. Messio, B. Bernu, and C. Lhuillier, Phys. Rev. Lett. 108, 207204 (2012).
 He and Chen (2015) Y.C. He and Y. Chen, Phys. Rev. Lett. 114, 037201 (2015).
 He et al. (2015) Y.C. He, S. Bhattacharjee, F. Pollmann, and R. Moessner, Phys. Rev. Lett. 115, 267209 (2015).
 Bieri et al. (2015) S. Bieri, L. Messio, B. Bernu, and C. Lhuillier, Phys. Rev. B 92, 060407 (2015).
 Iqbal et al. (2015) Y. Iqbal, H. O. Jeschke, J. Reuther, R. Valentí, I. I. Mazin, M. Greiter, and R. Thomale, Phys. Rev. B 92, 220404 (2015).
 McClarty et al. (2015) P. A. McClarty, O. Sikora, R. Moessner, K. Penc, F. Pollmann, and N. Shannon, Phys. Rev. B 92, 094418 (2015).
 Henelius et al. (2016) P. Henelius, T. Lin, M. Enjalran, Z. Hao, J. G. Rau, J. Altosaar, F. Flicker, T. Yavors’kii, and M. J. P. Gingras, Phys. Rev. B 93, 024402 (2016).
 Zhou et al. (2012) H. D. Zhou, J. G. Cheng, A. M. Hallas, C. R. Wiebe, G. Li, L. Balicas, J. S. Zhou, J. B. Goodenough, J. S. Gardner, and E. S. Choi, Phys. Rev. Lett. 108, 207206 (2012).
 Dun et al. (2014) Z. L. Dun, M. Lee, E. S. Choi, A. M. Hallas, C. R. Wiebe, J. S. Gardner, E. Arrighi, R. S. Freitas, A. M. ArevaloLopez, J. P. Attfield, H. D. Zhou, and J. G. Cheng, Phys. Rev. B 89, 064401 (2014).
 Wiebe and Hallas (2015) C. R. Wiebe and A. M. Hallas, APL Materials 3, 041519 (2015).
 Rau and Gingras (2015) J. G. Rau and M. J. P. Gingras, Phys. Rev. B 92, 144417 (2015).
 Yan et al. (2013) H. Yan, O. Benton, L. D. C. Jaubert, and N. Shannon, arXiv:1311.3501 (2013).
 Jaubert et al. (2015) L. D. C. Jaubert, O. Benton, J. G. Rau, J. Oitmaa, R. R. P. Singh, N. Shannon, and M. J. P. Gingras, Phys. Rev. Lett. 115, 267208 (2015).
 Robert et al. (2015) J. Robert, E. Lhotel, G. Remenyi, S. Sahling, I. Mirebeau, C. Decorse, B. Canals, and S. Petit, Phys. Rev. B 92, 064425 (2015).
 Nakano and Sakai (2013) H. Nakano and T. Sakai, Journal of the Physical Society of Japan 82, 083709 (2013).
 Siddharthan and Georges (2001) R. Siddharthan and A. Georges, Physical Review B 65, 014417 (2001).
 Richter et al. (2004) J. Richter, O. Derzhko, and J. Schulenburg, Physical Review Letters 93, 107206 (2004).
 Derzhko and Richter (2006) O. Derzhko and J. Richter, The European Physical Journal B 52, 23 (2006).
 Richter et al. (2009) J. Richter, J. Schulenburg, P. Tomczak, and D. Schmalfuss, Condensed Matter Physics 12, 507 (2009).
 Derzhko et al. (2013) O. Derzhko, J. Richter, O. Krupnitska, and T. Krokhmalskii, Phys. Rev. B 88, 094426 (2013).
 Nakano et al. (2014) H. Nakano, T. Sakai, and Y. Hasegawa, Journal of the Physical Society of Japan 83, 084709 (2014).
 Nakano et al. (2015) H. Nakano, Y. Hasegawa, and T. Sakai, Journal of the Physical Society of Japan 84, 114703 (2015).
 Ralko and Rousochatzakis (2015) A. Ralko and I. Rousochatzakis, Physical Review Letters 115, 167202 (2015).
 Tomczak and Richter (2003) P. Tomczak and J. Richter, Journal of Physics A: Mathematical and General 36, 5399 (2003).
 Glaetzle et al. (2014) A. W. Glaetzle, M. Dalmonte, R. Nath, I. Rousochatzakis, R. Moessner, and P. Zoller, Phys. Rev. X 4, 041037 (2014).
 Rousochatzakis et al. (2013) I. Rousochatzakis, R. Moessner, and J. van den Brink, Physical Review B 88, 195109 (2013).
 Derzhko et al. (2014) O. Derzhko, J. Richter, O. Krupnitska, and T. Krokhmalskii, Low Temperature Physics 40, 513 (2014), 1312.1111 .
 ChunFeng et al. (2006) S. ChunFeng, K. XiangMu, and Y. XunChang, Communications in Theoretical Physics 45, 555 (2006).
 Fisher (1959) M. E. Fisher, Phys. Rev. 113, 969 (1959).
 Strecka and Jascur (2015) J. Strecka and M. Jascur, Acta Physica Slovaca 65, 235 – 367 (2015).
 Lieb (1989) E. H. Lieb, Phys. Rev. Lett. 62, 1201 (1989).
 Harris et al. (1997) M. J. Harris, S. T. Bramwell, D. F. McMorrow, T. Zeiske, and K. W. Godfrey, Phys. Rev. Lett. 79, 2554 (1997).
 Villain (1979) J. Villain, Zeitschrift Für Physik BCondensed Matter 33, 31 (1979).
 (41) Besides the trivial timereversal symmetry.
 Elhajal et al. (2002) M. Elhajal, B. Canals, and C. Lacroix, Phys. Rev. B 66, 014422 (2002).
 Canals et al. (2008) B. Canals, M. Elhajal, and C. Lacroix, Phys. Rev. B 78, 214431 (2008).
 Chern (2010) G.W. Chern, arXiv (2010), arXiv:1008.3038 [condmat.strel] .
 McClarty et al. (2014) P. A. McClarty, P. Stasiak, and M. J. P. Gingras, Phys. Rev. B 89, 024425 (2014).
 Jaubert et al. (2013) L. D. C. Jaubert, M. J. Harris, T. Fennell, R. G. Melko, S. T. Bramwell, and P. C. W. Holdsworth, Phys. Rev. X 3, 011014 (2013).
 Isakov et al. (2004) S. V. Isakov, K. S. Raman, R. Moessner, and S. L. Sondhi, Phys. Rev. B 70, 104418 (2004).
 Ryzhkin (2005) I. A. Ryzhkin, Journal of Experimental and Theoretical Physics 101, 481 (2005).
 Conlon and Chalker (2010) P. H. Conlon and J. T. Chalker, Phys. Rev. B 81 (2010).
 Isoda (2008) M. Isoda, Journal of Physics: Condensed Matter 20, 315202 (2008).
 Li et al. (2010) W. Li, S.S. Gong, Y. Zhao, S.J. Ran, S. Gao, and G. Su, Phys. Rev. B 82, 134434 (2010).
 Macdonald et al. (2011) A. J. Macdonald, P. C. W. Holdsworth, and R. G. Melko, Journal of PhysicsCondensed Matter 23, 164208 (2011).
 Hablützel (1939) J. Hablützel, Helv. Phys. Acta 12, 489 (1939).
 Vaks et al. (1966) V. Vaks, A. Larkin, and Y. N. Ovchinnikov, Soviet Physics JETP 22, 820 (1966).
 Cladis (1988) P. Cladis, Molecular Crystals and Liquid Crystals 165, 85 (1988).
 Landau et al. (1971) D. P. Landau, B. E. Keen, B. Schneider, and W. P. Wolf, Phys. Rev. B 3, 2310 (1971).
 Wolf et al. (1972) W. P. Wolf, B. Schneider, D. P. Landau, and B. E. Keen, Phys. Rev. B 5, 4472 (1972).
 Wolf (2000) W. P. Wolf, Brazilian Journal of Physics 30, 794 (2000).
 Struck et al. (2013) J. Struck, M. Weinberg, C. Ölschläger, P. Windpassinger, J. Simonet, K. Sengstock, R. Höppner, P. Hauke, A. Eckardt, M. Lewenstein, and L. Mathey, Nature Physics 9, 738 (2013).
 Wang et al. (2006) R. F. Wang, C. Nisoli, R. S. Freitas, J. Li, W. McConville, B. J. Cooley, M. S. Lund, N. Samarth, C. Leighton, V. H. Crespi, and P. Schiffer, Nature 439, 303 (2006).
 Nisoli et al. (2013) C. Nisoli, R. Moessner, and P. Schiffer, Rev. Mod. Phys. 85, 1473 (2013).
 Kapaklis et al. (2012) V. Kapaklis, U. B. Arnalds, A. HarmanClarke, E. T. Papaioannou, M. Karimipour, P. Korelis, A. Taroni, P. C. W. Holdsworth, S. T. Bramwell, and B. Hjorvarsson, New Journal of Physics 14 (2012), 10.1088/13672630/14/3/035009.
 Farhan et al. (2013a) A. Farhan, P. M. Derlet, A. Kleibert, A. Balan, R. V. Chopdekar, M. Wyss, L. Anghinolfi, F. Nolting, and L. J. Heyderman, Nature Physics 9, 1 (2013a).
 Farhan et al. (2013b) A. Farhan, P. M. Derlet, A. Kleibert, A. Balan, R. V. Chopdekar, M. Wyss, J. Perron, A. Scholl, F. Nolting, and L. J. Heyderman, Phys. Rev. Lett. 111, 057204 (2013b).
 Morgan et al. (2013) J. P. Morgan, J. Akerman, A. Stein, C. Phatak, R. M. L. Evans, S. Langridge, and C. H. Marrows, Phys. Rev. B 87, 024405 (2013).
 Marrows (2013) C. Marrows, Nature Physics 9, 324 (2013).
 Anghinolfi et al. (2015) L. Anghinolfi, H. Luetkens, J. Perron, M. G. Flokstra, O. Sendetskyi, A. Suter, T. Prokscha, P. M. Derlet, S. L. Lee, and L. J. Heyderman, Nature Communications 6, 8278 (2015).
 Arnalds et al. (2016) U. B. Arnalds, J. Chico, H. Stopfel, V. Kapaklis, O. Bärenbold, M. A. Verschuuren, U. Wolff, V. Neu, A. Bergman, and B. Hjörvarsson, New Journal of Physics 18, 023008 (2016).
 Zhang et al. (2012) S. Zhang, J. Li, I. Gilbert, J. Bartell, M. J. Erickson, Y. Pan, P. E. Lammert, C. Nisoli, K. K. Kohli, R. Misra, V. H. Crespi, N. Samarth, C. Leighton, and P. Schiffer, Phys. Rev. Lett. 109, 087201 (2012).
 Chioar et al. (2014) I. A. Chioar, N. Rougemaille, A. Grimm, O. Fruchart, E. Wagner, M. Hehn, D. Lacour, F. Montaigne, and B. Canals, Phys. Rev. B 90, 064411 (2014).
 Rousochatzakis et al. (2012) I. Rousochatzakis, A. M. Läuchli, and R. Moessner, Phys. Rev. B 85, 104415 (2012).
 Rojas et al. (2012) M. Rojas, O. Rojas, and S. M. de Souza, Phys. Rev. E 86, 051116 (2012).
 Apel and Everts (2011) W. Apel and H.U. Everts, Journal of Statistical Mechanics: Theory and Experiment 2011, P09002 (2011).
 Benton and Shannon (2015) O. Benton and N. Shannon, Journal of the Physical Society of Japan 84, 104710 (2015).
 Hiroi et al. (2001) Z. Hiroi, M. Hanawa, N. Kobayashi, M. Nohara, H. Takagi, Y. Kato, and M. Takigawa, Journal of the Physical Society of Japan 70, 3377 (2001).
 Okamoto et al. (2013) Y. Okamoto, G. J. Nilsen, J. P. Attfield, and Z. Hiroi, Phys. Rev. Lett. 110, 097203 (2013).
 Kimura et al. (2014) K. Kimura, S. Nakatsuji, and T. Kimura, Phys. Rev. B 90, 060414 (2014).
 Husimi (1950) K. Husimi, Journal of Chemical Physics 18, 682 (1950).
 Chandra and Doucot (1994) P. Chandra and B. Doucot, Journal of Physics AMathematical and General 27, 1541 (1994).
 Yoshida et al. (2002) S. Yoshida, K. Nemoto, and K. Wada, Journal of the Physical Society of Japan 71, 948 (2002).
 Strečka and Ekiz (2015) J. Strečka and C. Ekiz, Phys. Rev. E 91, 052143 (2015).
 Foini et al. (2013) L. Foini, D. Levis, M. Tarzia, and L. Cugliandolo, Journal of Statistical Mechanics , P02026 (2013).
 Levis et al. (2013) D. Levis, L. Cugliandolo, L. Foini, and M. Tarzia, Phys. Rev. Lett. 110, 207206 (2013).
 Pauling (1935) L. Pauling, Journal of the American Chemical Society 57, 2680 (1935).
 Landau and Binder (2009) D. P. Landau and K. Binder, A guide to Monte Carlo simulation in statistical physics (Cambridge University Press, Berlin, 2009).
 Onsager (1944) L. Onsager, Physical Review 65, 117 (1944).