# Frustrating antiferromagnetic exchange interactions enhance specific valence-bond-pair motifs

###### Abstract

We present variational results for the ground state of the antiferromagnetic quantum Heisenberg model with frustrating next-nearest-neighbour interactions. The trial wave functions employed are of resonating-valence-bond type, elaborated to account for various geometric motifs of adjacent bond pairs. The calculation is specialized to a square-lattice cluster consisting of just sixteen sites, large enough that the system can accommodate nontrivial singlet dimer correlations but small enough that exhaustive enumeration of states in the total spin zero sector is still feasible. A symbolic computation approach allows us to generate an algebraic expression for the expectation value of any observable and hence to carry out the energy optimization exactly. While we have no measurements that could unambiguously identify a spin liquid state in the controversial region at intermediate frustration, we can say that the bond-bond correlation factors that emerge do not appear to be consistent with the existence of a columnar valence bond crystal. Furthermore, our results suggest that the magnetically disordered region may accommodate two distinct phases.

FrustrationDiep05 (); Lacroix11 (); Balents10 () is believed to be a key ingredient for stabilizing magnetically disordered statesAffleck87 (); Affleck88 (); Read89 () in quantum spin systems. One of the canonical and most widely discussed frustrated models is the so-called –, a spin-half quantum Heisenberg hamiltonian with nearest- and next-nearest-neighbour antiferromagnetic exchange interactions.Chandra88 () When defined on a bipartite lattice, the model possesses a weakly frustrated limit, in the vicinity of , in which the ground state has a well-defined Marshall sign structureMarshall55 (); Richter94 (); Beach09 (); Richter10 (); Zhang13 () and a tendency toward antiferromagnetic order (scrambled by quantum fluctuations in one dimensionMermin66 (); Hohenberg67 () but generally robust for higher-dimensional latticesReger89 (); Weihond91 (); Beach07a ()). As the relative coupling is tuned up from zero, the effects of the frustration become increasingly disruptive, to the point where they induce a quantum phase transition.

The – model on the linear chain has been extensively studied, and its behaviour over the full range of relative couplings is known.White96 () A crucial anchor for our understanding is the famous Majumdar-Ghosh point,Majumdar69 () at , where the ground state is a perfect spin-Peierls product state,Frahm97 () i.e., a fluctuationless crystal of dimerized singlet pairs or valence bonds.Rumer32 (); Pauling33 (); Hulthen38 () What the one-dimensional model clearly demonstrates is that frustration can cause a valence bond state that otherwise resonatesAnderson73 (); Fazekas74 (); Sutherland88 (); Liang88 () to freeze into a static bond pattern that breaks translational symmetry.

In two dimensions, however, there is no comparable, exactly solvable point and no compelling reason to believe that bond crystallization must occur. There do exist models with farther-neighbour interactions for which the columnar bond crystal is a genuine ground stateBose91 () (and others with staggered and herring bone dimer statesBose92 ()), but it is unclear whether such models are in any sense “close” to the – in parameter space. Indeed, Jiang, Yao, and Balents argue (forcefully, in Ref. Jiang11, ) that the columnar valence bond crystal is not a good candidate for the intermediate phase of the square-lattice version of the model. Rather, their density matrix renormalization group calculations suggest a spin liquid with robust gaps to the and excitations. Hu and coworkers have also presented evidence for a liquid state, but a gapless one.Hu13 () Their approach is to show that a Gutzwiller-projected fermionic trial wavefunction subject to additional Lanczos steps produces extremely good variational energies in the frustrated regime. We note that these two claims for a liquid state are just the latest salvos in the long-running dispute as to whetherGelfand89 (); Gelfand90 (); Singh90 (); Zhitomirsky96 (); Leung96 (); Kotov99 (); Singh99 (); Kotov00 (); Capriotti00 (); Saarloos00 (); Takano03 (); Mambrini06 (); Murg09 (); Reuther10 (); Reuther11 (); Yu12 () or notFigueirido90 (); Oguchi90 (); Locher90 (); Schulz92 (); Zhong93 (); Oitma96 (); Capriotti01 (); Zhang03 (); Capriotti03 (); Captiotti04a (); Capriotti04b (); Yunoki04 (); Wang11 (); Mezzacapo12 (); Li12 (); Wang13 () the ground state has broken translational symmetry.

This paper attempts to address the issue in the following, modest way. We posit that if a bond crystal is favoured over a spin liquid, there should be a signature in the form of explicit correlations between pairs of valence bonds—presumably, enhanced (suppressed) correlations between pairs whose arrangement is compatible (incompatible) with the global crystalline order. This is a reasonable assumption, since we know that an unbiased short-bond-only state has at most powerlaw dimer correlations.Albuquerque10 (); Tang11 (); Ju12 () Our approach is to construct expressive trial wavefunctions that include various pair correlation factors and then to look at the behaviour of the optimized variational parameters that emerge at a given level of frustration.

Specifically, we consider a resonating-valence-bond (RVB) wave function of the form

(1) |

The summation ranges over all singlet-product states that are purely bipartite.Beach06 (); Beach08 () On the lattice with periodic boundary conditions—the case to which we specialize—only two kinds of bond are possible: those that are nearest-neighbour [whose end-to-end vector is symmetry equivalent to ] and those that entangle spins a knight’s moveIske87 () apart [].

In Eq. (1), we allow for up to three continuous variational parameters, denoted , , and , each of which appears as a factor raised to an integer power. The exponent is the number of knight’s move bonds in , whereas and count particular motifs of two adjacent bonds, chosen from the six possibilities shown in Fig. 1. The motifs are equivalent under all lattice symmetries and thus depend only on the relative orientation of the two bonds. The counting is consistent with the convention for the correlated valence bond states established in Ref. Lin12, . We point out that the usual Liang-Doucot-Anderson amplitude product stateLiang88 () is recovered from Eq. (1) with the choice of values and . In accordance with the standard notation, refers to the independent amplitude for a single bond connecting sites and in opposite sublattices.note ()

We argue that, despite its small size, the lattice system provides a good caricature of the physics of frustration. The lattice is large enough to accommodate all possible pairs of adjacent nearest-neighbour bonds. Equally important, it accommodates the one crucial long bond that is responsible for the breakdown of the Marshall sign structure: we know that, initially at least, the knight’s move bond is the only one to turn negativeRichter94 (); Beach09 (); Zhang13 () in response to increasing frustration.

Computational approach. As two of us have pointed out elsewhere,Zhang13 () the system can be solved in symbolic form with the help of the Lehmer code.Lehmer60 (); Knuth73 () All observables

(2) |

can be computed in the form of a ratio of two mixed polynomials in , , and . For a given value of the coupling strength , the triplet of values is chosen to minimize the energy

(3) |

In the notation of Eq. (3), each configuration is a double dimer covering , and is the number of loops formed by the overlap of states and . The exponent is defined as , and there are corresponding definitions for and . The factors of 2 are loop fugacities that arise from overlapsSutherland88 () in the overcomplete basis, and is the loop estimatorBeach06 () of the hamiltonian evaluated for configuration .

The energy minimization is carried out by steepest descent (and results verified afterwards by a hierarchical global search). All that is required is knowledge of the local downhill direction

(4) |

where the gradient is expressed in the variational parameter coordinates. Since both and have rational polynomial form, the downhill vector does too; hence Eq. (4) is straightforward to compute, so long as we take care not to overrun the computer’s finite precision. Our strategy to cope with floating-point difficulties is two-fold: first, we employ the 80-bit Extended Precision Format provided in the x86 architecture; second, we ensure that evaluation of any function

(5) |

is treated as a triply nested polynomial

(6) |

with the floating-point operations carried out in the order dictated by Horner’s rule.Borwein95 () The components of the energy gradient can be cast into the same nested form:

(7) |

Because the location of the energy minimum evolves smoothly as the relative coupling is varied, it is easy to keep the computation under good numerical control as we sweep from low frustration to high. We simply carry out the variational calculation repeatedly in a sequence of smalls steps , at each stage seeding the search with the previous step’s results.

Numerical results. The uppermost panel of Fig. 2 reports the variational energy density of the single-motif and pure RVB (i.e., no-motif) wavefunctions. Treating the latter amounts to minimizing over alone, leaving fixed. In the remaining cases, we hold fixed and minimize over and with counting each of the single motifs in turn. We find that the addition of even one bond correlation factor can substantially improve the variational energy, especially for large frustration.

The middle panel of Fig. 2 shows the discrepancy between the energy density of the exact ground state and that of the best variational state determined at the current value of the relative exchange coupling. Several distinct regimes are revealed. In the range , the energy is lowest when the herring bone motif is slightly suppressed, its weight falling from 0.99 to 0.93 with increasing . For , the best energy is given by the staggered motif, whose weight climbs from 1.08 to 1.20. Right at the level-crossing point, , the three wavefunctions with herring bone, staggered, and line motifs are all energy degenerate.

Figure 3 illustrates the further improvement to the trial state when two motifs are taken into account. For , the best variational energy is achieved by strong suppression of the line () and herring bone () motifs. In the thin sliver , the best energy corresponds to a mild enhancement of the plaquette () and staggered () motifs. The first of the two ranges we have just discussed is situated in what would be the antiferromagnetic phase of the infinite system; the second is likely centred on the critical point where the Néel order is extinguished. Note that, in contrast to the results in Fig. 2, two optimization regimes rather than one are needed to span the magnetically disordered phase, which is believed to extend up to .Zhang13 () For , we find reduced plaquette () and increased line () correlations. But deeper in the disordered phase, over a range that abuts the antiferromagnet, the optimal state involves reductions of the plaquette () and herring bone () weights.

We now make some additional observations about the data in the bottom panels of Figs. 2 and 3. First, we highlight the fact that the behaviour of the optimized parameter depends only very weakly on the form of the variational wavefunction. Without exception, descends monotonically along a nearly consistent line of values and goes negative in the vicinity of , which confirms our previous understanding that the Marshall rule breaks down simply by way of changing sign. This appears to be a robust feature that survives the addition of any number of correlation factors. The second observation is that the identification of correlation trends is not quite as clean as we had anticipated. In the bottom panel of Fig. 3, the inconsistent behaviour of shows how the plaquette correlation can be favoured or disfavoured depending on which other bond-bond correlation factor it is combined with (in this case, each of the staggered, herring bond, and line). There is, however, consistency in the sense that the wavefunctions having the lowest overall energy do agree in predicting decreased plaquette correlations throughout the magnetically disordered phase.

Conclusions. Our variational results indicate that the ground state of the – model is very well approximated by a simple RVB state all the way up to . Beyond that level of frustration, explicit bond-bond correlations become increasingly important, and it is unlikely that any pure RVB descriptionLou07 (); Beach09 (); Li12 (); Zhang13 () of the disordered phase will be reliable.

Indeed we know of several groups that have anticipated this deficiency and are considering how to carry out calculations with beyond-pure-RVB states, particularly in the tensor network frameworkWang13 () (which, unlike stochastic sampling methods, permits evaluation of states with negative amplitudes).

Based on the small-scale (but exact) results reported here, we can offer some insights to those attempting future large-scale, variational calculations using correlated valence bond states. (i) The ramp and cross motifs are never competitive energetically in our scheme, so it may be sufficient to account for correlations between short bonds only. (ii) The sign of each configurational amplitude can be attributed to the individual bond amplitudes and probably even confined to just ; this means that correlation factors can be taken to be positive definite. (iii) In our system, the optimized values of the correlation factors are unique, but this may not be true generally. On larger lattices and when more than a few motifs are employed, there could be internal redundancy in the full set of variational parameters. (iv) The correlation weight trends do not always depend in an obvious way on the combination of motifs.

As for the nature of the intermediate phase, the correlations do not indicate a particular bond order (columnarKotov00 () or plaquette,Mambrini06 () say) in a straightforward way, but they do seem to be at odds with the columnar pattern. For instance, in the range we find that the weights associated with pairs of bonds in plaquette and line arrangements move in opposite directions (decreasing and increasing, respectively), when we would expect them both to increase if the state were tending toward a columnar crystal. (What this may point to instead is a resonating state with enhanced quasi-one-dimensional character, as suggested in Refs. Beach09, and Zhang13, .) Over the full range of the disordered phase , the plaquette correlation is disfavoured rather significantly.

Given the level crossing that appears at , we might also entertain the possibility that there are multiple disordered phasesSushkov01 () that intervene between the and Néel states. For example, a reasonable scenario (one consistent with our results) is that the leftmost phase transition is pinned at 0.4073, where the Marshall sign rule first breaks down; there, the antiferromagnet is superseded by a liquid phase that survives up to 0.5401, where a non-columnar bond crystalSaarloos00 (); Mambrini06 () then becomes stable; this is turn survives up to the onset of the antiferromagnet at 0.60. An interesting coincidence is that is the point up to which Néel order survives in the variational RVB calculation with no correlation factors and all individual bond amplitudes restricted to .Zhang13 () This may be a clue that the state in the range depends in some crucial way on the nontrivial sign structure, whereas the state in the range does not.

Acknowledgements. The authors offer their thanks to Anders Sandvik, Ling Wang, and Olexi Montunich for several useful discussions. Funding was provided by a Discovery grant from NSERC of Canada. KSDB also acknowledges the financial support and hospitality of the Department of Physics and IQIM at CalTech.

## References

- (1) Frustrated spin systems, edited by H. T. Diep editor, (World-Scientific, Singapore, 2005). ISBN 978-9812560919
- (2) Introduction to Frustrated Magnetism: Materials, Experiments, Theory, Springer Series in Solid-State Sciences, Vol. 164, edited by C. Lacroix, P. Mendels, and F. Mila (Springer, Berlin, 2011). ISBN 978-3642105883
- (3) L. Balents, Nature Phys. 464, 199 (2010).
- (4) I. Affleck, T. Kennedy, E. H. Lieb, and H. Tasaki, Phys. Rev. Lett. 59, 799 (1987); T. Kennedy, E. H. Lieb, H. Tasaki, J. Stat. Phys. 53, 383 (1988).
- (5) I. Affleck and J. B. Marston, Phys. Rev. B 37, 3774 (1988).
- (6) N. Read and S. Sachdev, Phys. Rev. Lett. 62, 1694 (1989).
- (7) P. Chandra and B. Doucot, Phys. Rev. B 38, 9335 (1988).
- (8) W. Marshall, Proc. R. Soc. Lond. A 232, 48 (1955).
- (9) J. Richter, N. B. Ivanov, and K. Retzlaff, Europhys. Lett. 25, 545 (1994).
- (10) K. S. D. Beach, Phys. Rev. B 79, 224431 (2009).
- (11) J. Richter and J. Schulenburg, Eur. Phys. J. B 73, 117 (2010).
- (12) X. Zhang and K. S. D. Beach, Phys. Rev. B 87, 094420 (2013).
- (13) N. D. Mermin and H. Wagner, Phys. Rev. Lett. 17, 1133 (1966).
- (14) P. C. Hohenberg, Phys. Rev. 158, 383 (1967).
- (15) J. D. Reger, J. A. Riera, and A. P. Young, J. Phys.: Condens. Matter 1, 1855 (1989).
- (16) Z. Weihong, J. Oitmaa, and C. J. Hamer, Phys. Rev. B 44, 11869 (1991).
- (17) K. S. D. Beach and A. W. Sandvik, Phys. Rev. Lett. 99, 047202 (2007).
- (18) S. R. White and I. Affleck, Phys. Rev. B 54, 9862 (1996).
- (19) C. K. Majumdar and D. Ghosh, J. Math. Phys. 10, 1388 (1969); C. K. Majumdar, J. Phys. C: Solid State Phys. 3, 911 (1970).
- (20) H. Frahm and J. Schliemann, Phys. Rev. B 56, 5359 (1997).
- (21) G. Rumer, Göttingen Nachr. Tech. 1932, 377 (1932).
- (22) L. Pauling, J. Chem. Phys. 1, 280 (1933).
- (23) L. Hulthén, Arkiv Mat. Astr. Fysik 26A, 1 (1938).
- (24) P. W. Anderson, Mater. Res. Bull. 8, 153 (1973).
- (25) P. Fazekas and P. W. Anderson, Philos. Mag. 30, 23 (1974).
- (26) B. Sutherland, Phys. Rev. B 37, 3786 (1988).
- (27) S. Liang, B. Doucot, and P. W. Anderson, Phys. Rev. Lett. 61, 365 (1988).
- (28) I. Bose and P. Mitra, Phys. Rev. B 44, 443 (1991).
- (29) I. Bose, Phys. Rev. B 45, 13072 (1992).
- (30) H.-C. Jiang, H. Yao, and L. Balents, Phys. Rev. B 86, 024424 (2012).
- (31) W.-J. Hu, F. Becca, A. Parola, and S. Sorella, Phys. Rev. B 88, 060402(R) (2013).
- (32) M. P. Gelfand, R. R. P. Singh, and D. A. Huse, Phys. Rev. B 40, 10801 (1989).
- (33) M. P. Gelfand, Phys. Rev. B 42, 8206 (1990).
- (34) R. R. P. Singh and R. Narayanan, Phys. Rev. Lett. 65, 1072 (1990).
- (35) M. E. Zhitomirsky and K. Ueda, Phys. Rev. B 54, 9007 (1996).
- (36) P. W. Leung and N. W. Lam, Phys. Rev. B 53, 2213 (1996).
- (37) V. N. Kotov, J. Oitmaa, O. P. Sushkov, and Z. Weihong, Phys. Rev. B 60, 14613 (1999).
- (38) R. R. P. Singh, Z. Weihong, C. J. Hamer, and J. Oitmaa, Phys. Rev. B 60, 7278 (1999).
- (39) V. N. Kotov and O. P. Sushkov, Phys. Rev. B 61, 11820 (2000).
- (40) L. Capriotti and S. Sorella, Phys. Rev. Lett. 84, 3173 (2000).
- (41) M. S. L. du Croo de Jongh, J. M. J. van Leeuwen, and W. van Saarloos, Phys. Rev. B 62, 14844 (2000).
- (42) K. Takano, Y. Kito, Y. Ōno, and K. Sano, Phys. Rev. Lett. 91, 197202 (2003).
- (43) M. Mambrini, A. Läuchli, D. Poilblanc, and F. Mila, Phys. Rev. B 74, 144422 (2006).
- (44) V. Murg, F. Verstraete, and J. I. Cirac, Phys. Rev. B 79, 195119 (2009).
- (45) J. Reuther and P. Wölfle, Phys. Rev. B 81, 144410 (2010).
- (46) J. Reuther, P. Wölfle, R. Darradi, W. Brenig, M. Arlego, and J. Richter, Phys. Rev. B 83, 064416 (2011).
- (47) J.-F. Yu and Y.-J. Kao, Phys. Rev. B 85, 094407 (2012).
- (48) F. Figueirido, A. Karlhede, S. Kivelson, S. Sondhi, M. Rocek, and D. S. Rokhsar, Phys. Rev. B 41, 4619 (1990).
- (49) T. Oguchi and H. Kitatani, J. Phys. Soc. Jpn. 59, 3322 (1990).
- (50) P. Locher, Phys. Rev. B 41, 2537 (1990).
- (51) H. J. Schulz and T. A. L. Ziman, Europhys. Lett. 18, 355 (1992).
- (52) Q. F. Zhong and S. Sorella, Europhys. Lett. 21, 629 (1993).
- (53) J. Oitmaa and Z. Weihong, Phys. Rev. B 54, 3022 (1996).
- (54) L. Capriotti, F. Becca, A. Parola, and S. Sorella, Phys. Rev. Lett. 87, 097201 (2001).
- (55) G.-M. Zhang, H. Hu, and L. Yu, Phys. Rev. Lett. 91, 067201 (2003).
- (56) L. Capriotti, F. Becca, A. Parola, and S. Sorella, Phys. Rev. B 67, 212402 (2003).
- (57) L. Capriotti, D. J. Scalapino, and S. R. White, Phys. Rev. Lett. 93, 177004 (2004).
- (58) L. Capriotti and S. Sachdev, Phys. Rev. Lett. 93, 257206 (2004).
- (59) S. Yunoki and S. Sorella, Phys. Rev. Lett. 92, 157003 (2004).
- (60) L. Wang, Z.-C. Gu, X.-G. Wen, and F. Verstraete, arxiv:1112.3331 (2011).
- (61) F. Mezzacapo, Phys. Rev. B 86, 045115 (2012).
- (62) T. Li, F. Becca, W. Hu, and S. Sorella, Phys. Rev. B 86, 075111 (2012).
- (63) L. Wang, D. Poilblanc, Z.-C. Gu, X.-G. Wen, and F. Verstraete, Phys. Rev. Lett. 111, 037202 (2013).
- (64) A. F. Albuquerque and F. Alet, Phys. Rev. B 82, 180408 (2010).
- (65) Y. Tang, A. W. Sandvik, and C. L. Henley, Phys. Rev. B 84, 174427 (2011).
- (66) H. Ju, A. B. Kallin, P. Fendley, M. B. Hastings, and R. G. Melko, Phys. Rev. B 85, 165121 (2012).
- (67) K. S. D. Beach and A. W. Sandvik, Nucl. Phys. B 750, 142 (2006).
- (68) K. S. D. Beach, M. Mambrini, and F. Alet, Phys. Rev. B 77, 146401 (2008).
- (69) P. L. Iske and W. J. Caspers, Physica A 142 360 (1987).
- (70) Y.-C. Lin, Y. Tang, J. Lou, A. W. Sandvik, Phys. Rev. B 86, 144405 (2012).
- (71) The full set of amplitudes irreducible by symmetry is , where the Manhattan length is odd and the coordinates sit in the ordered range .
- (72) D. H. Lehmer, Proc. Sympos. Appl. Math. 10, 179 (1960).
- (73) D. E. Knuth, The Art of Computer Programming, Volume 3: Sorting and Searching (Addison-Wesley, Reading, MA, 1973), p. 12. ISBN 0-201-03803-X
- (74) P. Borwein and T. Erdélyi, “Horner’s Rule.” §1.1.E.5 in Polynomials and Polynomial Inequalities, (Springer-Verlag, New York, 1995). ISBN 978-0387945095
- (75) J. Lou and A. W. Sandvik, Phys. Rev. B 76, 104432 (2007).
- (76) O. P. Sushkov, J. Oitmaa, and Z. Weihong , Phys. Rev. B 63, 104420 (2001).