Phase diagram of the 1/3-filled extended Hubbard model on the Kagome lattice
We study the phase diagram of the extended Hubbard model on the kagome lattice at 1/3 filling. By combining a configuration interaction approach to an unrestricted Hartree-Fock, we construct an effective Hamiltonian which takes the correlations back on top of the mean-field solution. We obtain a rich phase diagram with, in particular, the presence of two original phases. The first one consists of local polarized droplets of metal standing on the hexagons of the lattice, and an enlarged kagome charge order, inversely polarized, on the remaining sites. The second, obeying a local ice-rule type constraint on the triangles of the kagome lattice, is driven by an antiferromagnetically coupling of spins and is constituted of disconnected six-spin singlet rings. The nature and stability of these phases at large interactions are studied via trial wave functions and perturbation theory.
pacs:71.10.Hf, 73.20.Qt, 71.30.+h, 74.70.Kn
Lattice systems with geometrical frustration are at the center of intense theoretical and experimental research activity of current condensed matter, thanks to the extremely rich exotic phases encountered. One of the most important features of frustrated systems is the presence, in the classical limit, of highly degenerate ground states. These states, under the presence of quantum fluctuations, yield unconventional phases generic either in fermionic or bosonic systems such as spin liquidsMisguich (); Hermele (), valence bond physicsRokhsar (); Moessner (), pinball phasesPinball1 (); Pinball2 (), fractionalized defectsRuegg2 (), and even supersolidsRalko0 (); Trousselet (). Hence, frustration is at the crossroads of many disciplines, ranging from quantum antiferromagnets to optical lattices, via itinerant electrons and bosonic systems.
Intensively studied in quantum magnetism for exotic phases and spin liquid ground states, the role of the frustration is more and more considered in fermionic systems Pinball1 (); Pinball2 (); Ruegg2 (); Wen (); Pollmann (). Among the typical two-dimensional lattices possessing a geometrical frustration, one can mention the triangular, checkerboardRalko0 (); Wessel (), Cairo-pentagonalRalko3 (); Ioanis (), and kagome Ralko1 (); Ralko2 (); Isakov01 (); Indergrand02 (). The latter is the most frustrated two-dimensional lattice and, as a function of the electron filling, possesses several induced local constraint phases. At half filling, under strong Hubbard repulsion, several possible ground states have been proposed for the Heisenberg modelbook (), ranging from valence bond solidsMarston (); Balents (); Singh (); Evenbly (); Poilblanc () to various spin liquids such as a U(1) algebraicRan (), a triplet-gappedJiang01 () and a topological state with anyonic excitationsYan (); Depenbrock (); Jiang02 (). At intermediate interaction and for various fillings, original electronic features have recently been pointed out, such as peculiar Mott transitionOhashi (); Laeuchli (), anomalous quantum Hall effectsOhgushi (); LeHur (), Fermi surface instabilities Kiesel (), competing ordersWang (), and even anyonic states Ludwig (), to name a few.
Surprisingly, at lower filling, only a few works have been done so far. In particular, the -filled system takes a special role since, in the non-interacting limit, the Fermi sea is at the Dirac cones at the corner of the Brillouin zone (the K point) and the ground state (GS) is a semi-metal (SM). With interactions, let us cite studies on topological phases Wen (), Mott-insulator transitions for spinless fermions Nishimoto () or fractionalization of defects with an extra bow-tie term Ruegg2 (). Under the presence of the Hubbard interaction, spin charge density waves with antiferromagnetic spin stripes are stabilizedWen (), but this phase was enforced to preserve the translations. With nearest-neighbor interaction, a local constraint emerges, often called the ice-rule, very reminiscent of those observed in spin-ices. Thus, the ground state manifold is macroscopically degenerate with a finite entropy at zero temperature. This has to be connected with the physics of hard-core bosons which arises in corner-sharing frustrating unit lattices such as the checkerboard Ralko0 (); Indergrand02 (), the pyrochlore and the Cairo-pentagonal lattices Ralko3 (); Ioanis () for example.
On such a manifold, the quantum fluctuations select a resonating plaquette phase by an order-to-disorder mechanism through a perturbation theory of the third order similar to that of Resonating Valence Bond (RVB) physics, and the local ice-rule constraint enforces the system to have exactly one doubly occupied site per triangle Ruegg2 (). Despite the presence of all the ingredients –frustration, local constraint, ground state degeneracy, etc– needed to observe a very rich and exotic physics, the phase diagram of this model is still lacking in the literature, even at the mean-field level.
In this paper, we provide a theoretical study of the phase diagram of the extended Hubbard model at 1/3 filling at low and intermediate on-site and nearest neighbor interactions. In Sec. II are presented the model and methods, as well as the physical quantities used to determine the phases. In particular, we briefly introduce the configuration interaction (CI) method Malli (); Friedman (); Louis2 (); Li (); Bunge (); Sambataro () that we combine to an Unrestricted Hartree-Fock in order to partially bring back the correlations lost in the mean-field decoupling. Section III gives an overview of the thermodynamic limit phase diagram and a brief description of the encountered phases, emphasizing two original ones, the pinned metal droplet phase (PMD) and a spin charge density wave (SCDW). The identification of the phases, the extrapolation to the thermodynamic limit, and the stability of the phases upon the correlations by the CI are discussed in Sec. IV. Finally, the structure of the phases beyond the UHF and CI is discussed in Sec. V. We present a variational theory that helps in the understanding of the mechanism stabilizing the two peculiar phases, one being the presence of local polarized droplets of metal isolated by pins, and the other consisting of disconnected hexagons of six-spin singlets.
Ii Model and methods
The model considered in this paper is the extended Hubbard model on the kagome lattice:
where is a spinor of fermion (electron) creation operators at site , is the number of electrons with spin and is the total density. is the on-site Hubbard term, is the nearest neighbor Coulomb repulsion and is the hopping term. Throughout the paper, the lowest band is completely filled, corresponding to the 1/3 filling. As reported in the literature, at small interactions and , the system is a semi-metal (SM), while at large (keeping small) a resonating plaquette charge density wave (CDW) is stabilized Wen (); Ruegg2 (); Nishimoto (); Wen66 (); Pollmann (). Indeed, in this large limit, the system is governed by a local constraint enforcing each triangle of the latticeWen67 () to have exactly only one doubly occupied site and the others empty. This ensures the ground state, at the electrostatic limit (), to be macroscopically degenerate. Then, a ring exchange of particles on a hexagon of order lifts this degeneracy and selects the CDW. How those two phases behave upon intermediate values of is one of the questions addressed in this work.
The kagome lattice with periodic boundary conditions is defined by its linear size , and has three sites per unit cell, called , and in Fig.1. The number of sites is given by . Most of the computations of the ground states performed in this paper have been done by (i) using an Unrestricted Hartree-Fock (UHF) and (ii) partially bringing back on the top of this mean-field solution the correlations by using a Configuration Interaction on single and double excitations (SDCI) of the UHF Slater determinantMalli (); Friedman (); Louis2 (); Li (); Bunge (); Sambataro (). As a reminder, in the standard UHF approach, the interaction terms are decoupled as a Hartree term and a Fock term . The sets and , for and , are computed self-consistently, as well as the single electron basis wave functions. The configuration interaction (CI) method Malli (); Friedman (); Louis2 (); Li (); Bunge (); Sambataro (), then consists in constructing an effective Hamiltonian by projecting the microscopic Hamiltonian in a subspace made of Slater determinants , constructed in this single electron basis obtained by the UHF. If all Slater determinants are kept in the basis, then the full Hilbert space is spanned and is nothing else but . This is usually referred to as the full CI (FCI) in the literatureMalli (). Obviously, in this case, the limitations are the same as those for the exact diagonalizations (EDs), the Hilbert space is increasing exponentially, and only small cluster sizes are available. An advantage to working with the Slater basis however, is that a simple energy criterion can be exploited. Indeed, at small and intermediate values of and , the correlations are expected to be well taken into account by single and double excitations from the Fermi sea. Hence, it is reasonable to truncate the basis to those processes and to restrict our effective Hamiltonian to this subspace. This is the so called single-double CI (SDCI). Note that this has been intensively studied in quantum chemistry where it has been shown that, for certain molecules, up to 80% of the correlation energy was restored. Here, in addition to our results, we show in Sec. IV-B an improvement of about 60% for certain phases.
In order to characterize the phase diagram, we define a set of order parameters. The structure factor , where indicates spontaneous translational symmetry breaking of the ground state. Since at very large Coulomb repulsion the system is governed by the local ice-rule constraints, we define the local triangle operators with and . for a phase obeying an ice-rule of two particles per triangle with exactly one doubly occupied site, and if the ice-rule of the type of two particles with one empty site is realized. Both order parameters are zero in the homogeneous phase, and finite otherwise. To study the internal spin ordered phases, we define the operator , which is finite if there are no vacant sites.
Iii Overview of the phase diagram
The phase diagram of the 1/3-filled extended Hubbard model on the kagome lattice extrapolated to the thermodynamic limit from the UHF results up to site clusters is depicted in Fig.2.
Four phases have been obtained in the plane, the semi-metal (SM), a charge density wave (CDW), a spin charge density wave (SCDW), and a pinned metal droplet phase (PMD). As already mentioned, the first two phases have already been reported elsewhereNishimoto (); Wen66 (); Wen67 (). For the CDW, a perturbative effective model of order maps the present model onto a quantum dimer model on the hexagonal latticeNishimoto (); Wen66 (); Wen67 (); Pollmann (). Indeed, in the classical limit and , the ground state is macroscopically degenerate, and the configurations that minimize the energy obey a rule of two electrons per triangle with one doubly occupied site. A finite lifts this degeneracy and in the limit the system is effectively described by a hardcore quantum dimer model (QDM) on the honeycomb lattice Pollmann (); Nishimoto (); Wen66 () whose ground state consists of resonating plaquettes. The mean-field treatment cannot capture the resonating plaquettes of the QDM, but its solution can be thought of as its electrostatic counterpartWen (). This can be seen in the UHF local densities in Fig.4. A second order transition between the SM and the CDW is obtained. To the best of our knowledge, the two other phases reported in this work, namely the PMD and SCDW phases, have not been reported so farPollmann2 (). Starting from the CDW phase in the large regime, by increasing with , it is more favorable for the system to forbid double occupancies. This is happening through a first order phase transition. In contrast to the recent paper of Wen et al. [Wen, ], since we do not restrict our phases to preserve translational symmetries, we obtain a unit-cell phase with six-spin antiferromagnetic hexagonal rings. All spin-singlet rings are disconnected by empty sites, as depicted in Fig.4-a (see the text for more details). We call this phase the spin charge density wave (SCDW). Its phase transition with the SM is of second order. Eventually, at lower interaction, the antiferromagnetic rings are not favored anymore, and the system prefers to separate up and down spins via another first order transition. One spin sector crystallizes on one third of the sites in order to form an enlarged kagome structure, the other sector forming disconnected half-filled hexagonal metal rings on the remaining sites. This phase, called the pinned metal droplets (PMDs) phase, since the hexagons (metal droplets) are disconnected by localized electrons (pins), is stabilized thanks to the kinetic energy of those rings, as discussed in Sec. V. Except for the semi-metal, all the phases reported in this phase diagram break the translational symmetry, as a finite indicates (see Fig.3). We provide in Sec. IV the complete analysis about the nature of these phases by UHF, and then we describe the SDCI method and apply it to confirm their stability upon the correlations. We have followed the following strategy. First, we identify the phases by the UHF. Then, we check their stability deep in their domains and close to the phase boundaries by constructing a configuration interaction effective Hamiltonian in the single and double excitations (SDCI) Slater basis (see the next section for more details). We also verify that no new phases arise. Finally, in Sec. V, we unravel the mechanism of our phases by constructing a trial wave function and a perturbative theory based on our UFH and SDCI results at large and limit. This validates the present phase diagram.
Iv Identification of the phases and Thermodynamic Limit
iv.1 UHF thermodynamic limit
In this section, we first present the thermodynamic limit (TL) phase diagram, displayed in Fig.2, extracted from finite-size scalings of the UHF results. Then, we confirm the existence and the stability of all the phases beyond the mean-field level thanks to the single and double configuration interaction (SDCI). We have performed finite size scalings on the four order parameters, , , and on clusters up to sites. To illustrate our results, we consider the scan line displayed in Fig.2, which passes through all the phases. The results are displayed in Fig.3. As one can see, size effects are important (inset), as expected for the kagome lattice, but thanks to the large size clusters available, precise size scalings are performed. Our analysis reveals two second order phase transitions, between the SCDW and the SM and between the SM and the CDW. The finite order parameter (blue) at the TL indicates the spontaneous symmetry breaking of translations at the K point and thus a unit cell of nine sites. The quantity (gray) is zero as soon as on each triangle at least one zero magnetization site is present or finite otherwise. We observe that the PMD phase is the only one with a finite TL value, and all sites of each triangle are occupied. Finally, (red, green) provide insights about the internal charge (magnetization) structure of each triangle. In Table 1 is reported which order parameter is compatible with which phase.
All of them have been identified unambiguously within the UHF approach as seen in Fig.3. The corresponding density maps of the PMD and the SCDW, for the charge and the magnetization, are very close to the ones obtained by the SDCI (see below) depicted in Fig.4.
iv.2 Configuration interaction: beyond the Hartree-Fock
In order to get further insights about the nature of the phases and their stability beyond the UHF, we have developed a configuration interaction method which aims to bring back the correlations on top of a mean-field single particle solution, here obtained by the UHF (see Sec.II). In Table 2 are reported the ground state energies obtained by ED (FCI), SDCI and UHF on a cluster, the largest available by ED at the 1/3 filling. Keeping all single and double excitations above the Fermi sea corresponds to a total of 3489 Slater determinants. The parameters are chosen in the three non trivial phases, the CDW, the SCDW and the PMD. For each of them, the SDCI clearly improves the GS energies, hence showing that correlations have been incorporated.
In our model, if the phases obtained by UHF are not stable under the correlations (for example the UHF is known to overestimate ferromagnetism), we expect the SDCI to stabilize other solutions.
We were able to distinguish the four phases identified at the mean-field level. The snapshots of the local densities of the PMD and the SCDW phases are depicted in Fig.4. By turning on the correlations on top of the UHF solution, we found that all the phases are stable up to the largest cluster available by the CI, in our case the unit cells (108 sites). Note that for such a large cluster, the corresponding basis contains 263226945 determinants, which is not reachable by our computers. In order to reduce this size, it is possible to keep among the single and double excitations only those for which the excitation energy does not exceed a chosen cutoff . By adjusting correctly this cutoff, one can find a compromise between the improvement of the GS wave function and the size of the basis we want to deal with. For the site cluster, we have kept 1850220 Slater determinants for energy improvement of the order of on the GS energies. On the site cluster, for which no cutoff is required (97039 Slater determinants), the UHF energy is and for the SDCI in the SCDW phase at and [an improvement of ]. In the PMD phase at and , the UHF gives against for the SDCI . In Fig.4 are depicted the local charge and spin densities obtained by SDCI on the 108 site cluster, for the PMD and the SCDW phases. For information, in our computations, the change in the density values between the SDCI and the UHF never exceeded 3-4%, even close to the phase boundaries. This indicates that (i) these phases are stable upon the correlations, and (ii) the UHF phase diagram at the TL of Fig.2 captures the essential physics. Let us conclude this section by a remark concerning the phase boundaries of Fig.2. In this work, we have employed the SDCI to validate the presence of the main phases beyond the UHF and to verify whether or not new phases could be stabilized. Hence, we have focused deep in their domains, and made checking close to the main transitions. The question of the exact location of the boundaries within this approach requires an optimized methodology in selecting the active determinants, in order to reduce the size of the basis while keeping the correlations strongLi (); Bunge (); Sambataro (). This is left for another work.
V Trial wave function at large U and V
Let us now try to understand the characteristics of the PMD and the SCDW phases by considering the large and limit.
At large , we find a phase which is the separation between pins (one electron) on each site of an enlarged kagome lattice and balls of three polarized electrons on each hexagon, that we expect to be in a local metal state (delocalized electrons on six sites). Indeed, as a first approximation, one could think that these hexagons are isolated by the pins, as suggested by the real space density maps of the PMD depicted in Figs. 4(c) and 4(f). To validate this picture, we have calculated the energy of the trial wavefunction (TWF) where the first product is over single polarized electron states standing on the sites of the enlarged kagome lattice and the second is over all the remaining hexagons, being the ground state of an extended Hubbard model of three inversely polarized electrons on a ring:
By computing the ground state energy and the kinetic operator (the sites are labeled clockwise from zero to five), an inflexion point appears about . For , the ground state is delocalized (metal). At , the ground state on each hexagon is a simple metal of three electrons, and in terms of local momenta on the six-site ring, the wave function is trivially given by the product state of the single electron Bloch states , where , and , with energies and respectively. By applying the perturbation theory on at large , we find a correction to the energy as , , and for the pins , for a total trial energy . In Fig.5 are depicted the GS energies at obtained by UHF and the trial energy obtained from the TWF with metallic hexagons; their agreement is excellent.
We turn now tn the SCDW. In the large limit, at , there are no double occupancies and the system, for any finite , enters an ice-rule constraint of two electrons per triangle. Its ground state manifold is macroscopically degenerate and consists of closed loops of various lengths. When the quantum fluctuations are turned on, the UHF and the SDCI select the smallest possible loops as observed in Fig.4-d. From a perturbation theory, there are two types of processes which can lift this degeneracy: an antiferromagnetic coupling of spins of order at the second order, and resonating plaquette terms (ring exchange of three electrons on a hexagon) in at the third order, similar to those of the QDMNishimoto () at and . Unfortunately, neither the UHF (mean-field) nor the SDCI (up to double excitations) allow us to distinguish which process is dominant in functions of and . All we can say is that, at large enough and , the Heisenberg mechanism will be obviously favored. But let us go back to Eq. (2) and . The inflexion point in appearing around indicates that ring exchange terms become dominant for higher and thus the PMD is destroyed at such large values of . At the mean-field level, the UHF snapshot would be similar to that of Fig.4-d, but now the resonating hexagons are no longer the same as the Heisenberg loop 6 discussed above. Hence, two scenarios are possible at the transition of the PMD: (i) the system enters a phase with resonating plaquettes separated by pins based on our trial wave function, and then a further transition to a spin state driven by Heisenberg couplings is stabilized (with no pins), or (ii) the Heisenberg state is selected right away. Our UHF and SDCI results present a phase transition about up to , which is much smaller than , at which the internal structure of the hexagon is changing. This seems to indicate that scenario (ii) is favored, or in this case, based on the results of Fig.4-d, a ground state made of six-spin singlet hexagons of energy per site is selected among the degenerate configurations. This is what we call the spin charge density wave (SCDW). Finally, let us compare our results with previous work. In Wen et al. Wen, , since no translational symmetry breaking is allowed, a stripe phase with disconnected antiferromagnetic chains is obtained. In this case, the quantum fluctuations of the perturbation theory lowers the classical energy by . The hexagonal spin chain is then more likely when the symmetry constraint is removed, as confirmed by the snapshots of the SCDW depicted in Fig. 4(c) and 4(d).
Vi Concluding remarks
In this paper, we have studied an extended Hubbard model of spinful electrons on the kagome lattice at 1/3 filling. In addition to two already known phases, the semi-metal (SM) and the charge density wave (CDW), we have established the presence of two original phases, the pinned metal droplets (PMDs) and the spin charge density wave (SCDW) stabilized at intermediate and .
We have first used an unrestricted Hartree-Fock approach in order to extract the phase diagram at the TL. We have then constructed an effective configuration interaction Hamiltonian by projecting the microscopic model onto a basis of Slater determinants containing the UHF wave function, and all single and double excitations (SDCI). This allows for bringing back partial correlations on top of the UHF solution. For various sets of parameters, the ground-state energy (and related GS) is strongly improved (up to 60% of the correlation energy). We have confirmed that (i) no new phases are stabilized by the correlations and (ii) the two original phases, namely the PMD and the SCDW, are very robust against the correlations.
Finally, to get more insights about the nature of two new phases, we have performed a perturbation theory on a trial wave function. We have shown that the properties of the pinned metal droplets (PMDs) are perfectly captured by the UHF and explained by our trial wave function. It consists of polarized pins (electrons) standing on the sites of an enlarged kagome lattice, separating droplets of metal living on the remaining hexagons. The second important phase, the spin charge density wave (SCDW), is driven by Heisenberg interactions and leads to loop-6 ring isolated spin states. Note that we have discussed the possibility of a third phase that could arise in the system, close to the PMD but with resonating plaquettes instead of metal droplets on the remaining hexagons. In some sense, this would correspond to an internal crossover of the hexagon states, the pins remaining unchanged. However, within our approaches, this phase does not seem to be stabilized at the intermediate parameters considered in this work. Instead, the PMD and the six-spin singlet SCDW are the more stable.
We recently became aware of related work by F. Pollmann et. al (Ref. [Pollmann2, ]), where the authors have considered an effective model of Eq. (1) at large and limit. We thank F. Pollmann for stimulating discussions about their results corroborating the existence of the phases reported in this manuscript.
Acknowledgements – A.R. is grateful to F. Mila for insightful remarks about the stability of the phases. We thank S. Fratini for enlightening discussions and his careful rereading of the manuscript. A.R. thanks G. Bouzerar and M.-B. Lepetit for their valuable comments about the methods used in this paper. This work is supported by the French National Research Agency through Grants No. ANR-2010-BLANC-0406-0 NQPTP and No. ANR-12-JS04-0003-01 SUBRISSYME.
- (1) G. Misguich, D. Serban, and V. Pasquier, Phys. Rev. Lett. 89, 137202 (2002).
- (2) M. Hermele, Y. Ran, P. A. Lee, and X.-G. Wen, Phys. Rev. B77, 224413 (2008).
- (3) D.S. Rokhsar and S.A. Kivelson, Phys. Rev. Lett. 61, 2376 (1988).
- (4) R. Moessner and S.L. Sondhi, Phys. Rev. B63, 224401 (2001).
- (5) J. Merino, A. Ralko, and S. Fratini, Phys. Rev. Lett. 111, 126403 (2013).
- (6) F. Trousselet, A. Ralko, and A.M. Oleś, Phys. Rev. B86, 014432 (2012).
- (7) A. Ruëgg, and G. A. Fiete, Phys. Rev. B83, 165118 (2011).
- (8) A. Ralko, F. Trousselet and D. Poilblanc, Phys. Rev. Lett. 104, 127203 (2010).
- (9) F. Trousselet, P. Rudea-Fonseca, and A. Ralko, Phys. Rev. B89, 085104 (2014).
- (10) J. Wen, A. Rüegg, C.-C.J. Wang, and G.A. Fiete, Phys. Rev. B82, 075125 (2010).
- (11) F. Pollmann, P. Fulde, and K. Shtengel, Phys. Rev. Lett. 100, 136404 (2008).
- (12) S. Wessel, Phys. Rev. B86, 140501(R) (2012).
- (13) A. Ralko, Phys. Rev. B84, 184434 (2011).
- (14) I. Rousochatzakis, A. M. Laüchli, and R. Moessner, Phys. Rev. B85, 104415 (2012).
- (15) O. Cepas, and A. Ralko, Phys. Rev. B84, 020413(R) (2011).
- (16) D. Poilblanc, and A. Ralko, Phys. Rev. B82, 174424 (2010).
- (17) S. V. Isakov, S. Wessel, R. G. Melko, K. Sengupta, and Y. B. Kim Phys. Rev. Lett. 97, 147202 (2006).
- (18) M. Indergrand, A. Lauchli, S. Capponi, and M. Sigrist, Phys. Rev. B74, 064429 (2006).
- (19) P. Mendels and A. Wills, in Introduction to Frustrated Magnetism, edited by C. Lacroix, P. Mendels and F. Mila (Springer Berlin Heidelberg, 2011), vol. 164 of Springer Series in Solid-State Sciences, pp. 207-238, ISBN 978-3-642-10588-3, and references therein.
- (20) J. B. Marston, C. Zeng, J. Appl. Phys. 69, 5962 (1991).
- (21) L. Balents, Nature 464, 199 (2010).
- (22) R. R. P. Singh and D. A. Huse, Phys. Rev. B76, 180407 (2007).
- (23) G. Evenbly and G. Vidal, Phys. Rev. Lett. 104, 187203 (2010).
- (24) D. Poilblanc, M. Mambrini, D. Schwandt, Phys. Rev. B 81, 180402(R) (2010).
- (25) Y. Ran, M. Hermele, P. A. Lee, and X.-G. Wen, Phys. Rev. Lett. 98, 117205 (2007).
- (26) H. C. Jiang, Z. Y. Weng, and D. N. Sheng, Phys. Rev. Lett. 101,117203 (2008).
- (27) H.-C. Jiang, Z. Wang, and L. Balents, Nature Physics 8, 902 (2012).
- (28) S. Yan, D. A. Huse, and S. R. White, Science 332, 1173 (2011).
- (29) S. Depenbrock, I. P. McCulloch, and U. Schollwoeck, Phys. Rev. Lett. 109, 067201 (2012).
- (30) T. Ohashi, N. Kawakami, and H. Tsunetsugu, Phys. Rev. Lett. 97, 066401 (2006).
- (31) A. Läuchli and D. Poilblanc, Phys. Rev. Lett. 92, 236404 (2004).
- (32) K. Ohgushi, S. Murakami, and N. Nagaosa, Phys. Rev. B62, R6065(R) (2000).
- (33) A. Petrescu, A.A. Houck, and K. Le Hur, Phys. Rev. A86, 053804 (2012).
- (34) M.L. Kiesel, C. Platt, and R. Thomale, Phys. Rev. Lett. 110, 126405 (2013).
- (35) W.-S. Wang, Z.-Z. Li, Y.-Y. Xiang, and Q.-H. Wang, Phys. Rev. B87, 115135 (2013).
- (36) B. Bauer, L. Cincio, B. P. Keller, M. Dolfi, G. Vidal, S. Trebst, and A. W. W. Ludwig, cond-mat:1401.3017 (2014).
- (37) S. Nishimoto, M. Nakamura, A. O’Brien, and P. Fulde, Phys. Rev. Lett. 104, 196401 (2010).
- (38) Relativistic and Electron Correlation Effects in Molecules and Solids, Edited by G.L. Malli, Plenum Press, New York (1994).
- (39) B. Friedman and G. Levine, Phys. Rev. B55, 9558 (1997).
- (40) E. Louis, F. Guinea, M.P. López Sancho, and J.A. Vergés, Phys. Rev. B59, 14005 (1999).
- (41) S. Li, J. Ma, and Y. Jiang, Int. Nat. J. Quant. Chem. 78 153 (2000).
- (42) C.F. Bunge and R. Carbó-Dorca, J. Chem. Phys. 125 014108 (2006).
- (43) M. Sambataro, D. Gambacurta, and L. Lo Monaco, Phys. Rev. B83, 045102 (2011).
- (44) A. O’Brien, F. Pollmann, and P. Fulde, Phys. Rev. B81, 235115 (2010).
- (45) R. Moessner, S. L. Sondhi, and P. Chandra, Phys. Rev. B64, 144416 (2001).
- (46) F. Pollmann, K. Roychowdhury, C. Hotta, and K. Penc, arxiv:1402.4932 (2014).