Efficient MPS methods for extracting spectral information on rings and cylinders
Maarten Van Damme, Robijn Vanhove, Jutho Haegeman, Frank Verstraete, Laurens Vanderstraeten
EEfficient MPS methods for extracting spectral information on rings and cylinders
Maarten Van Damme, Robijn Vanhove, Jutho Haegeman, Frank Verstraete, and Laurens Vanderstraeten Department of Physics and Astronomy, University of Ghent, Krijgslaan 281, 9000 Gent, Belgium
Based on the MPS formalism, we introduce an ansatz for capturing excited states in finite systemswith open boundary conditions, providing a very efficient method for computing, e.g., the spectralgap of quantum spin chains. This method can be straightforwardly implemented on top of anexisting DMRG or MPS ground-state code. Although this approach is built on open-boundaryMPS, we also apply it to systems with periodic boundary conditions. Despite the explicit breakingof translation symmetry by the MPS representation, we show that momentum emerges as a goodquantum number, and can be exploited for labeling excitations on top of MPS ground states. Weapply our method to the critical Ising chain on a ring and the classical Potts model on a cylinder.Finally, we apply the same idea to compute excitation spectra for 2-D quantum systems on infinitecylinders. Again, despite the explicit breaking of translation symmetry in the periodic direction,we recover momentum as a good quantum number for labeling excitations. We apply this methodto the 2-D transverse-field Ising model and the half-filled Hubbard model; for the latter, we obtainaccurate results for, e.g., the hole dispersion for cylinder circumferences up to eight sites.
I. INTRODUCTION
Matrix product states (MPS) [1, 2], or the density-matrixrenormalization group (DMRG) [3, 4], provide an ef-ficient formalism for simulating one-dimensional (1-D)and quasi 1-D quantum lattice systems with very highprecision. Although MPS are introduced most elegantlyon a periodic system [5–7], MPS and DMRG algorithmsare traditionally formulated on finite systems with openboundary conditions, because this setting allows for anefficient calculation of expectation values and the fix-ing of the gauge degrees of freedom in canonical forms.The disadvantage of open boundaries is that the trans-lation symmetry of the model is explicitly broken, butthe formulation of MPS directly in the thermodynamiclimit [8–12] has made it possible to restore translationalsymmetries without sacrificing the efficiency or canonicalforms of MPS. In this setting, translation symmetry ofthe MPS ground state has been used as a basis for formu-lating the MPS version [13] of the Feynman-Bijl ansatz[14] or the single-mode approximation [15–18], which de-scribes excited-state wavefunctions with definite momen-tum quantum number with high precision for generic 1-Dlattice models [19–21].Yet, systems with periodic boundary conditions areimportant in, at least, two contexts. In both cases, peri-odic boundary conditions are used because the finite sizeeffects are much smaller for systems with periodic as op-posed to open boundary conditions. First, in the studyof 1-D critical models the low-energy spectrum of a finiteperiodic system provides a very clear fingerprint of theconformal field theory (CFT) that captures the infraredproperties of the model [22]. Therefore, computing thespectrum on a periodic system is paramount for identify-ing the effective CFT for a given critical model. Second,the finite size effects imply that the most efficient way ofstudying two-dimensional (2-D) quantum systems withMPS consists of imposing periodic boundary conditionsin one direction, i.e. studying the model on a cylindricalgeometry around which the MPS snakes or spirals. In the case of 1-D periodic systems, there have beena number of proposals for alleviating the computationalcost of periodic MPS [23–27] and DMRG [28, 29] sim-ulations. In addition, a variational ansatz for captur-ing elementary excitations on top of a periodic MPS wasproposed [30] and refined [27, 31], which captures manylow-lying energy levels that reflect the CFT finite-sizespectrum of critical spin chains[27]. Here, it is crucialthat translation invariance is explicitly conserved in theground state – by choosing a uniform MPS ansatz – suchthat the momentum is a good quantum number for la-beling the excited states.When systems with a cylindrical geometry are studiedusing a real space MPS ansatz [32], translation symme-try around the cylinder is explicitly broken by the MPSstructure. While one can hope that the symmetry isrestored in the ground state for sufficiently large bonddimension, transversal momentum cannot explicitly beused as a good quantum number in the MPS ansatz totarget excitations. In the case of fermionic systems, onecan switch to the momentum basis in the periodic di-rection such that translation symmetry can be preservedexplicitly [33, 34]. In the case of spin systems, the trans-formation to the transversal momentum basis cannot beimplemented as a canonical transformation and has assuch not been considered, as it would not exhibit a ten-sor product structure.In this paper, we show that we can recover momen-tum as a reliable quantum number to label excitationseven when the associated translation symmetry is ex-plicitly broken by the MPS structure. In particular, wewill use open-boundary MPS and consider its tangentspace as a variational ansatz for excitations, analogouslyto the excitation ansatz that has proven successful forinfinite systems. First, we set the stage by applying thisansatz to an actual finite 1-D spin system, namely tostudy the magnon excitations in the spin-1 Heisenbergmodel, and compare the results to traditional DMRG-based approaches for targeting excited states. Next, weapply the same formalism to 1-D systems with periodic a r X i v : . [ c ond - m a t . s t r- e l ] F e b boundary conditions, i.e. living on a ring, where we findthat momentum is recovered with good accuracy. Thisis useful in particular for studying critical systems, andwe obtain very accurate CFT spectra for Ising and Pottsmodels. Finally, we translate the same idea to capturingexcitations in infinite cylinders. Again, despite the ex-plicit breaking of translation symmetry around the cylin-der, the transversal momentum can be recovered accu-rately and we discuss a strategy to target specific mo-mentum sectors. We illustrate this approach for the 2-Dtransverse-field Ising model and the half-filled 2-D Hub-bard model. All algorithms in this paper can be found in MPSKit.jl [35], an open-source software package writtenin the scientific programming language Julia.
II. QUASIPARTICLES ON A FINITE SYSTEM
There are a variety of techniques to variationally targetexcited states in traditional DMRG or MPS simulations.The easiest is when the excitation has a non-trivial quan-tum number, such that it could directly be obtained asthe solution of a ground state problem in that symme-try sector [4]. If multiple states in the same sector arerequired, they can be simultaneously captured by an ex-tended DMRG procedure [4, 36], but the bond dimensiongrows rapidly with the number of targeted states. In theMPS representation, one would instead sequentially findhigher excited states by imposing orthogonality on pre-viously found states [37]. Finally, in the case of criticalsystems, it was realized recently that just changing onetensor in the middle of the chain is often sufficient tocapture the low-lying excitations [38]. In this section, weintroduce a variational ansatz for directly capturing theexcited states of finite quantum chains and show thatwe can outperform these standard approaches, both inefficiency and accuracy.
A. Method
Let us consider a finite one-dimensional spin chain oflength N with open boundary conditions, described bya generic model Hamiltonian H , for which the groundstate is described by an MPS | Ψ( A . . . A N ) (cid:105) = A A . . . A N . (1)This MPS representation has a residual gauge freedom.In particular, we can always choose the tensors left andright of a given site to be respectively left and right isome-tries, leading to the representation | Ψ( A . . . A N ) (cid:105) = A l . . . A ci . . . A rN A l . (2)with A li and A ri the left and right isometric MPS tensors,and A ci the center-site tensor. The standard sweepingalgorithm [1] can be used to find an optimal MPS ap-proximation for the ground state. Once we have found a good MPS approximation of theground state, we can build excitations on top. Inspiredby the success of the quasiparticle ansatz in the thermo-dynamic limit, we propose a very similar ansatz | Φ( B , . . . , B N ) (cid:105) = (cid:88) i A l . . . A rN B i . . . . (3)In this ansatz, we take a sum of N terms, where in eachterm we change one tensor. Note that one can representthis state as a finite MPS with two times the bond di-mension of the ground state, but there are only O ( N χ )free variational parameters.The quasiparticle ansatz exhibits a gauge freedom, asevery tensor B i can be modified as B (cid:48) i = B i + X i A ri − A li X i +1 and we would still describe the same state. Wecan follow the ideas from infinite MPS [12], and fix thisgauge freedom by imposing that B i lives in the null-spaceof ( A li ) † . The tensor B i is then parametrized as B i = V i X i (4)where X i contains the actual degrees of freedom and V i spans the null space of ( A li ) † , i.e. it contains the orthonor-mal columns to complement A li to a full unitary matrix.As a consequence, the quasiparticle state is automati-cally orthogonal to the ground-state MPS. Our ansatznow becomes | Φ( X . . . X N (cid:105) = (cid:88) i | Φ i ( X i ) (cid:105) , | Φ i ( X i ) (cid:105) = A l . . . V i X i . . . A rN . (5)The overlap of two distinct quasiparticle states simplifiesto the simple euclidean inner product of the tensors X i (cid:104) Φ( X . . . X N ) | Φ( X (cid:48) . . . X (cid:48) N ) (cid:105) = (cid:88) i ( (cid:126)X i ) † X (cid:48) i = ( (cid:126)X (cid:48) ) † (cid:126)X (6)This is important, as minimizing the energy within themanifold of quasiparticle excitations now becomes a sim-ple eigenvalue problem (cid:88) j ( H eff ) ij (cid:126)X j = ω (cid:126)X i , (7)with the effective hamiltonian matrix( (cid:126)X i ) † ( H eff ) ij (cid:126)X j = (cid:104) Φ i ( X i ) | H | Φ j ( X j ) (cid:105) . (8)The action of H eff on a set of tensor X i can be computedeasily, and the eigenvalue problem can be solved by aniterative Krylov method, typically the Lanczos method. Here we denote (cid:126)X i for the vectorized version of the tensor X , and (cid:126)X for the concatenation of the vectors (cid:126)X i into a large vector. FIG. 1. The energy density of the four lowest-lying spin-1magnon excitations on top of the ground state of the spin-1Heisenberg chain with 100 sites. The ground state has bonddimension D = 64. B. Magnon in the spin-1 chain
In infinite systems, the excitation ansatz leads to apicture of gapped excitations as dressed quasiparticlesagainst a correlated background state [39]; in that set-ting the ansatz describes a traveling wave or Bloch wavewith a definite momentum. On a finite system with openboundary conditions, however, we would expect to finda standing-wave configuration of this very same dressedparticle.As an illustration of this scenario, we simulate the exci-tation spectrum of the spin-1 Heisenberg chain on a finitesystem with open boundary conditions. In order to focuson the bulk excitations, we place spin-1/2s at the ends,thus eliminating the gapless edge modes[40]. In Fig. 1we plot the energy density of the first five excited states,showing indeed the different standing-wave patterns ofthe magnon behaving as a particle-in-a-box [41].We can compare this method for extracting excitedstates with the conventional MPS methods. As explainedabove, the standard method requires an entire sweepingoptimization for every excited state, and one needs higherbond dimensions to faithfully capture the excitation. Inthat respect, the quasiparticle ansatz is numerically muchcheaper as it only requires to solve a single eigenvalueproblem for a given number of excited states. Moreover,the excitation ansatz reaches the same level of accuracyfor the energy, as we show explicitly in Table I. One ex-pects the excitation ansatz to fail for higher excitations,as soon as they start to involve multiple particles. It thenmakes sense to go to a hybrid setup where we switch toDMRG after the first few lowest eigenvectors. This willnot be necessary however when doing simulations of asystem close to criticality and for which the correlationlength of the MPS is larger than the system size (i.e. thefinite size scaling regime): in that case, the local tensorsin the quasiparticle ansatz have a global effect, and seemto be able to represent multiparticle excitations [38].
Energy QP Energy DMRG Variance QP Variance DMRG0.4165739 0.4165864 3.1e-6 7.05e-50.4165739 0.4165865 3.1e-6 7.09e-50.4165739 0.4165864 3.1e-6 7.05e-50.4344129 0.4344322 3.2e-6 0.000110.4344129 0.4344321 3.2e-6 0.0001090.4344129 0.4344321 3.2e-6 0.0001090.462799 0.4628202 3.3e-6 0.00012070.462799 0.4628204 3.3e-6 0.00012090.462799 0.4628205 3.3e-6 0.00012090.5001103 0.5001329 3.3e-6 0.000129TABLE I. Comparison of the energies and energy vari-ances obtained using the quasiparticle (QP) ansatz and usingDMRG, for the lowest-lying excitations in the spin-1 chainwith N sites (also see Fig. 1). For the DMRG excitations, themaximum number of sweeps was set at 100. III. CRITICAL SYSTEMS ON A RING
As highlighted above, MPS techniques are significantlyless efficient for systems with periodic boundary condi-tions, in part due do the inferior scaling in bond dimen-sion when using an MPS with periodic boundary condi-tions. Still, both ground states and excited states can betargeted with high precision using periodic MPS.An alternative approach is to re-use the above tech-nique for open boundary systems, but with a periodichamiltonian or transfer matrix. Although it would re-quire a quadratically larger bond dimension to representa periodic MPS as an open-boundary MPS, we will showthat we obtain quantitatively good results at reasonablebond dimensions. In particular, we can compute the ex-pectation value of the translation operator, which for theMPS ground state [Eq. 1] is given by (cid:104) Ψ( A . . . ) | T | Ψ( A . . . ) (cid:105) = ¯ A ¯ A ¯ A N A A A N . . .. . . , (9)and show that the momentum of the ground states, aswell as the excited states , are correct up to very highprecision.In Ref. 27 it was shown that the excitation ansatz forperiodic MPS is able to reproduce a surprising number ofenergy levels in the finite-size CFT spectrum of critical Eigenvalues occur with different degeneracies and the corre-sponding eigenvectors will be momentum superpositions. Wethen have to diagonalize the translation operator within thisdegenerate-energy subspace to extract the momentum labels. Itis also possible to work the other way around, and impose acertain momentum. That is precisely what we do in the nextsection, for quantum systems on a cylinder.
FIG. 2. The lowest-lying energy levels for the quantum Isingchain on a circle with circumferences N = 20 (left) and N =50 (right). The blue dots are results from the MPS excitationansatz with D = 50 and the red crosses are results from exactdiagonalization. The energies were shifted and rescaled suchthat the ground state is at e = 0 and the gap ∆ (cid:15) = 1.We see that the MPS with open boundary conditions restorestranslation invariance, and that degeneracies agree with exactresults. spin chains – even the multi-particle excitations are wellcaptured by this ansatz. In order to make this possible,the bond dimension of the ground-state MPS needs tobe large enough such that we are in the finite-size scalingregime [25]: through the virtual level of the MPS one canchange the state over the whole system by modifying asingle tensor locally. This same effect was observed inRef. 38 for systems with open boundary conditions. Mo-tivated by this effect, we now apply our open-boundaryMPS excitation ansatz to critical models on a ring. A. Quantum Ising chain
Let us first look at the simplest critical 1-D model, thecritical Ising chain. As shown in Fig. 2, our results agreewell with exact diagonalization. Despite the open bound-ary conditions, the MPS restores translation invarianceand we retrieve momentum eigenstates. Because ourmethod does not scale exponentially in system size, wecan push the simulation to much higher system sizes.This is important for systems with a large local dimen-sion, where one often cannot reach large enough sizes toperform finite-size scaling.
B. Classical 2-D Potts model
We can also deal with statistical-mechanical problemson an infinite cylinder, by finding the leading eigenvec-tors of the transfer matrix [42] in the periodic direction.The largest few eigenvalues contain information aboutthe free energy and the dominant correlation functions
FIG. 3. The rescaled spectrum of the transfer matrix of the3-state Potts model with system size N = 28, as obtainedwith the excitation ansatz with bond dimension D = 150. Wetake the negative logarithm of the transfer-matrix eigenvalues, f = − log λ , such that f corresponds to a free energy. We haverescaled the values such that fixed-point free energy is set at f = 0 and the first excited state at f = 2 / in the system. It is possible to find the dominant eigen-vector by modifying any of the well known ground-statealgorithms to look for the largest magnitude eigenvector,and the above algorithm for finding the excited states ofthe transfer matrix are straightforwardly adapted.Here we show the finite-size spectrum of the critical 3-state Potts model. In Fig. 3 we show the spectrum of thetransfer matrix with circumference N = 28 – a systemsize that is not feasible with exact diagonalization – atbond dimension D = 150. The rescaling was done suchthat the ground state sits at free energy f = 0 and thefirst excited state at f = 2 /
15 (the smallest non-trivialscaling dimension for the Potts CFT); all other eigenval-ues are approaching the CFT prediction. In Tab. II wecompare the obtained energies with both the exact re-sults on smaller systems and the CFT prediction for theinfinite-size limit. An alternative rescaling could be per-formed, where the free energy at momentum 2 is rescaledsuch that it lies at free energy f = 2 (the first excitedstate of the identity tower), exploiting the momentuminformation explicitly. Such a rescaling is independent ofthe CFT [22], but looking at Tab. II at N = 28, it makesalmost no difference compared to the rescaling that wehave used. IV. SYSTEMS ON A CYLINDER
We can generalize this approach to the setting of MPS ap-proximations for cylindrical systems, a setup that is veryoften used for simulating 2-D systems. Here, the MPSis wrapped around the infinite cylinder in a snake-likefashion, so the translation symmetry in the transversedirection is explicitly broken by the MPS representation.However, just like for the one-dimensional rings above,
10 12 14 28 CFT0.0 0.0 0.0 0.0 0.00.13333 0.13333 0.13333 0.13333 0.133330.13333 0.13333 0.13333 0.13333 0.133330.83376 0.82863 0.82499 0.8139 0.81.148 1.1428 1.1398 1.134 1.13331.148 1.1428 1.1398 1.134 1.13331.148 1.1428 1.1398 1.134 1.13331.148 1.1428 1.1398 1.134 1.13331.3213 1.3225 1.3235 1.3277 1.33331.3213 1.3225 1.3235 1.3277 1.33331.7462 1.75 1.7537 1.7704 1.81.7462 1.75 1.7537 1.7704 1.81.8995 1.8763 1.8618 1.8272 1.81.8995 1.8763 1.8618 1.8272 1.82.0532 2.0336 2.0224 2.0016 2.02.0532 2.0336 2.0224 2.0017 2.0TABLE II. Comparison of the rescaled free energies of the3-state Potts model at different system sizes with the CFTprediction. Results for system sizes N = (10 , ,
14) wereobtained with exact diagonalization, while results for systemsize N = 28 where obtained with our MPS-based excitationansatz with bond dimension D = 150. we can hope that this translation symmetry is restoredfor large enough bond dimensions and that we can ex-ploit this to create quasiparticle excitations with fixedtransverse momentum. In the case of infinite cylinders,translation symmetry along the cylinder can be imposedstraightforwardly, so in this way we can have access tothe momentum quantum numbers in both directions. A. Method
The ground state of an infinite cylinder with a circum-ference of N sites can be represented by an infinite MPSwith an N -site unit cell | Ψ( { A i } ) (cid:105) = A . . . A N A . . . A N . . .. . . . (10)This MPS can be found using variational ground-statesearches [12, 43] or the infinite DMRG algorithm [10].We normalize the state such that the leading eigenvalueof the N -site transfer matrix is one, λ max (cid:32) A . . . ¯ A A N ¯ A N . . . (cid:33) = 1 . (11) The analysis is easily generalized to larger unit cells.
There are two translation operators acting on this state.The first one ( T x ) corresponds to translation along thecylinder, and shifts the full unit cell T x | Ψ( { A i } ) (cid:105) = A A N A . . . A N . . . . (12)The second one ( T y ) corresponds to translation aroundthe cylinder and is a transformation within the unit cell, T y | Ψ( { A i } ) (cid:105) = A A N A A N . . . . . . . (13)The MPS is invariant under T x by construction, buttranslation invariance under T y is less trivial. The ex-pectation value of this operator (cid:104) Ψ( { A i } ) | T y | Ψ( { A i } ) (cid:105) , (14)is a quantity that scales exponentially with the numberof unit cells, so the characteristic quantity is the leadingeigenvalue of the mixed N -site transfer matrix, µ = λ max ¯ A ¯ A ¯ A N A A A N . . .. . . , (15)so that, formally, (cid:104) Ψ( { A i } ) | T y | Ψ( { A i } ) (cid:105) ∼ µ N x with N x the diverging number of unit cells. For future reference,we can associate left and right fixed points to this mixedtransfer matrix as ¯ A ¯ A ¯ A N A A A N l . . .. . . = µ l (16) ¯ A ¯ A ¯ A N A A A N r. . .. . . = µ r . (17)Now the MPS has well-defined transverse momentum ifthe eigenvalue µ lies on the unit circle with the angle amultiple of 2 π/N .We can make an excitation on top of this MPS withthe ansatz | Φ p x ( B ) (cid:105) = (cid:88) n e ip x n T nx A A N A A N B (18)with B = B . . . A N + . . . · · · + A B N . . . . (19)This is the multi-site version [20] of the excitation ansatz[13] for infinite MPS. Similar to Eq. 4, we can fix the re-dundant gauge degrees of freedom in the B tensors suchthat the overlap between two of these excited states re-duces to the simple euclidean inner product with a δ -function normalization for the momentum, (cid:104) Φ p (cid:48) x ( B (cid:48) ) | Φ p x ( B ) (cid:105) = 2 πδ ( p x − p (cid:48) x ) N (cid:88) i =1 ( (cid:126)B (cid:48) i ) † (cid:126)B i . (20)Consequently, an ordinary eigenvalue problem H p x , eff (cid:126)B = ω (cid:126)B, (21)with (cid:104) Φ p (cid:48) x ( B (cid:48) ) | H | Φ p x ( B ) (cid:105) = 2 πδ ( p x − p (cid:48) x )( (cid:126)B ) † H p x , eff (cid:126)B, (22)finds the optimal B tensors for representing the lowest-energy excitation in the system.This ansatz clearly is an eigenstate of the T x operator, T x | Φ p x ( B ) (cid:105) = e ip x | Φ p x ( B ) (cid:105) , (23)but the momentum around the cylinder is not straight-forward. Indeed, we can compute the normalized ex-pectation value for the T y operator by dividing by theground-state expectation value12 πδ ( p x − p (cid:48) x ) (cid:104) Φ p (cid:48) x ( B ) | T y | Φ p x ( B ) (cid:105)(cid:104) Ψ( A ) | T y | Ψ( A ) (cid:105) = 1 µ ¯ BB rl + e − ip x µ Bl ¯ A ¯ A N ¯ BA A N r + . . . , (24)The infinite sums in this expression converge so that thisexpression yields a finite number; we can hence computethe transverse momentum of the excited states that wefind variationally. We can, however, directly target exci-tations that are approximate eigenvectors of the T y op-erators. Indeed, if we define the effective T y operatoras (cid:104) Φ p x ( B (cid:48) ) | T y | Φ p x ( B ) (cid:105)(cid:104) Ψ( A ) | T y | Ψ( A ) (cid:105) = 2 πδ ( p x − p (cid:48) x )( (cid:126)B (cid:48) ) † T eff (cid:126)B, (25)we solve the eigenvalue problem (cid:0) H p x , eff − α e − ip y T eff (cid:1) B p x ,p y = λB p x ,p y (26)for the eigenvalue λ = ω − α with the most negativereal part. When α is sufficiently large and translation It would be more elegant to consider the operator H p x , eff − α e − ip y T eff − α e ip y T † eff , (27)such that the eigenvalue problem remains hermitian. In prac-tice, however, we find that adding the hermitian conjugate isnot needed for the stability of the eigenvalue problem and onlyincreases the computational cost. FIG. 4. The momentum per rung of the MPS ground state[Eq. 15] for the square-lattice Ising model on an N = 12cylinder as a function of bond dimension for three values ofthe field, λ = 2 . λ = 2 . λ = 3 . symmetry is sufficiently well captured by the MPS, thiseigenvalue should indeed be real and correspond to aneigenstate with transversal momentum p y and energyeigenvalue ω . B. Benchmarks
We first illustrate this approach on the two-dimensionaltransverse-field Ising model, defined by the Hamiltonian H Ising = (cid:88) (cid:104) ij (cid:105) S zi S zj + λ (cid:88) i S xi . (28)We choose to run simulations at two values in thesymmetry-broken phase, λ = 2 . λ = 2 .
9, and at thequantum critical point λ = 3 . N = 12 at different bonddimensions, showing that the expectation value of T y perrung [Eq. 15] nicely converges to one. The convergence isslower as one approaches the critical point, which pointsto the fact that a larger bond dimension is needed nearcriticality.Next we study the excitations, where we can takedifferent cuts of the two-dimensional Brillouin zone byimposing the transverse momentum as in Eq. (26), seeFig. 5. The p y = 0 cut shows the dispersion becom-ing gapless at the critical point, whereas the p y = π isgapped. We push the bond dimension to get a goodestimate for the gapless point, showing that the finitecircumference induces a non-zero gap even for the criti-cal value of the field; only in the infinite 2-D plane themodel becomes truly gapless. As we show in the bottompanel of Fig. 5, the finite-size gap decreases for increasingcylinder width. FIG. 5. The dispersion relation of the square-lattice Isingmodel on a cylinder of circumference N = 12 for three valuesof the field, λ = 2 . λ = 2 . λ = 3 . p y = 0, nicely capturing the spectrum becoming gap-less at the critical point. The bottom panel is for p y = π ,where the dispersion changes very little as the field is tuned.Here, all simulations were done with an MPS bond dimension D = 256.FIG. 6. The bottom panel shows the gap at λ = 3 . N = 8 (blue) and N = 12 (red). We can convergethe value of the gap with increasing bond dimension in bothcases, although the larger cylinder needs higher bond dimen-sion. The non-zero gap is due to the finite circumference ofthe cylinder, and clearly goes down as N increases. FIG. 7. The dispersion relation of the Hubbard model with t = 1 and U = 12 on a cylinder of circumference N = 4, inthree cuts of the Brillouin zone p y = 0 (blue), p y = π/ p y = π (orange). The top panel shows the dispersion ofthe magnon excitation (charge-neutral spin-1 excitation) inthree momentum cuts, and the bottom panel shows the holedispersion (charged spin-1/2 excitation). Here, the total bonddimension is around D = 3500, where the largest SU(2) ⊗ U(1)symmetry block is D max = 154.FIG. 8. The minimum of the dispersion relation at the Spoint for the N = 4 (blue) and the N = 8 (red) cylinder as afunction of truncation threshold of the MPS. Now we apply our method to the two-dimensional Hub-bard model on the square lattice at half filling, withHamiltonian H Hubbard = − (cid:88) σ = {↑ , ↓} (cid:88) (cid:104) ij (cid:105) (cid:16) c † σ,i c σ,j + c † σ,j c σ,i (cid:17) + U (cid:88) i (cid:16) c †↑ ,i c ↑ ,i (cid:17) (cid:16) c †↓ ,i c ↓ ,i (cid:17) , (29)where c σ,i and c † σ,i are fermionic creation and annihila-tion operators with spin σ at site i . In our MPS repre-sentation for the ground state on the cylinder, we fix thefilling by implementing the U(1) symmetry for the elec-tron charge, use SU(2) spin-rotation symmetry (which isleft unbroken on a cylinder with finite circumference) andimplement a graded Z symmetry [45] for encorporatingfermionic statistics. We work at U = 12, deep in theMott-insulating phase.We optimize the ground state on infinite cylinders ofcircumference N = 4 and N = 8 using the vumps algo-rithm [43], and look at the excitation spectrum on topof the ground state. Using the U(1) and SU(2) symme-tries in the MPS representation, we can fix the chargeand spin quantum numbers ( q c , q s ) in the excited stateansatz.First, we study the magnon excitations with quantumnumbers (0 ,
1) in the top panel of Fig. 7. We observea strong minimum in the X-point at ( π, π ), in agree-ment with the fact that the model effectively behaves as aHeisenberg antiferromagnet in the charge-neutral sector.We also see a minimum at (0 , π, π )converges quickly with bond dimension, which points tothe fact that the magnon is a well-defined quasiparticlefor which our ansatz is ideally suited. Correspondingly,we observe that the energy of the two-magnon state at(0 ,
0) converges a lot slower (not shown).Next, we consider the charged excitations with quan-tum numbers ( − , ). The hole dispersion can be directlyobserved in angle-resolved photoemission spectroscopy,and has gained a lot of (renewed) theoretical attentiondue to recent cold-atom experiments [46, 47]. In particu-lar, the nature of the magnetic polaron as a bound stateof chargon and spinon [48–50] has been investigated insome detail. In the bottom panel of Fig. 7 we have plot-ted three momentum cuts for the N = 4 cyinder, and weobserve a minimum in the S point at ( π/ , π/ t - J model [50]. The p y = 0 and p y = π cuts alsoshow pronounced minima, and a level crossing at higherenergies.In Fig. 8 we plot this minimum of the hole disper-sion, i.e. the charge gap, as a function of Schmidt-valuethreshold in the MPS ground state for the N = 4 and The Schmidt-value threshold is the approximate value of the the N = 8 cylinders. We observe that convergence israther slow, which points to the composite nature of theexcitation: The quasiparticle ansatz needs a large bonddimension for describing the extended bound-state com-plex. As Fig. 8 shows, the charge gap is not completelyconverged for the N = 8 cylinder with our highest bonddimension ( D ≈ V. CONCLUSIONS
In this paper we have applied open-boundary MPS al-gorithms to systems with periodic boundary conditions.In particular, although translation symmetry is brokenby the MPS representation, we have shown that momen-tum emerges as a good quantum number. This can beexploited for labeling excitations on top of MPS groundstates.We have also introduced an ansatz for capturing ex-cited states in finite systems with open boundary condi-tions, which is a very efficient method for computing e.g.the spectral gap, and which can be straightforwardly im-plemented on top of an existing DMRG or MPS ground-state code. We expect that this approach can be of greatuse to the DMRG/MPS community.Our method was shown to work very well for obtainingCFT finite-size spectra of critical models. We expectthis will prove useful for characterizing the effective CFTfor 1-D quantum and 2-D classical models, as well as tocompute entanglement spectra of projected entangled-pair states in order to characterize the (chiral) edge fieldtheory [51].Applying the framework to 2-D cylinders has shownthat we can obtain excitation spectra of challenging mod-els such as spin liquids or (doped) Hubbard models. In-deed, whereas simulating time evolution is often pro-hibitively expensive for these models and excitation spec-tra are, therefore, rather infeasible to compute, the exci-tation ansatz can be applied with a cost that is similarto a ground-state simulation. In that respect, it wouldbe extremely interesting to look at fractionalized excita-tions such as spinons or holons: In principle, the excita-tion ansatz can be straightforwardly generalized to alsocapture these topological excitations [20, 52]. Another in-teresting question is how this MPS excitation ansatz onthe cylinder compares to the excitation ansatz for pro-jected entangled-pair states [53–55], which is formulateddirectly on the infinite plane. smallest Schmidt value in all symmetry sectors of the MPS; it isthe analog of the truncation threshold in (i)DMRG simulations,which does not have a clear meaning in variational optimizationschemes without truncation steps.
VI. ACKNOWLEDGEMENTS
The authors would like to thank Ji-Yao Chen, Boris Pon-sioen, Philippe Corboz, Michael Knap and Frank Poll- mann for inspiring discussions. This work was supportedby the Research Foundation Flanders (G0E1520N,G0E1820N) and the ERC grants QUTE (647905) andERQUAF (715861). [1] U. Schollw¨ock, The density-matrix renormalization groupin the age of matrix product states, Annals of Physics , 96 (2011).[2] J. I. Cirac, D. Perez-Garcia, N. Schuch, and F. Ver-straete, Matrix product states and projected entan-gled pair states: Concepts, symmetries, and theorems,arXiv:2011.12127 (2020).[3] S. R. White, Density matrix formulation for quantumrenormalization groups, Phys. Rev. Lett. , 2863 (1992).[4] S. R. White, Density-matrix algorithms for quantumrenormalization groups, Phys. Rev. B , 10345 (1993).[5] S. Rommer and S. ¨Ostlund, Class of ansatz wave func-tions for one-dimensional spin systems and their relationto the density matrix renormalization group, Phys. Rev.B , 2164 (1997).[6] F. Verstraete, D. Porras, and J. I. Cirac, Density matrixrenormalization group and periodic boundary conditions:A quantum information perspective, Phys. Rev. Lett. ,227205 (2004).[7] D. Perez-Garcia, F. Verstraete, M. M. Wolf, andJ. I. Cirac, Matrix product state representations,arXiv:quant-ph/0608197 (2006).[8] S. ¨Ostlund and S. Rommer, Thermodynamic limit of den-sity matrix renormalization, Phys. Rev. Lett. , 3537(1995).[9] G. Vidal, Classical simulation of infinite-size quantumlattice systems in one spatial dimension, Phys. Rev. Lett. , 070201 (2007).[10] I. P. McCulloch, Infinite size density matrix renormaliza-tion group, revisited, arXiv:0804.2509 (2008).[11] J. Haegeman, J. I. Cirac, T. J. Osborne, I. Piˇzorn, H. Ver-schelde, and F. Verstraete, Time-dependent variationalprinciple for quantum lattices, Phys. Rev. Lett. ,070601 (2011).[12] L. Vanderstraeten, J. Haegeman, and F. Verstraete,Tangent-space methods for uniform matrix productstates, SciPost Phys. Lect. Notes , 7 (2019).[13] J. Haegeman, B. Pirvu, D. J. Weir, J. I. Cirac, T. J. Os-borne, H. Verschelde, and F. Verstraete, Variational ma-trix product ansatz for dispersion relations, Phys. Rev.B , 100408 (2012).[14] R. P. Feynman, Atomic theory of the λ transition in he-lium, Phys. Rev. , 1291 (1953).[15] S. M. Girvin, A. H. MacDonald, and P. M. Platzman,Collective-excitation gap in the fractional quantum halleffect, Phys. Rev. Lett. , 581 (1985).[16] D. P. Arovas, A. Auerbach, and F. D. M. Haldane, Ex-tended heisenberg models of antiferromagnetism: Analo-gies to the fractional quantum hall effect, Phys. Rev.Lett. , 531 (1988).[17] M. Takahashi, Excitation spectra of s=1 antiferromag-netic chains, Phys. Rev. B , 3045 (1994).[18] E. S. Sørensen and I. Affleck, Equal-time correlations inhaldane-gap antiferromagnets, Phys. Rev. B , 15771(1994). [19] A. K. Bera, B. Lake, F. H. L. Essler, L. Vanderstraeten,C. Hubig, U. Schollw¨ock, A. T. M. N. Islam, A. Schnei-dewind, and D. L. Quintero-Castro, Spinon confinementin a quasi-one-dimensional anisotropic heisenberg mag-net, Phys. Rev. B , 054423 (2017).[20] V. Zauner-Stauber, L. Vanderstraeten, J. Haegeman,I. P. McCulloch, and F. Verstraete, Topological nature ofspinons and holons: Elementary excitations from matrixproduct states with conserved symmetries, Phys. Rev. B , 235155 (2018).[21] L. Vanderstraeten, M. Van Damme, H. P. B¨uchler, andF. Verstraete, Quasiparticles in quantum spin chains withlong-range interactions, Phys. Rev. Lett. , 090603(2018).[22] P. Francesco, P. Mathieu, and D. S´en´echal, Conformalfield theory (Springer Science & Business Media, 2012).[23] D. Porras, F. Verstraete, and J. I. Cirac, Renormaliza-tion algorithm for the calculation of spectra of interactingquantum systems, Phys. Rev. B , 014410 (2006).[24] B. Pirvu, F. Verstraete, and G. Vidal, Exploiting trans-lational invariance in matrix product state simulationsof spin chains with periodic boundary conditions, Phys.Rev. B , 125104 (2011).[25] B. Pirvu, G. Vidal, F. Verstraete, and L. Tagliacozzo,Matrix product states for critical spin chains: Finite-size versus finite-entanglement scaling, Phys. Rev. B ,075117 (2012).[26] D. Draxler, J. Haegeman, F. Verstraete, and M. Rizzi,Continuous matrix product states with periodic bound-ary conditions and an application to atomtronics, Phys.Rev. B , 045145 (2017).[27] Y. Zou, A. Milsted, and G. Vidal, Conformal dataand renormalization group flow in critical quantum spinchains using periodic uniform matrix product states,Phys. Rev. Lett. , 230402 (2018).[28] P. Pippan, S. R. White, and H. G. Evertz, Efficientmatrix-product state method for periodic boundary con-ditions, Phys. Rev. B , 081103 (2010).[29] D. Rossini, V. Giovannetti, and R. Fazio, Stiffness in 1dmatrix product states with periodic boundary conditions,Journal of Statistical Mechanics: Theory and Experi-ment , P05021 (2011).[30] B. Pirvu, J. Haegeman, and F. Verstraete, Matrix prod-uct state based algorithm for determining dispersion re-lations of quantum spin chains with periodic boundaryconditions, Phys. Rev. B , 035130 (2012).[31] W.-L. Tu, H.-K. Wu, N. Schuch, N. Kawashima, andJ.-Y. Chen, Generating function for tensor network di-agrammatic summation, arXiv:2101.03935 (2021).[32] E. M. Stoudenmire and S. R. White, Studying two-dimensional systems with the density matrix renormal-ization group, Annual Review of Condensed MatterPhysics , 111 (2012).[33] J. Motruk, M. P. Zaletel, R. S. K. Mong, and F. Poll-mann, Density matrix renormalization group on a cylin- der in mixed real and momentum space, Phys. Rev. B , 155139 (2016).[34] G. Ehlers, S. R. White, and R. M. Noack, Hybrid-space density matrix renormalization group study of thedoped two-dimensional hubbard model, Phys. Rev. B ,125125 (2017).[35] M. Van Damme, G. Roose, M. Hauru, and J. Haegeman,Mpskit.jl (2020).[36] U. Schollw¨ock, The density-matrix renormalizationgroup, Rev. Mod. Phys. , 259 (2005).[37] I. P. McCulloch, From density-matrix renormalizationgroup to matrix product states, Journal of Statistical Me-chanics: Theory and Experiment , 10014 (2007).[38] N. Chepiga and F. Mila, Excitation spectrum and densitymatrix renormalization group iterations, Phys. Rev. B , 054425 (2017).[39] L. Vanderstraeten, F. Verstraete, and J. Haegeman, Scat-tering particles in quantum spin chains, Phys. Rev. B ,125136 (2015).[40] S. R. White and D. A. Huse, Numerical renormalization-group study of low-lying eigenstates of the antiferromag-netic s=1 heisenberg chain, Phys. Rev. B , 3844 (1993).[41] E. S. Sørensen and I. Affleck, Large-scale numerical evi-dence for bose condensation in the s=1 antiferromagneticchain in a strong field, Phys. Rev. Lett. , 1633 (1993).[42] J. Haegeman and F. Verstraete, Diagonalizing transfermatrices and matrix product operators: A medley of ex-act and computational methods, Annual Review of Con-densed Matter Physics , 355 (2017).[43] V. Zauner-Stauber, L. Vanderstraeten, M. T. Fishman,F. Verstraete, and J. Haegeman, Variational optimizationalgorithms for uniform matrix product states, Phys. Rev.B , 045145 (2018).[44] H. W. J. Bl¨ote and Y. Deng, Cluster monte carlo sim-ulation of the transverse ising model, Phys. Rev. E ,066110 (2002).[45] N. Bultinck, D. J. Williamson, J. Haegeman, andF. Verstraete, Fermionic matrix product states and one-dimensional topological phases, Phys. Rev. B , 075108(2017).[46] C. S. Chiu, G. Ji, A. Bohrdt, M. Xu, M. Knap, E. Demler,F. Grusdt, M. Greiner, and D. Greif, String patterns inthe doped hubbard model, Science , 251 (2019).[47] A. Bohrdt, C. S. Chiu, G. Ji, M. Xu, D. Greif, M. Greiner,E. Demler, F. Grusdt, and M. Knap, Classifying snap-shots of the doped Hubbard model with machine learn-ing, Nature Physics , 921 (2019).[48] P. B´eran, D. Poilblanc, and R. Laughlin, Evidence forcomposite nature of quasiparticles in the 2d t-j model,Nuclear Physics B , 707 (1996).[49] R. B. Laughlin, Evidence for quasiparticle decay in pho-toemission from underdoped cuprates, Phys. Rev. Lett. , 1726 (1997).[50] A. Bohrdt, E. Demler, F. Pollmann, M. Knap, andF. Grusdt, Parton theory of angle-resolved photoemis-sion spectroscopy spectra in antiferromagnetic mott in-sulators, Phys. Rev. B , 035139 (2020).[51] D. Poilblanc, N. Schuch, and I. Affleck, SU(2) chiraledge modes of a critical spin liquid, Phys. Rev. B ,174414 (2016).[52] L. Vanderstraeten, E. Wybo, N. Chepiga, F. Verstraete,and F. Mila, Spinon confinement and deconfinement inspin-1 chains, Phys. Rev. B , 115138 (2020). [53] L. Vanderstraeten, M. Mari¨en, F. Verstraete, andJ. Haegeman, Excitations and the tangent space of pro-jected entangled-pair states, Phys. Rev. B , 201111(2015).[54] L. Vanderstraeten, J. Haegeman, and F. Verstraete, Sim-ulating excitation spectra with projected entangled-pairstates, Phys. Rev. B , 165121 (2019).[55] B. Ponsioen and P. Corboz, Excitations with projectedentangled pair states using the corner transfer matrixmethod, Phys. Rev. B , 195109 (2020). Appendix A: Pseudocode for excitation ansatz withopen-boundary MPS
In this subsection we elaborate on the implementationdetails of the quasiparticle ansatz for finite systems. Afinite MPS can be gauged in such a way that every tensorleft and right from a given site is respectively left/rightisometric | Ψ( A . . . A N (cid:105) = A l . . . A ci . . . A rN A l . (A1)In what follows we will work in the MPO-representationof the Hamiltonian [1, 37]. We will need the partiallycontracted environments of the state, much like in exist-ing DMRG codes: ρ ln +1 = ¯ A ln A ln ρ ln O , ρ rn − = A rn ¯ A rn ρ rn O . (A2)We will also require V i , the null space of A li , such that ¯ A ln V n = 0 . (A3)All degrees of freedom can then be absorbed in tensors X i , as we have seen earlier B i = V i X i , | Φ i ( X i ) (cid:105) = A l . . . V i X i . . . A rN . (A4)The last two quantities we need to define are ρ B,ln and ρ B,rn , satisfying: ρ B,ln +1 = ¯ A ln B n ρ ln O + ¯ A ln A rn ρ B,ln O (A5)and ρ B,rn − = B n ¯ A rn ρ rn O + A ln ¯ A rn ρ B,rn O . (A6)1All these quantities together make H eff ( (cid:126)X ) take on apleasant form. For every site i , one has to calculate thefollowing tensor contractions to get the derivative withrespect to the tensor B (cid:48) i : T = A ri ρ B,li O ρ ri + B i ρ li O ρ ri + A li ρ li O ρ B,ri . (A7)The fundamental degrees of freedom are the X (cid:48) i tensors,which can then be extracted by projecting down on V i : X (cid:48) i = ¯ V i T . (A8)The pseudocode for H eff ( (cid:126)X ) can be found in Algorithm1. Eigenvectors can be found using any iterative eigen-solver. This algorithm is very simple, straightforward toimplement and requires similar contractions as in ‘usual’DMRG codes. Not much changes when going to the ther-modynamic limit, except the calculation for ρ B,ln , ρ B,rn involves infinite sums [12].
Algorithm 1
Pseudocode for the action of H eff ( (cid:126)X ) Inputs ( (cid:126)X, (cid:126)ρ l , (cid:126)ρ r ) (cid:126)X (cid:48) ← (cid:126)ρ B,l ← (cid:126)ρ B,r ← (cid:46) initialization for i ∈ [1 , len ( (cid:126)X )] do (cid:46) environments calculate ρ B,li +1 (A5) calculate ρ B,rend − i (A6) for i ∈ [1 , len ( (cid:126)X )] do calculate T (A7) calculate X (cid:48) i (A8) return (cid:126)X(cid:126)X