The spin-1/2 Kagome XXZ model in a field:
between lattice nematic and solid orders
We study numerically the spin-1/2 XXZ model in a field on an infinite Kagome lattice. We use different algorithms based on infinite Projected Entangled Pair States (iPEPS) for this, namely: (i) an approach with simplex tensors and 9-site unit cell, and (ii) an approach based on coarse-graining three spins in the Kagome lattice and mapping it to a square-lattice model with local and nearest-neighbor interactions, with usual PEPS tensors, 6- and 12-site unit cells. Similarly to our previous calculation at the SU(2)-symmetric point (Heisenberg Hamiltonian), for any anisotropy from the Ising limit to the XY limit, we also observe the emergence of magnetization plateaus as a function of the magnetic field, at using 6- 9- and 12-site PEPS unit cells, and at and using a 9-site PEPS unit cell, the later set-up being able to accommodate solid order. We also find that, at , (lattice) nematic and VBC-order states are degenerate within the accuracy of the 9-site simplex-method, for all anisotropy. The 6- and 12-site coarse-grained PEPS methods produce almost-degenerate nematic and VBC-Solid orders. We also find that, within our accuracy, the 6-site coarse-grained PEPS method gives slightly lower energies, which can be explained by the larger amount of entanglement this approach can handle, even in the cases where the PEPS unit-cell is not commensurate with the expected ground state unit cell. Furthermore, we do not observe chiral spin liquid behaviors at and close to the XY point, as has been recently proposed. Our results are the first tensor network investigations of the XXZ model in a field, and reveal the subtle competition between nearby magnetic orders in numerical simulations of frustrated quantum antiferromagnets, as well as the delicate interplay between energy optimization and symmetry in tensor network numerical simulations.
pacs:03.65.Ud, 71.27.+a, 75.10.-b
Frustrated quantum antiferromagnets pose some of the hardest numerical challenges in condensed matter physics. This is so because, given frustration, quantum Monte Carlo is hampered by the infamous sign problem. In this context, the Kagome lattice is one of the most famous examples of frustrated lattices given the presence of triangles, where competing antiferromagnetic interactions between nearest-neighbour spins lead to many eigenstates with very different structure but with very similar and low energy, all of them competing to be the true ground state of the system. ÊThis makes the Kagome lattice an ideal experimental and theoretical setup to have very large quantum fluctuations at low temperatures, specially for spin-1/2 systems with Heisenberg-like interactions, but also for larger spins fru ().
An important model in this context is the the spin-1/2 XXZ model in a field on the Kagome lattice. In such model, an anisotropy angle is introduced that distinguishes the interactions between nearest neighbours from the and ones. For generic values of anisotropy and field the model has symmetry, and one recovers the full symmetry of the Heisenberg model at the isotropic point when the field is absent. The anisotropy allows us to study the crossover from the XY model () to the Ising model (). As such, the model exhibits frustration, and in the presence of a magnetic field it gives rise to quantized plateaus in the longitudinal (i.e. along the field) magnetization , being the most prominent at , Êbut also possible at and plateaus (). The existence, nature and properties of such plateaus have been a hot topic of discussion, with some works proposing that in the XY point (i.e. for purely and interactions) the plateau may give rise to a chiral topological quantum spin liquid chirxy (). In fact, even the existence of the plateau at the XY point is also a matter of controversy. Moreover, it is not clear if the plateau undergoes a phase transition as one tunes the anisotropy between the XY and the Ising points.
In this paper we study the zero-temperature phase diagram of the above model with several methods based on infinite Projected Entangled Pair States (iPEPS) ipeps (), similarly to previous studies of the Heisenberg point () heisipeps (). In particular, we use:
As we shall see, the use of different unit cells produces different results in the phase diagrams, depending on whether the magnetic structure of the phase at some plateaus is commensurate with the used cell. Moreover, we shall also see that methods for which the structure of some plateaus is not commensurate may produce, however, slightly lower variational energies because they are able to handle a larger amount of entanglement. In particular, we do a detailed analysis of this situation for the plateau. Our different methods reveal a subtle competition between symmetry and energy in the model: while the 9-site unit cell seems better-suited for the expected symmetry properties of the ground state in some parameter regimes such as the and plateaus, the approach based on the 6-site unit cell, even if not favourable from the symmetry perspective, turns out to be able to handle more entanglement in the wavefunction and may thus produce slightly better variational energies 111In fact, the property that symmetric tensor networks do not necessarily produce better variational energies is well-known for tensor network methodspess ().. We also observe a competition between different types of order: (i) Valence Bond Crystal solid (VBC-Solid) order breaking lattice translation symmetry and (ii) nematic order only breaking -rotation symmetry heisipeps0 (). More specifically, we see that the 9-site simplex-method produces degenerate nematic and VBC-Solid states 222Notice that in Ref.heisipeps () this order was called “Solid 1”, in order to distinguish it from other orders with a unit cell at other values of the magnetization., within our accuracy, for all values of the anisotropy up to the Ising point, whereas the 6- and 12-site coarse-grained PEPS methods produce almost degenerate nematic and VBC-Solid states, Êwith a slightly lower (but almost degenerate) energy than the 9-site approach (see Fig.1 for schematic diagrams of these orders). We also do not find any of the characteristic signatures of a chiral spin liquid phase in the XY point chirxy (). This paper is also the first tensor network study of the XXZ model in a field on a Kagome lattice.
Ii Model and methods
ii.1 The Hamiltonian
Ê Here we consider the antiferromagnetic XXZ model for spin-1/2 in the presence of an external magnetic field along the direction on the Kagome lattice. Its Hamiltonian is given by
where is the spin-1/2 th operator at site , is the magnetic field, is the anisotropy angle, and the interaction is for nearest-neighbour spins. For generic values of and this Hamiltonian is -symmetric under spin rotations around the -axis. Different interesting points can be accessed by tuning parameter , namely:
The Heisenberg point for .
The XY point for .
The Ising point for .
ÊAfter doing an overall analysis of the phase diagram, we shall focus on the incompressible phase at reduced magnetization where, in the classical Ising limit, a macroscopic degeneracy occurs between all configurations with two spin up and one spin down (for ) in every triangle wannier (). When a small XY (ferromagnetic or antiferromagnetic) anisotropy is added, i.e. , Eq. (1) can be mapped to an effective quantum dimer model on the dual hexagonal lattice cabra2005 () whose ground-state was argued to be a VBC-Solid moessner2001 (). This result is supported by direct Quantum Monte Carlo simulations on the (non-frustrated) hard-core boson model isakov2006 () equivalent to the ferromagnetic XXZ chain given by . Here we wish to investigate what happens away from the Ising point in this case.
ii.2 Numerical approaches
Following our approach in Ref.heisipeps (), we use imaginary-time evolution in order to obtain approximations of the ground state of the system in the thermodynamic limit. To implement this, here we use several approaches based on infinite-PEPS (or iPEPS), also similarly to what we did in Ref.heisipeps (). On the one hand, we use an approach based on introducing simplex tensors, also called “Projected Entangled Simplex States” (PESS) pess (), where we use a 9-site unit cell and the so-called “simple update” su (). On the other hand, we also use an approach where three spins in the Kagome lattice are mapped to a single coarse-grained 8-dimensional spin on the square lattice while keeping all interactions local, see Ref.heisipeps (), and then use a standard PEPS algorithm for the square lattice with 2-site and 4-site square-lattice unit cells, amounting to 6-site or 12-site unit cells in the original Kagome lattice. We have implemented this second approach both with the “simple update” and also the so-called “fast full update” ffu (), but we saw no significant difference in our results with these two update schemes, so we stick to the simple update because it is more efficient. The approach based on the coarse-grained Kagome lattice is less efficient but, however, has more variational parameters, Êwhich allows to handle a larger amount of entanglement in the ansatz wavefunction. For more specific details about these approaches, we refer the reader to Ref.heisipeps ().
In all cases, our refining parameter is the PEPS bond dimension , which controls the amount of entanglement in the tensor network, and which we consider up to . Depending on the choice of algorithm, translation symmetry may be broken in different ways, and also the number of variational parameters may be different. For instance, for the PESS approach with a 9-site unit cell, translation symmetry may be broken in, e.g., a superstructure. On the other hand, the coarse-grained lattice approach with 6-site and 12-site unit cells may break translation invariance in other ways, such as a structure. Concerning the number of variational parameters, the coarse-grained approach easily has many more than the PESS approach (e.g. a factor of 4 for and the 6-site unit cell), which means that the calculations need more CPU time, but they Êare also able to handle more entanglement in the wavefunction. As in Ref.heisipeps (), the characterization of the phases is possible by checking the local magnetizations as a function of the field, as well as the link energy terms , with as in Eq.(1). Expectation values and environment contractions are computed using Corner Transfer Matrices, see, e.g., Refs.ctm (); tn (). Moreover, in order to check the possible chirality of the wave-function, we compute the entanglement spectra of the state wrapped around half an infinite cylinder, by following the procedure in, e.g, Ref.chirpeps ().
In Fig.2 and Fig.3 we show the different plateaus that arise as a function of the external magnetic field and the anisotropy angle , as computed with the 9-site PESS and the 6-site coarse-grained PEPS. We find that the 9-site PESS reproduces different plateaus at and , whereas the 6-site coarse-grained PEPS reproduces only the , due to the incommensurability of this unit cell with the expected magnetic structure of the other plateaus. Concerning the ground state energy, the approach based on the 6-site cell seems to produce slightly lower values than the 9-site cell. For instance, at the plateau region for we find that for , the ground state energies per site are (where in fact there is no plateau) and . This slight energetic advantage may be due to the larger number of variational parameters involved in the 6-site wavefunction for a given . However, variations beyond the third significant digit should be taken with care, since the relative difference is still very small specially in the context of the several approximations involved in both algorithms, and which we control to the best of our possibilities.
For concreteness, in Fig.4 we show the magnetization curve, i.e. as a function of the external field obtained with the 6-site algorithm, for different values of the anisotropy angle . We observe the presence of a very prominent magnetization plateau at , for all the scanned values of , including the XY, Heisenberg and Ising points. Away from the Ising point we also see a jump towards the saturation value Êfor several values of theta jump () – a more fine-grained analysis (not displayed) shows that this is indeed the case –. The width of the plateau increases as we tune from the XY towards the Ising point, see Fig.5. Results for the 12-site algorithm are essentially identical. Notice also that the 9-site PESS algorithm also finds a plateau at the XY point. 333Note that it was suggested in Ref.cabra2005 () that the plateau width vanishes at .
In Fig.6 we show the ground state energy per site of in the plateau phase as a function of the anisotropy angle , as computed with our three methods. We see that the 6-site approach is the one that produces a lower variational energy. The difference with the 9-site PESS is around the third digit, roughly 1%, whereas the difference with the 12-site PEPS is even smaller, around 0.1%. As said before, these small differences should be taken with caution, given the different approximations involved in the methods used. The fact that the 9-site PESS approach produces higher energy may be a consequence of the smaller number of variational parameters in this ansatz. The quasi-degeneracy between the 6-site and the 12-site PEPS energies shows that a doubling of the unit cell is preferred over a quadrupling. We have also observed that the 9-site PESS produces degenerate nematic and VBC-Solid states (within our energy resolution), similar to what we found in Ref.heisipeps () at the Heisenberg point. This degeneracy is very slightly lifted in the 6- and 12-site PEPS approaches, in favor of the nematic order for XY anisotropy, in the vicinity of the Heisenberg point (with a typically relative energy difference of at the Heisenberg point) and for a small Ising anisotropy, and in favor of VBC-Solid order close to the Ising point. This is true even for the largest bond dimensions that we could reach. It may be possible that the energy difference between these states becomes narrower in the large-size limit, thus explaining this observation. The stability of one phase w.r.t. the other is assessed when systematically obtaining the same one in different runs of our algorithms with different initial conditions.
We next perform an analysis of the plateau Êwith the different algorithms. In Fig.7 we show diagrammatically the longitudinal magnetizations at every site as well as the link energy terms at every different link, for the 6-site PEPS (left) and the 9-site PESS (right) calculations, for points in the plateau and four different anisotropy angles . In our simulations we have seen that, for the 6-site PEPS approach, for the phase corresponds, according to the classification in Ref.heisipeps (), to a nematic state, where rotation symmetry is broken down to (3-fold degenerate). For we obtain a VBC-Solid phase (again according to the classification in Ref.heisipeps ()), where translation symmetry is broken down to a cell, and the discrete symmetries (lattice 6-fold rotation) and (reflexion about one axis) are fully broken. Note that the Heisenberg point lies then within the nematic region. Comparison with the results obtained from the 9-site PESS approach is useful. In the right panel of Fig. 7 we show different orders obtained in this way. In practice we have seen that, for every possible , nematic and VBC-solid orders are degenerate in energy within our accuracy. The structures in the figure are examples of what we obtain in some runs of the algorithm. By running the algorithm with different random initial conditions, we get either a nematic or VBC-Solid state, both with the same energy. Such a degeneracy between nematic and VBC-Solid seems to be lifted slightly in the 6- and 12-site PEPS approaches.
Finally, in order to see possible signs of chirality, we have computed the entanglement spectrum of the PEPS obtained from the 6-site PESS, for half an infinite cylinder of width 12, within the plateau at the XY point. Such state has been conjectured to be chiral chirxy (). In our simulation, however, the entanglement spectrum is perfectly symmetric under time-reversal and therefore not chiral, see Fig.8 444ÊAs a word of caution, it is well-known that finding a gapped chiral topological state with PEPS may actually be quite difficult because of fundamental constraints of the PEPS tensor network, and therefore it is also difficult to get an intuition about from a numerical optimization.
Ê Here we have studied the zero-temperature phase diagram of the Kagome XXZ model in a field and in the thermodynamic limit using different tensor network approaches. We find different plateau structures for the longitudinal magnetization as a function of the field and the anisotropy, depending on the unit cell being used. We also find that a 6-site coarse-grained PEPS approach seems to produce slightly lower ground state energies, probably because the larger number of available parameters which allows for more entanglement in the ansatz wavefunction, and in spite of the fact that the unit cell cannot accommodate the superstructures found in other recent work huerga2016 (). Such small energy differences should be considered with caution, though, given the different approximation schemes involved in each method. From our results, we cannot discern if our 6-site ansatz, which contains more adjustable parameters, slightly bias the results towards a 6-site superstructure, or if this superstructure is indeed a genuine feature of the system. In any case, our results show a tight competition between nematic and VBC-Solid orders, which we find for all values of the anisotropy, Êas well as the delicate interplay between the implementation of (lattice) symmetries and the optimization of the energy in tensor network algorithms.
Acknowledgements.A. K. and R. O. acknowledge financial support from JGU, as well as the University of Toulouse and the CCBPP, where part of this work was performed. Discussions with Y. Iqbal and M. Ziegler are also acknowledged.
- (1) S. Yan, D. A. Huse, and S. R. White, Science, 332 (6034):1173-1176 (2011); S. Depenbrock, I. P. McCulloch, and U. Schollwöck, Phys. Rev. Lett. 109, 067201 (2012); Y. Iqbal, F. Becca, S. Sorella, and D. Poilblanc, Phys. Rev. B 87, 060405(R) (2013); D. Poilblanc and N. Schuch, Phys. Rev. B 87, 140407 (2013); D. Poilblanc, N. Schuch, D. Pérez-García, and J. I. Cirac, Phys. Rev. B 86, 014404 (2012); K. Hida, J. Phys. Soc. Jpn, 70:3673 (2001); A. Honecker, J. Schulenburg, and J. Richter, J. Phys.: Condens. Matter 16, S749-S758 (2004); Y. Okamoto, M. Tokunaga, H. Yoshida, A. Matsuo, K. Kindo, and Z. Hiroi, Phys. Rev. B 83, 180407(R) (2011); Krishna Kumar, Kai Sun, and Eduardo Fradkin, Phys. Rev. B 92, 094433 (2015); K. Hida, J. Phys. Soc. Jpn. 69(12):4003-4007, (2000); T. Liu, W. Li, A. Weichselbaum, J. von Delft, and G. Su, Phys. Rev. B 91, 060403(R) (2015); H. J. Changlani and A. M. Läuchli, Phys. Rev. B 91, 100407 (2015); W. Li, A. Weichselbaum, J. von Delft, and H.-Hao Tu, Phys. Rev. B 91, 224414 (2015); D. Ixert, T. Tischler, and K. P. Schmidt, Phys. Rev. B 92, 174422 (2015); O. Götze, D. J. J. Farnell, R. F. Bishop, P. H. Y. Li, and J. Richter, Phys. Rev. B 84, 224428 (2011); W. Li, S. Yang, M. Cheng, Z.-Xin Liu, and H.-Hao Tu, Phys. Rev. B 89, 174411 (2014); S. Nishimoto and M. Nakamura, Phys. Rev. B 92, 140412(R) (2015); K. Okuta, S. Hara, H. Sato, Y. Narumi, and K. Kindo, J. Phys. Soc. Jpn. 80(6):063703 (2011); T. Heng Han, J. S. Helton, S. Chu, D. G. Nocera, J. A. Rodriguez-Rivera, C. Broholm, and Y. S. Lee, Nature, 492: 406-410 (2012); D. E. Freedman, R. Chisnell, T. M. McQueen, Y. S. Lee, C. Payen, and D. G. Nocera, Chem. Commun., 48:64-66 (2012); S. Hara, H. Sato, Y. Narumi, J. Phys. Soc. Jpn. 81(7):073707 (2012). J. B. Marston and C. Zeng, J. Appl. Phys. 69 (8):5962-5964, (1991); K. Yang, L. K. Warman, and S. M. Girvin, Phys. Rev. Lett. 70, 2641 (1993); M. B. Hastings, Phys. Rev. B 63, 014413 (2000); F. Wang and A. Vishwanath, Phys. Rev. B 74, 174423 (2006); Y. Ran, M. Hermele, P. A. Lee, and Xiao-Gang Wen, Phys. Rev. Lett. 98, 117205 (2007); R. R. P. Singh and D. A. Huse, Phys. Rev. B 77, 144415 (2008); G. Evenbly and G. Vidal, Phys. Rev. Lett. 104, 187203 (2010). J.-Wei Mei, J.-Yao Chen, H. He, X.-Gang Wen, arXiv:1606.09639.
- (2) S. Capponi, O. Derzhko, A. Honecker, A. M. Lauchli, and J. Richter, Phys. Rev. B 88, 144416 (2013); S. Nishimoto, N. Shibata, and C. Hotta, Nat. Commun. 4, 2287 (2013).
- (3) K. Kumar, K. Sun, and Eduardo Fradkin, Phys. Rev. B 90, 174409 (2014).
- (4) J. Jordan, R. Orús, and G. Vidal, F. Verstraete, and J. I. Cirac, Phys. Rev. Lett. 101, 250602 (2008); J. Jordan, R. Orús, and G. Vidal, Phys. Rev. B 79, 174515 (2009).
- (5) T. Picot, M. Ziegler, R. Orús and D. Poilblanc, Phys. Rev. B 93, 060407 (2016).
- (6) Z. Y. Xie, J. Chen, J. F. Yu, X. Kong, B. Normand, and T. Xiang, Phys. Rev. X 4, 011025 (2014); H. J. Liao, Z. Y. Xie, J. Chen, Z. Y. Liu, H. D. Xie, R. Z. Huang, B. Normand, and T. Xiang, arXiv:1610.04727.
- (7) Thibaut Picot and Didier Poilblanc, Phys. Rev. B 91, 064415 (2015).
- (8) Ê G. H. Wannier, Phys. Rev. 79, 357 (1950).
- (9) D. C. Cabra, M. D. Grynberg, P. C. W. Holdsworth, A. Honecker, P. Pujol, J. Richter, D. Schmalfuß, and J. Schulenburg, Phys. Rev. B 71, 144420 (2005).
- (10) R. Moessner, S.L. Sondhi, P. Chandra, Phys. Rev. B 64, 144416 (2001).
- (11) S. V. Isakov, S. Wessel, R. G. Melko, K. Sengupta, Yong Baek Kim, Phys. Rev. Lett. 97, 147202 (2006).
- (12) H. C. Jiang, Z. Y. Weng, and T. Xiang, Phys. Rev. Lett. 101, 090603 (2008).
- (13) H. N. Phien, J. A. Bengua, H. D. Tuan, P. Corboz, and R. Orús, Phys. Rev. B 92, 035142 (2015).
- (14) R. Orús and G. Vidal, Phys. Rev. B 80, 094403 (2009);
- (15) F. Verstraete, J. I. Cirac, and V. Murg, Adv. Phys. 57,143 (2008); J. I. Cirac and F. Verstraete, J. Phys. A: Math. Theor. 42, 504004 (2009); R. Augusiak, F. M. Cucchietti, and M. Lewenstein, in Modern Theories of Many-Particle Systems in Condensed Matter Physics, Lect. Not. Phys. 843, 245-294 (2012); J. Eisert, Modeling and Simulation 3, 520 (2013); N. Schuch, QIP, Lecture Notes of the 44th IFF Spring School (2013); R. Orús, Eur. Phys. J. B 87, 280 (2014); R. Orús, Ann. Phys.-New York 349 117158 (2014).
- (16) J. I. Cirac, D. Poilblanc, N. Schuch, and F. Verstraete, Phys. Rev. B 83, 245134 (2011).
- (17) J. Schulenburg, A. Honecker, J. Schnack, J. Richter, and H.-J. Schmidt, Phys. Rev. Lett. 88, 167207 (2002).
- (18) R. Orús, T.-Chieh Wei, O. Buerschaper, A. García-Saez, Phys. Rev. Lett. 113, 257202 (2014).
- (19) Daniel Huerga, Sylvain Capponi, Jorge Dukelsky, Gerardo Ortiz, Phys. Rev. B 94, 165124 (2016).