Topology and criticality in the resonating Affleck-Kennedy-Lieb-Tasaki loop spin liquid states
We exploit a natural Projected Entangled-Pair State (PEPS) representation for the resonating Affleck-Kennedy-Lieb-Tasaki loop (RAL) state. By taking advantage of PEPS-based analytical and numerical methods, we characterize the RAL states on various two-dimensional lattices. On square and honeycomb lattices, these states are critical since the dimer-dimer correlations decay as a power law. On kagome lattice, the RAL state has exponentially decaying correlation functions, supporting the scenario of a gapped spin liquid. We provide further evidence that the RAL state on the kagome lattice is a spin liquid, by identifying the four topological sectors and computing the topological entropy. Furthermore, we construct a one-parameter family of PEPS states interpolating between the RAL state and a short-range Resonating Valence Bond state and find a critical point, consistent with the fact that the two states belong to two different phases. We also perform a variational study of the spin-1 kagome Heisenberg model using this one-parameter PEPS.
The quest for physical systems harboring quantum spin liquid states has been a long-standing challenge in condensed matter physics Balents-2010 (). These states do not exhibit any symmetry breaking and thus admit no descriptions in terms of local order parameters. In contrast to conventional ordered states, spin liquids may have topological order Wen-1990 (); Wen-Niu-1990 () and support exotic fractional excitations. Among various theoretical constructions, a class of gapped spin liquids has been well-studied and recently receives strong numerical support White-2011 (); Schollwoeck-2012 (); Jiang-2012a (); Wang-2011 () as candidate ground states of some physical spin Hamiltonians on frustrated lattices. Their ground-state wave functions contain long-range entanglement and the low-energy effective theory is conveniently described by an emergent gauge theory Kitaev-2003 (); Senthil-2000 (), which also offers an intuitive picture for the highly-entangled ground state: a soup of fluctuating electric field lines Levin-2005 (). Moreover, if the spin system has symmetry, the gauge charges and fluxes in different spin liquid states may transform as linear (integer spin) or projective (half-integer) representations under the action of the symmetry group, i.e., they may belong to different classes of Symmetry Enriched Topological (SET) phases Essin-2013 (); Mesaros-2013 (); Hung-2013 (); Lu-2013 (); Fidkowski-unpub ().
In this paper, we characterize a type of spin- spin liquids Liu-2010 (); FWang-2011 (); Serbyn-2011 (); Cenke-2012 (); Bieri-2012 (); ZXLiu-2012 (); Tu-2013 (); ZCai-2007 (). The representative wave function of these states is a superposition of strongly fluctuating, fully packed loops (see Fig. 1), where each loop carries an Affleck-Kennedy-Lieb-Tasaki (AKLT) state Affleck-1987 (). We call this state a resonating AKLT-loop (RAL) state following Ref. [HYao-2010, ], where it was considered as an example of nontrivial SET phases. We provide an explicit Projected Entangled-Paired State (PEPS) representation of the state Verstraete-2004 (); Verstraete-2006 (), which allows us to exploit efficient PEPS-based analytical and numerical techniques to characterize these wave functions on various lattices. We demonstrate that on bipartite lattices (e.g. square and honeycomb lattices) these states have algebraically decaying dimer-dimer correlations and exponentially decaying spin-spin and quadrupole-quadrupole correlations, indicating a gapless spin liquid. On the non-bipartite kagome lattice, all correlation functions we have tested (spin, dimer and quadrupole) decay exponentially, whose correlation lengths do not increase with the system size. Additionally, we compute the topological entanglement entropy Kitaev-2006 (); Levin-2006 () of these wave functions on a long cylinder and obtain a universal value . All these results are indicative of a gapped spin liquid.
The AKLT state is a paradigmatic example of symmetry-protected topological phase in one dimension Gu-2009 (); Pollmann-2010 (); Schuch-2011 (), featuring degenerate spin- excitations on boundaries protected by the spin rotation symmetry. With no surprise, the RAL state also exhibits interesting symmetry properties. Namely, terminations of the loops always carry spin- degrees of freedom. In contrast, the spin-1 short-range RVB state, unlike its spin-1/2 counterpartAnderson-1973 (); Rokhsar-1988 (), does not support deconfined spin- excitations. We discuss the manifestation of this interesting symmetry enrichment in the construction of the degenerate topological sectors on an infinitely long cylinder. As a demonstration, we construct a one-parameter family of PEPS states interpolating between the RAL state and the RVB state. We find a critical point along the path, where the dimer-dimer correlation function decays algebraically. We also use these states as variational wave functions for the spin-1 kagome Heisenberg model and find an upper bound for the ground-state energy.
The rest of the paper is organized as follows. In Sec. II, we introduce the PEPS representations of the RAL states and related PEPS algorithms for extracting physical quantities. In Sec. III, we show that the RAL states are critical on square and honeycomb lattices. In Sec. IV, we exploit a combination of analytical and numerical techniques to provide evidence that the RAL state on a kagome lattice is a gapped spin liquid with topological order. In Sec. V, a one-parameter PEPS based on the RAL state are utilized for a variational study of the spin-1 Heisenberg model. Lastly, the Sec. VI is devoted to a summary.
Ii PEPS representation and algorithms
In this section, we introduce the PEPS representation of the RAL state and describe the PEPS algorithms utilized to extract physical quantities.
ii.1 PEPS representation for RAL states
Let us introduce the RAL state on a square lattice. In this case, we associate each lattice site with four virtual particles, each of which has three basis states, , , and , representing a spin-0 vacancy and a spin- doublet, respectively. In the PEPS language that we shall extensively use below, the virtual Hilbert space is denoted by the representation (i.e. direct sum of spin 0 and spin 1/2), whose dimension is called bond dimension, denoted by . Thus, we have in the present case. The RAL state can be constructed as
where is a maximally entangled bond state connecting the virtual particles between neighboring sites and , defined by , and is a local projector acting on site , which maps two virtual spin-1/2 states, out of the four virtual particles, onto the physical spin-1 Hilbert space [see Fig. 2(a) for an illustration]. More explicitly, can be written as
where, for instance, is given by
Here and denote the physical and virtual Hilbert spaces, respectively. is the Clebsch-Gordan (CG) coefficient symmetrizing two spin-1/2 particles into a physical spin-1, with nonvanishing coefficients being and .
To uniquely define the wave function, one needs to specify the orientation of the virtual singlets up to a gauge choice. For square and honeycomb lattices, we use the sign conventions shown in Fig. 2(a) and (b), where the singlets are oriented according to the arrows. For the kagome lattice, we consider two inequivalent sign conventions shown in Fig. 2(c) and (d).
Up to now, we have introduced the RAL state (1) in the PEPS language. If we substitute (2) into (1) and expand the product , the resonating AKLT-loop picture of (1) shown in Fig. 1 also becomes transparent, as each configuration contains different patterns of fully packed spin-1 AKLT loops covering the whole lattice. Since the projector in (2) gives the same weight to all allowed ways of combining two spin-1/2’s into a physical spin-1, the wave function (1) can be viewed as an equal weight superposition (up to a sign depending on the sign convention) of all possible AKLT-loop configurations, which was first proposed in Ref. [HYao-2010, ].
At this point, we note that our PEPS representation naturally generalizes to more complicated resonating loop states, where the loops carry other 1D matrix-product states (e.g. AKLT state Tu-2008 ()).
ii.2 RAL states on the kagome lattice
A straightforward construction of the RAL state on the kagome lattice proceeds along the description in the previous subsection. As shown in Fig. 3(a), there are four virtual particles on each vertex, in the representation , and the on-site projection operator is exactly the same as Eq. (3).
For the RAL state on the kagome lattice, this PEPS representation has a redundancy which can be revealed by the following procedure: we group the two virtual particles belonging to the same triangle (see e.g. in Fig. 3(a)) and block them into a single virtual particle (see in Fig. 3(b)). It turns out that not all virtual degrees of freedom survive after this blocking procedure. After removing this redundancy, we find that the new virtual particle is in the representation and has dimension (see Fig. 3(b)), instead of from a naive counting of the dimension of the tensor product . Let us note that this procedure is exact and the removal of redundancy allows us to reduce computational costs in our numerical calculations.
Furthermore, in our framework it is straightforward to introduce a one-parameter family of PEPS, which interpolates between the RAL and the spin-1 RVB states. This is done by extending the virtual representation [ in Fig. 3 (a)] from to . Thus, the virtual bonds in (1) are modified as , where denotes the three states in the virtual spin-1 space. Accordingly, the projector in (2) has to be modified as
where is a local projector which maps one virtual spin-1 state, out of the four virtual spin-1 particles, onto the physical spin-1 Hilbert space. Its explicit form is given by
where, for instance, is defined as
Here store the trivial CG coefficients . In this one-parameter PEPS, we recover the RAL state (1) for and the RVB state for . Hence we have a continuous family of “mixed” RAL states parametrized by the weight controlling the density of the spin-1 valence bonds in the mixed loop-dimer configuration. Like in the pure RAL state, we can combine two virtual particles into a single virtual particle (see Fig. 3). As a result, the virtual particle is also in the representation and thus the resulting PEPS has bond dimension .
ii.3 PEPS algorithms
In our numerical simulations, we mainly deal with the RAL states defined on (i) a cylinder with finite circumference and infinite length and (ii) infinite-size 2D lattices. For the former, when the circumference is relatively small ( for RAL state), we perform exact contractions for extracting physical observables (e.g. correlation functions) and entanglement properties. In the latter case, we employ the infinite-PEPS (iPEPS) algorithm Jordan-2008 (); Orus-2009 () for approximate contractions with high precision.
Firstly, let us briefly explain the exact contraction method on a cylinder. In order to evaluate the norm of the wave function, or the expectation values of concerned physical quantities, we need to calculate the inner-product of two PEPS wave functions by contracting the double-layer tensor network from both open ends along the (infinitely long) horizontal direction. In each step, the boundary vector is contracted with a single column of transfer operators, and the size of the boundary vectors does not grow up after a step of contraction. Therefore, this process can be repeated until the vectors, as well as the concerned observables, converge. In practice, it takes iteration steps to converge. However, when treating topological PEPS, we note that the initial boundary vector has to be chosen with special care, since it plays an important role in selecting the topological sectors, as we shall show in later sections.
When both directions of the lattice are infinite, exact contraction is no longer feasible. In this case, we utilize the iPEPS technique Jordan-2008 (); Orus-2009 () to contract the tensor network with high accuracy. The basic idea is to (approximately) express the boundary vector in the form of a translation-invariant boundary matrix-product state (MPS), and to contract the boundary MPS with the transfer operator (in a matrix-product operator form, see Fig. 4) of the double-layer tensor network. The size of the boundary MPS, i.e., the dimension of its virtual bonds (see Fig. 4), increases exponentially with contraction steps. Therefore, in order to make the contraction process sustainable, it is necessary to truncate the virtual bond space of the boundary MPS after a few steps. This introduces the truncation bond dimension, denoted by . We adopt both the standard canonicalizationVidal-2008 () and bi-canonicalizationXiang-2013 () truncation techniques, and compare the results with the extrapolated value obtained in exact contractions to double check the reliability and accuracy of our numerical calculations.
Iii Criticality on bipartite lattices
In this section, we consider the RAL states on two bipartite lattices, i.e., square and honeycomb lattices. For the honeycomb RAL state [see Fig. 2 (b) for its PEPS form], we block two nearest-neighbor sites together and thus obtain effectively a square-lattice PEPS. When computing two-point correlation functions, the distances are measured along a line on the effective square lattice, which corresponds to a zig-zag line on the original honeycomb lattice.
We employ the iPEPS method to accurately evaluate the correlation functions. When doing this, we have checked the convergence of the results with different truncation parameter . The calculated correlation functions are shown in Fig. 5.
The spin-spin correlation function and the quadrupole-quadrupole correlation function , where , both decay exponentially as can be seen in Fig. 5(a,c), suggesting the absence of Néel or quadrupolar long-range order. However, the dimer-dimer correlation function, , decays algebraically on both lattices, as shown in Fig. 5(b,d). Therefore, we conclude that the RAL states are critical on square and honeycomb lattices HYao-2010 (). We expect that it is a consequence of the fully-packed arrangement of the loops, similar to those critical loop models in classicial statistical mechanics Nienhuis-1982 (); Nienhuis-1994 ().
Iv The RAL states on a kagome lattice
iv.1 topological order
Now we turn to the non-bipartite kagome lattice and provide numerical evidence that the kagome RAL state is gapped and possesses topological order. In this Section, we only present results of the RAL state with sign convention shown in Fig. 2(c). For the other sign convention as in Fig. 2(d), the results lead to the same conclusion and thus are not present here.
We first compute the spin-spin, dimer-dimer and quadrupole-quadrupole correlation functions. These results are shown in Fig. 6(a). We note that, by adopting the simplex construction in Fig. 3(b), one obtains an effective honeycomb-lattice representation of the RAL states, which can then be transformed to a square-lattice PEPS as we did in the previous Section. The correlation functions are thus measured along a line of this coarse-grained square lattice. The correlation functions in Fig. 6(a) all exhibit exponentially decaying behaviors, supporting that the kagome RAL state is a gapped spin liquid.
To characterize the RAL state on the kagome lattice, we compute the entanglement entropy on an infinite cylinder following the scheme described in Ref. [Jiang-2012, ], utilizing the PEPS technique in Ref. [Cirac-2011, ]. The cylinder is assumed to be periodic (open) along the vertical (horizontal) direction. We trace over the states on half of the cylinder to obtain the reduced density matrix and the entanglement entropy. Before presenting the numerical results, we first analyze the topological sectors on the cylinder to construct the so-called minimally entangled states YZhang-2012 (); YZhang-2011 (); Poilblanc-2012 (); Schuch-2012 (); Schuch-2013 (); Poilblanc-2013 (). Since the wave function closely resembles the string picture of a topological phase (or a gauge theory), we use the language of the gauge theory throughout our discussion.
We can easily identify four topological sectors distinguished by two non-local operators, the magnetic flux along the horizontal (vertical) direction denoted by , which counts the parity of the number of AKLT valence bonds along a horizontal (vertical) cut. First we set both ends of the cylinder open (i.e. AKLT lines can not terminate on the ends). The two sectors can be chosen as eigenstates of with eigenvalues , as shown in Fig. 7(a,b). Shifting the virtual valence bonds around the circumference of the cylinder toggles between the two sectors. This operator is exactly the Wilson loop denoted by , i.e., create a pair of spinons by breaking a virtual valence bond, wind one spinon out of the two around the cylinder and annihilate them again on the other side). In practice, we insert a line of diagonal -matrices to the PEPS (i.e., thread a magnetic flux along the cylinder, see Fig. 8 for illustration), and then construct two eigenstates of : which are minimally entangled YZhang-2012 () and ready for entanglement entropy calculation.
In order to construct the other two sectors, we are enforced to introduce virtual polarized spin- spinons [between which there is an open string of AKLT state, see Fig. 7(c,d)] on the two boundaries of the cylinder. They act as two charges where the electric field lines can terminate. It is easy to see that a vertical path necessarily cuts through odd number of virtual spin-1/2 valence bonds with this particular boundary conditions comment:sym-breaking (). We refer to the sectors with this boundary condition as the odd topological sector (), while the previous two with open boundaries as the even ones (). In the even sector, the reduced density matrix only has nonzero eigenvalues corresponding to eigenvectors with integer total spin; while in the odd sector, because of the virtual spin-1/2 spinons on the boundaries, all Schmidt eigenstates have half-odd-integer total spin, i.e., each level in the entanglement spectrum has an even-fold degeneracy (a detailed analysis on entanglement spectra will be present latter in section IV.4), which we have confirmed numerically comment:sym-breaking (). This is a manifestation of the fact that the charge forms a representation of the symmetry group while the flux excitation carries no nontrivial projective representation (i.e. integer spin).
In Fig. 6(b), we present the entanglement entropies of the two minimally entangled states in even and odd sectors. The entanglement entropy as a function of the circumference can be evaluated using both converged (left and right) boundary vectors obtained from the exact contraction Cirac-2011 (). The topological entropy is extracted by fitting to and we find (close to ) in both sectors, indicating that the kagome RAL state has a topological order.
Fig. 6(c) shows the variational energies of the RAL states in all four topological sectors with respect to the spin-1 antiferromagnetic Heisenberg model on the kagome lattice. When increasing , the variational energies of the wave functions in four different topological sectors (, ) converge rapidly, all getting close to the value obtained by iPEPS method, suggesting that they are degenerate in the thermodynamic limit.
iv.2 Characterization of the spin-1 RVB state
So far we have constrained ourselves to the pure RAL state ( in Eq. 4), and excluded the spin-1 valence bonds (equivalent to a short AKLT loop with length 2). However, it is interesting to consider the other limit () where only the spin-1 short-range valence bonds are allowed.
In the case of , we have a spin-1 short-range RVB state, which can be represented by a PEPS with bond dimension , i.e, the virtual spin is in the representation . The vision line can be defined as , and the odd topological sector is constructed by attaching a pair of virtual spin-1 on both ends.
In Fig. 9, we present the correlation functions (including the variational energies of the kagome Heisenberg model from nearest-neighbor spin-spin correlations) and the entanglement entropy of the spin-1 RVB state on a kagome lattice. Fig. 9(a) shows that both the spin-spin and dimer-dimer correlation functions decay exponentially, suggesting that the spin-1 kagome RVB state is fully gapped. In Fig. 9(b), we calculate the variational energies of the spin-1 kagome Heisenberg model in four topological sectors. All of them converge to the iPEPS energy result rapidly when the circumference increases. Fig. 9(c) shows the von Neumann entanglement entropies as a function of the circumference , from which we extract the topological entanglement entropy (close to ). This is a clear indication of the underlying topological order.
iv.3 Interpolating between the RAL and RVB states
In the previous subsections, we have shown that both the pure RAL () and the spin-1 RVB () states on a kagome lattice have topological order. As mentioned in Sec. II.2, the PEPS defined by Eq. 4 can describe the states with mixed loop and dimer configurations on equal footing, and we can tune the parameter to interpolate between the pure RAL and spin-1 RVB states.
The following question arises naturally: can this path adiabatically connect the spin-1 RAL and RVB states? Theoretically, although the two limiting states have identical intrinsic topological order, they behave very differently under spin symmetry. To be more concrete, the RVB state has no deconfined half-integer spin excitations. Since the path under consideration preserves the symmetry all the way from to , we expect there must be at least one singular point where the state becomes critical. In Fig. 10 we show that, although the spin-spin correlators always decay exponentially [see Fig. 10(a)], the dimer-dimer correlation function decays algebraically at around , as depicted in Fig. 10(b).
For this one-parameter family of PEPS , one can always find a local, positive semidefinite Hamiltonian such that , i.e., the RAL state is the exact ground state of this co-called parent Hamiltonian Schuch-2010 (). The existence of a critical point implies that the path in the space of parent Hamiltonians should be gapless at . However, as the dimension of virtual Hilbert space for this PEPS is , one needs to block a large number of sites so that a local projector annihilating the PEPS exists, which constitutes the parent Hamiltonian. Due to its complicated form, we do not present the analytical form of the parent Hamiltonian here.
iv.4 Entanglement spectra
Entanglement spectrum (ES) is defined as the eigenvalues of the operator , where is the reduced density matrix. The full ES has been suggested to reveal entanglement properties, providing more useful information than the entanglement entropy Haldane-2008 ().
In this section, we provide the entanglement spectra of the pure RAL, spin-1 RVB and the mixed RAL states in different topological sectors. The construction of four topological sectors for the pure RAL and spin-1 RVB states have been discussed in Secs. IV.1 and IV.2.
First we consider the ES of the pure RAL () and spin-1 RVB () states, shown in Fig. 11. We have taken advantange of the spin-rotational symmetry and the translation symmetry in the vertical (circumference) direction of the cylinder and label each level in the spectrum by its spin (a spin- multiplet represents states) and momentum (in vertical direction). For both states, the ES in the even sectors are filled with integer-spin multiplets and do not have any recognizable features regarding the degeneracy. However, the RAL and RVB states have distinct ES in the odd sectors. The ES levels of the RAL state all have half-integer spins, giving rise to even-fold degeneracies in the spectrum. For the spin-1 RVB state, the ES in the odd sector are still filled with integer-spin multiplets.
One interesting question regarding the entanglement spectrum shown in Fig. 11 is whether the “low energy” part of the ES is critical and described by a conformal field theory. However, as the largest system size we can achieve for computing the ES is , which is still too small for a finite-size scaling of the gap, we do not have a definite answer to this question and leave it for a future study.
For the mixed RAL states with , it is not entirely clear how to construct the topological sectors. Due to the phase transition at , we expect that the ES for in the odd sector is qualitatively different from those for . Since the mixed RAL state for is adiabatically connected to the pure RAL state, we expect that the four topological sectors also evolve adiabatically in this parameter region. Formally we can construct a series of wave functions, for both and , with the same boundary conditions (i.e. put a pair of oppositely polarized spinons on two open ends of the cylinder comment:sym-breaking ()). We calculate the ES of these states and the results are shown in Fig. 12. All levels in the ES have even-fold degeneracies for (RAL region), while they are featureless regarding the degeneracy for (RVB region).
V Variational study of the spin-1 kagome Heisenberg model
In this section, we use the RAL states as variational wave functions for the spin-1 Heisenberg model
on the kagome lattice, whose ground-state properties remain elusive to the best of our knowledge.
We start with the pure RAL state as a parameter-free variational wave function for spin-1 kagome Heisenberg model. We first compare with the exact diagonalization (ED) result on a 18-site torus with unit cells (see Fig. 13), which is the largest size considered in the ED studyHida (). We place the PEPS tensors of RAL states on this small torus geometry, and evaluate the energy expectation value by exact contraction. The variational energy of the pure RAL state with the sign convention given in Fig. 2(c) is per site, closer to the ED value than that of the hexagonal singlet solid variational state ( per site) Hida (); Hida_HSS ().
Next, we extend the variational study to the one-parameter family of PEPS wave functions, i.e., the mixed RAL states, which interpolate the pure RAL and spin-1 RVB states (controlled by the parameter ). The results from the iPEPS calculations in the thermodynamic limit is shown in Fig. 14. The best variational energy (per site) is achieved at with the sign convention given in Fig. 2(d).
To conclude, we have systematically investigated a family of resonating AKLT-loop states on square, honeycomb, and kagome lattices. Using a natural PEPS representation, we have shown that the RAL states are critical on square and honeycomb lattices, while on kagome lattice it is a gapped spin liquid. We also discussed the realization of the spin-rotation symmetry and clarified its manifestation through explicitly constructing the topological sectors and evaluating the corresponding entanglement spectra on infinitely long cylinders. We considered a one-parameter family of PEPS interpolating between the RAL and RVB states which have distinct symmetry realizations. A critical point has been identified along this interpolation path. Lastly, we have used the RAL states to obtain the best-to-date variational energy for the spin- Heisenberg model on a kagome lattice.
We thank Ignacio Cirac, Didier Poilblanc, Yi-Zhuang You, Zi Cai, Hong Yao and Fang-Zhou Liu for helpful discussions. WL acknowledges the hospitality of the Max-Planck Institute for Quantum Optics, where part of the work has been performed. HHT gratefully acknowledges Tai-Kai Ng and K. T. Law for hospitality during his visit in Hong Kong University of Science and Technology. This work has been supported by the DFG through SFB-TR12, the EU project SIQS, the DFG Cluster of Excellence NIM, and NSFC 11204149.
Note added.– Upon finalizing the manuscript we noticed a recent preprint Huang-2013 () on closely related topics.
- (1) L. Balents, Nature 464, 199 (2010).
- (2) X.-G. Wen, Int. J. Mod. Phys. B 4, 239 (1990).
- (3) X.-G. Wen and Q. Niu, Phys. Rev. B 41, 9377 (1990).
- (4) S. Yan, D. A. Huse, and S. White, Science 332, 1173 (2011).
- (5) S. Depenbrock, I. P. McCulloch, and U. Schollwöck, Phys. Rev. Lett. 109, 067201 (2012).
- (6) H.-C. Jiang, H. Yao and L. Balents, Phys. Rev. B 86, 024424 (2012).
- (7) L. Wang, Z.-C. Gu, F. Verstraete and X.-G. Wen, arXiv:1112.3331.
- (8) A. Kitaev, Ann. Phys. 303, 2 (2003).
- (9) T. Senthil and M. Fisher, Phys. Rev. B 62, 7850 (2000).
- (10) M. A. Levin and X.-G. Wen, Phys. Rev. B 71, 045110 (2005).
- (11) A. M. Essin and M. Hermele, Phys. Rev. B 87, 104406 (2013).
- (12) A. Mesaros and Y. Ran, Phys. Rev. B 87, 155115 (2013).
- (13) L.-Y. Hung and X.-G. Wen, Phys. Rev. B 87, 165107 (2013).
- (14) Y.-M. Lu and A. Vishwanath, arXiv:1302.2634.
- (15) L. Fidkowski, N. Lindner and A. Kitaev, unpublished.
- (16) Z.-X. Liu, Y. Zhou, and T.-K. Ng, Phys. Rev. B 81, 224417 (2010); 82, 144422 (2010).
- (17) M. Serbyn, T. Senthil, and P. A. Lee, Phys. Rev. B 84, 180403(R) (2011).
- (18) F. Wang and C. Xu, arXiv:1110.4091.
- (19) C. Xu, F. Wang, Y. Qi, L. Balents, and M. P. A. Fisher, Phys. Rev. Lett. 108, 087204 (2012).
- (20) S. Bieri, M. Serbyn, T. Senthil, and P. A. Lee, Phys. Rev. B 86, 224409 (2012).
- (21) Z.-X. Liu, Y. Zhou, H.-H. Tu, X.-G. Wen, and T.-K. Ng, Phys. Rev. B 85, 195144 (2012).
- (22) H.-H. Tu, Phys. Rev. B 87, 041103(R) (2013).
- (23) Z. Cai, S. Chen, S. P. Kou, and Y. P. Wang, Phys. Rev. B 76, 054443 (2007); Z. Cai, S. Chen, and Y. P. Wang, J. Phys.: Condens. Matter, 21, 456009 (2009).
- (24) I. Affleck, T. Kennedy, E. H. Lieb, and H. Tasaki, Phys. Rev. Lett. 59, 799 (1987); Commun. Math. Phys. 115, 477, (1988).
- (25) H. Yao, L. Fu, and X.-L. Qi, arXiv:1012.4470.
- (26) F. Verstraete and J. I. Cirac, cond-mat/0407066; F. Verstraete and J. I. Cirac, Phys. Rev. A 70, 060302 (2004).
- (27) F. Verstraete, M. M. Wolf, D. Perez-Garcia, and J. I. Cirac, Phys. Rev. Lett. 96, 220601 (2006).
- (28) A. Kitaev and J. Preskill, Phys. Rev. Lett. 96, 110404 (2006).
- (29) M. Levin and X.-G. Wen, Phys. Rev. Lett. 96, 110405 (2006).
- (30) Z.-C. Gu and X.-G. Wen, Phys. Rev. B 80, 155131 (2009).
- (31) F. Pollmann, A. M. Turner, E. Berg, and M. Oshikawa, Phys. Rev. B 81, 064439 (2010).
- (32) N. Schuch, D. Perez-Garcia, and J. I. Cirac, Phys. Rev. B 84, 165139 (2011).
- (33) P. W. Anderson, Mater. Res. Bull. 8, 153 (1973); Science 235,1196 (1987).
- (34) D. S. Rokhsar and S. A. Kivelson, Phys. Rev. Lett. 61, 2376 (1988).
- (35) H.-H. Tu, G.-M. Zhang, and T. Xiang, Phys. Rev. B 78, 094404 (2008); J. Phys. A 41, 415201 (2008); H.-H. Tu, G.-M. Zhang, T. Xiang, Z.-X. Liu, and T.-K. Ng, Phys. Rev. B 80, 014401 (2009).
- (36) J. Jordan, R. Orus, G. Vidal, F. Verstraete, and J. I. Cirac, Phys. Rev. Lett. 101, 250602 (2008).
- (37) R. Orus and G. Vidal, Phys. Rev. B 80 094403 (2009); R. Orus, Phys. Rev. B 85, 205117 (2012).
- (38) R. Orus and G. Vidal, Phys. Rev. B 78, 155117 (2008).
- (39) Because the matrix-product transfer operator of the bilayer tensor network is non-Hermitian, the left- and right-boundary MPS are different, which can be simultaneously canonicalized and truncated with the bi-canonicalization technique. J. Chen and T. Xiang, private communication.
- (40) B. Nienhuis, Phys. Rev. Lett. 49, 1062 (1982).
- (41) H. W. J. Blöte and B. Nienhuis, Phys. Rev. Lett. 72, 1372 (1994).
- (42) H.-C. Jiang, Z. Wang, and L. Balents, Nature Phys. 8, 902 (2012).
- (43) J. I. Cirac, D. Poilblanc, N. Schuch, and F. Verstraete, Phys. Rev. B 83, 245134 (2011).
- (44) Y. Zhang, T. Grover, A. Turner, M. Oshikawa, and A. Vishwanath, Phys. Rev. B 85, 235151 (2012).
- (45) Y. Zhang, T. Grover, and A. Vishwanath, Phys. Rev. B 84, 075128 (2011).
- (46) D. Poilblanc, N. Schuch, D. Perez-Garcia, and J. I. Cirac, Phys. Rev. B 86, 014404 (2012).
- (47) N. Schuch, D. Poilblanc, J. I. Cirac, and D. Perez-Garcia, Phys. Rev. B 86, 115108 (2012).
- (48) N. Schuch, D. Poilblanc, J. I. Cirac, and D. Perez-Garcia, Phys. Rev. Lett. 111, 090501 (2013).
- (49) D. Poilblanc and N. Schuch, Phys. Rev. B 87, 140407(R) (2013).
- (50) One may wonder that by putting polarized spinons on the boundaries we weakly break the symmetry through boundary conditions. However, since the boundaries are infinitely far away, this effect on the bulk density matrix is negligibly small. In addition, in case the half-integer spin sector is not a degenerate topological sector, the spinon-pair boundary does not effectively change the parity of entanglement spectrum in the bulk [for example, see ES in the spin-1 RVB region () shown in Fig. 12].
- (51) N. Schuch, J. I. Cirac, and D. Perez-Garcia, Ann. Phys. (Amsterdam) 325, 2153 (2010).
- (52) H. Li and F. D. M. Haldane, Phys. Rev. Lett. 101, 010504 (2008).
- (53) K. Hida, J. Phys. Soc. Jpn. 69, 4003 (2000).
- (54) K. Hida, J. Phys. Soc. Jpn. 70, 3673, (2001); J. Phys. Soc. Jpn. 71, 3021 (2002).
- (55) C.-Y. Huang, X. Chen, and F. Pollmann, arXiv:1312.3093.