Efficient Basis Formulation for (1+1)-Dimensional SU(2) Lattice Gauge Theory: Spectral calculations with matrix product states

We propose an explicit formulation of the physical subspace for a (1+1)-dimensional SU(2) lattice gauge theory, where the gauge degrees of freedom are integrated out. Our formulation is completely general, and might be potentially suited for the design of future quantum simulators. Additionally, it allows for addressing the theory numerically with matrix product states. We apply this technique to explore the spectral properties of the model and the effect of truncating the gauge degrees of freedom to a small finite dimension. In particular, we determine the scaling exponents for the vector mass. Furthermore, we also compute the entanglement entropy in the ground state and study its scaling towards the continuum limit.


I. INTRODUCTION
Gauge theories play a central role in our understanding of modern particle physics, with the Standard Model being one of the most prominent examples. In the usual Lagrangian or Hamiltonian formulation, the local gauge symmetry is ensured by introducing additional degrees of freedom (d.o.f.) in the form of a gauge field. However, this also leads to redundant d.o.f. in the theory. As the physical observables are strictly gauge invariant, the only relevant subspace is the one spanned by the gauge-invariant states, which is in general much smaller than the full Hilbert space of the theory. Because of the absence of transversal directions for the special case of (1+1) dimensions, the gauge fields are not genuinely independent d.o.f. Therefore, it is possible to remove them by integrating the Gauss law. This long-known fact, to the best of our knowledge, has been explicitly exploited in practice only for the Abelian case of the Schwinger model [1][2][3].
Despite these encouraging prospects, there are also some limitations. In particular, the Hilbert spaces for the gauge d.o.f. are typically infinite dimensional. Hence, in cases for which the gauge fields cannot be integrated out, they typically have to be truncated to a finite dimension to allow for a TN approach or a potential implementation in a quantum simulator. Previous works therefore resorted to the truncation methods from Refs. [9,50] to achieve a finite dimension while simultaneously preserving gauge invariance. A different type of finitedimensional gauge models explored in that context are quantum link models [51][52][53], where the gauge d.o.f. are replaced by discrete spins. However, these truncated models do not necessarily correspond to the continuum theory in the limit of vanishing lattice spacing, or might not have a continuum limit at all [54].
Here, we address these questions for a (1+1)dimensional SU(2) lattice gauge theory. In a first step, we show how, starting from a color-neutral basis developed in Ref. [55], the gauge d.o.f. can be integrated out on a lattice with open boundary conditions (OBC). The resulting formulation allows for truncating the colorelectric flux at an arbitrary value. These truncated models can be efficiently addressed with TN. In principle, since the maximum flux on a finite lattice with OBC is upper bounded, this enables an exact treatment of the arXiv:1707.06434v2 [hep-lat] 24 Jan 2018 model. In practice, to solve the model with matrix product states (MPS), we limit the number of color-flux sectors and compute the low-lying spectrum for this family of truncated SU(2) gauge models and investigate truncation effects in a systematic manner. In particular, here we explore such effects in the closing of the mass gaps as we approach criticality, and in the entanglement entropy of the ground state. Recent developments in the context of the holographic principle have suggested a deep connection between entanglement and emergent geometry [56,57], and have rekindled interest in understanding the peculiarities of entanglement in gauge theories [58][59][60][61][62]. Being especially well suited to compute entanglement entropies, TNs allow us to study how the truncation of the flux alters the entanglement of the vacuum in the approach to the continuum.
In our study, we are interested in various aspects. On the one hand, the basis we develop efficiently describes the physical subspace and can in principle be used with other analytical or numerical methods. Using TNs, we demonstrate the suitability of our formulation for addressing LGT problems as they appear in high-energy and also condensed-matter physics. Because of the vastly reduced number of basis states compared to the full basis, we are able to explore much larger values for the maximum color-electric flux than has been achieved in previous studies of the model [25,28]. On the other hand, the questions we explore are also relevant for quantum simulation of gauge theories. The corresponding, Abelianized, Hamiltonian in our basis is nonlocal, similar to the one recently realized in trapped ions for the Schwinger model [49]. Hence, our formulation might have potential applications for the design of future quantum simulators. Some proposals for quantum simulation of gauge models [9,[33][34][35][36][37][38][39][40][41][42][43][44][45][46]48] rely on the representation of the gauge variables by finite-dimensional d.o.f. For those cases, a truncation of the color flux will be required, and our formulation provides a tool for the systematic study of the effect on various observables.
The rest of the paper is organized as follows. In Sec. II, we introduce the model we are studying. After a brief review of the color-neutral basis developed in Ref. [55], we present our new formulation for systems with OBC where the gauge field is integrated out in Sec. III. Furthermore, we explain how this formulation readily allows for a truncation of the link Hilbert spaces in a gaugeinvariant manner to a finite dimension. In Sec. IV, we briefly review the MPS methods we are applying and present our results for the low-lying spectrum and the entanglement properties of the ground state while approaching the continuum limit. Finally, we conclude in Sec. V.

II. MODEL
The model we are studying is a (1+1)-dimensional SU(2) lattice gauge theory. We use a Hamiltonian lattice formulation with Kogut-Susskind staggered fermions [63] in the temporal gauge, given by (1) In the expression above, ψ k † is a single-component fermionic field creating a fermion of color on site k, U k acts on the gauge link between sites k and k + 1, and J 2 k gives the color-electric energy on the link. The parameter g is the coupling constant, a the lattice spacing, and m the bare fermion mass.
The operators U k are SU(2) matrices in the fundamental representation and can be interpreted as rotation matrices. Hence, the Hilbert space for each gauge link is analogous to a quantum rigid rotor with total angular momentum j, which can be described in two reference frames, the body-fixed system and the space-fixed (inertial) frame of reference [63]. Consequently, the links can be labeled by the angular momentum z components of the rotor , , one corresponding to the body-fixed coordinate system and one corresponding to the space-fixed coordinate system and the total angular momentum j (the same in both reference frames). The angular momentum operators L τ , τ ∈ {x, y, z} (for the body-fixed reference frame) and R τ (for the space-fixed reference frame) can be interpreted as the left and right electric field on a link and they are related to the color-electric flux energy as J 2 k = τ L τ k L τ k = τ R τ k R τ k . Hence, the operator for the color-electric flux in this basis is simply a total angular momentum operator with eigenvalues j(j+1), j = 0, 1/2, 1, . . . . From the considerations above, a suitable basis for addressing the Hamiltonian Eq. (1) is given by |n 1 , n 2 ⊗ |j ⊗ |n 1 , n 2 ⊗ . . . where n is the fermionic occupation number for color .
The physical states |φ of the system have to satisfy the Gauss law G τ k |φ = 0, ∀k, τ , where are the generators for gauge transformations. In the formula above, Q τ k = 2 =1 1 2 ψ k † σ τ ψ k are the components of the non-Abelian charge at site k and σ τ are the usual Pauli matrices. For the following spectral calculations it is convenient to use a dimensionless formulation of Hamiltonian Eq.
The adimensional parameters of the problem in units of the coupling g are x = 1/(ag) 2 and µ = 2 √ xm/g. In the strong coupling limit, the hopping term can be neglected and the Hamiltonian can be solved analytically. The gauge-invariant ground state is then simply given by the lattice analog of the Dirac sea corresponding to odd sites occupied by a fermion of each color, empty even sites, and vanishing color flux on the links [55]: In the formula above, the numbers in bold face represent the fermionic occupation numbers while |000 represents a link carrying no flux.
The basis considered in this paragraph still contains all the information about the color d.o.f. and, in particular, states which are not color singlets. However, the eigenstates for any physical observable consist of color-neutral superpositions of basis states satisfying the Gauss law. As we show in the next paragraph, restricting oneself to the physically relevant subspace of these color-neutral superpositions allows for significantly reducing these superfluous d.o.f.

III. INTEGRATING OUT THE GAUGE FIELD
A. Color-neutral basis A first step towards a physical basis was made by Hamer in the context of a strong coupling expansion of the model [55]. Here, we briefly review the basis formulation developed there. As shown in Refs. [55,64], the physically relevant states can be generated by applying the operator V from Eq. (3) repeatedly to a certain colorneutral initial state having the desired quantum numbers. This operator has no uncontracted color indices, thus, it can locally only generate or annihilate excitations consisting of color-neutral superpositions of quark-antiquark (antiquark-quark) pairs connected by a color-flux string, as illustrated in Figs. 1(a) and 1(b). The resulting superposition has a well-defined value of j on the links and fermionic occupation number n k = n 1 k + n 2 k . In particular, applying V to such a color singlet characterized by n k , n k+1 , and j k results, in general, again in a superposition of different color singlets with n k , n k+1 and j k . In Ref. [55], all possible transitions were worked out and the matrix elements of the operator V for each of those vertices are shown in Fig. 1(c).
Looking at Eq. (4), one can easily see that the states generated in that manner are eigenstates of the mass term, as it depends only on the total occupation number, and of the color-electric energy, as it only depends on j. Thus, W 0 acts identically on all those states. Consequently, instead of working in the basis containing the full color information, we can restrict ourselves to a basis formed by those color-singlet states, characterized by the fermionic occupation number of each site and by the color-electric flux j carried by a link, |φ = |n ⊗ |j ⊗ |n ⊗ |j · · · .
Here, n ∈ {0, 1, 2}, as we are not distinguishing between fermions of different colors anymore. Gauss law in this color-singlet basis simply reduces to the fact that the electric flux j k+1 can differ only from the one on the previous link by one quantum, if the site is occupied by a single fermion: In addition to reducing the d.o.f. significantly compared to the full basis, the color-singlet basis also offers the possibility to trivially truncate the color-electric flux at a certain value of j max in a gauge-invariant manner. Taking into account only states with j ≤ j max results in a truncated model with Hilbert spaces of dimension d link = 2j max + 1 for the gauge links. Compared to the full basis, where for j max = 1/2, 1, 3/2, 2 one would have link Hilbert spaces of dimension 5, 14, 30, 55, one only has to deal with spaces of dimension d link = 2, 3, 4, 5.

B. Removing the gauge fields
The color-singlet basis reviewed in the previous paragraph still contains redundant information which can be (partially) removed. While the discussion above still applies to both open and periodic boundary conditions, we restrict ourselves from now on to OBC.
Realizing that due to Eq. (5) the flux j k+1 depends only on j k and n k+1 , one can reconstruct the colorelectric flux at every link from the value j 0 at the left boundary, by recursive application of the Gauss law. However, the case of a site occupied by a single fermion is ambiguous, as the SU(2) non-Abelian Gauss law allows a change of the color-electric flux by ±1/2. To lift this ambiguity, we use two states |1 + , |1 − to describe the singly occupied site, which encode if the electric flux is increasing (|1 + ) or decreasing (|1 − ) with respect to the link to the left. As a result, the basis for a single fermionic site is again four dimensional and consists of the set of states {|0 , |1 − , |1 + , |2 }. The expense of increasing the basis for the fermionic sites by one allows us to integrate out the gauge links from the Hamiltonian, as the j k can now be reconstructed solely from the fermionic content viâ where j 0 is the color-flux value at the left boundary. Hence, a suitable basis for a system with N sites is given by with |α k ∈ {|0 , |1 − , |1 + , |2 }. Similar to the (Abelian) Schwinger model, in this basis the configuration of the sites uniquely determines the content of the gauge links, thus, effectively Abelianizing the model. Additionally, one can immediately see that this construction leads to long-range interactions in Hamer's color-electric energy term. Moreover, as some of the matrix elements for the hopping term in Fig. 1(c) also depend on the colorelectric flux, the hopping term becomes nonlocal, too (details about the Hamiltonian in this basis are given in Appendix B). It is instructive to study the dimension of the physical subspace in our basis. Without further constraint, it still contains unphysical states implying negative values of j k ; e.g., |φ = |1 − ... implies a value of j 1 = −1/2 [65]. For vanishing background field, j 0 = 0, the case on which we focus in our numerical calculations, physical basis states are characterized by a simple condition: the number l k,− of sites with |1 − up to a site k can never exceed the corresponding number l k,+ of sites with |1 + , l k,+ ≥ l k,− ∀k = 1, . . . , N . The dimension of the physically relevant subspace fulfilling this condition is given by 4 N (1 − N k=1 C k /4 k ), where C k = (2k)!/(k + 1)!k! is the Catalan number (for details, see Appendix A). Compared to formulations for the physical subspace for the U(1) case with dynamical fermions [3,8], we observe that the number of basis states in our formulation is exactly the square.
A simple isometry maps the states of the reduced basis from Eq. (7) to the full one. This transformation, explicitly shown in Appendix F, sequentially reconstructs the color flux on each link from the fermion content, and prepares a suitable combination of states in the full basis such that the state is a color singlet. Thus, the map can be written as a sequence of isometries M loc , schematically shown in Fig. 2, that act from left to right and take as input one fermionic site and the corresponding incoming flux link, and expand the basis to include the outgoing link, too (see Appendix F for more details). From a quantum information point of view, this is simply a quantum circuit of depth equal to the system size [66].
Notice that the reduced basis formulation is completely general and contains the entire information about the physical subspace. Hence, it lends itself to any analytical or numerical method. Moreover, one could, in principle, treat arbitrary gauge groups SU(N c ) with N c ≥ 2 in a similar fashion. After obtaining the corresponding matrix elements for the vertices, the gauge field for colorsinglet states can likewise be encoded in the fermionic sites.

IV. RESULTS
While the formulation derived in Sec. II can be approached with any numerical method, it is particularly suited for MPS. In order to demonstrate that, we compute the ground state and the vector state of the theory at vanishing background field, j 0 = 0, in the sector of vanishing total charge, j N = 0. Moreover, in contrast to Monte Carlo methods, the MPS approach also allows for access to the entanglement entropy, and we can also study the scaling of the von Neumann entanglement entropy in the ground state while approaching the continuum limit.

A. Numerical methods
For our numerical simulations we use the MPS Ansatz with OBC [67] which is given by for a system of N sites. In the formula above the states {|i k d k=1 } form a basis for the d-dimensional Hilbert space on site k, A i k k ∈ C D×D for 1 < k < N and A i1 The bond dimension of the MPS, D, determines the number of variational parameters in the Ansatz and limits the amount of entanglement that can be present in the state (for detailed reviews about MPS, see, e.g., Refs. [68][69][70]).
Using standard methods [71], we can determine MPS approximations for the ground-state as well as for lowlying excitations. The ground-state approximation is obtained by variationally minimizing the energy with respect to the tensor A i k k , while keeping the others fixed and iterating these updates from left to right and back, until the relative change of the energy per sweep is below a certain tolerance ε. The optimal tensor in every step is found by solving an eigenvalue problem [72,73] for an effective Hamiltonian describing the interaction of site k with its environment. Excited states can be computed by projecting the Hamiltonian onto a subspace orthogonal to each of the previously computed states and then using the same variational method with the projected Hamiltonian [15]. In the continuum, the vector candidate is the lowest-lying zero momentum mesonic excitation of the ground state with charge conjugation quantum number −1 and parity −1. On the finite lattices with OBC we are working with, however, charge conjugation as well as the momentum are no longer good quantum numbers due to the broken translational invariance. Nonetheless, the remnants of these symmetries allow us to properly identify the vector state (see Appendix D for details).
In addition, to be able to address the problem with MPS, we have to express the Hamiltonian as a matrix product operator (MPO) [74]. This can be done exactly with the bond dimension of the MPO representation growing linearly with the maximum color-electric flux j max present in the system (for details, see Appendix B). In the sector j 0 = 0 = j N , on which we focus in our calculations, j max is upper bounded by N/2 × 1/2 for a system with N sites. For the system sizes we are interested in, this would lead to a very large computational effort and, hence, we truncate j max to smaller values. In particular, this allows us to explore the effects of truncating the gauge d.o.f. to a finite dimension as might be necessary for a potential future quantum simulator [9,38,39,42]. Taking advantage of the fact that our basis formulation allows for an efficient truncation, we can easily reach maximum values for the color-electric flux far beyond j max = 1/2 amenable in previous numerical studies with TNs [25,28].
In order to avoid any influence of the unphysical states, one could in principle implement the symmetries directly at the level of the tensors [7,8]. In our calculations, we are targeting only the low-lying spectrum; hence, we choose a simpler approach and remove the unphysical states by adding an energy penalty to the Hamiltonian. Moreover, as we are only interested in the vector excitations, we also remove possible baryonic states from the low-lying spectrum with an additional penalty term. A third penalty ensures that we are in the sector of vanishing total charge (see Appendix C for more details on the penalty terms).

B. Low-lying spectrum
In order to demonstrate the power of the basis developed in Sec. II for TN calculations, we compute the ground state and the vector state of the model for a range of masses, m/g ∈ [0.1, 1.6]. To probe for truncation effects, we explore a family of models with maximum colorelectric flux j max = 1/2, 1, 3/2, 2. Moreover, we consider for each combination of (m/g, j max ) several system sizes N ∈ [100, 200] with lattice spacings corresponding to x ∈ [50, 150] to be able to extrapolate to the continuum. Compared to a conventional lattice calculation, we have an additional source of error due to the limited bond dimension that can be reached in the numerical simulations. To control this error, we repeat the calculation for each combination of parameters (m/g, j max , x, N ) for several bond dimensions D ∈ [50,200]. To estimate the exact ground-state energy values, E 0 (N, x), and vector mass gap values, ∆ vec (N, , we first extrapolate our data to the limit 1/D → 0, as illustrated in Figs. 3(a) and 3(b) (details about the extrapolation procedure are given in Appendix E). Subsequently, we can proceed in a standard manner and estimate the continuum values by first extrapolating to the thermodynamic limit and then to the limit of vanishing lattice spacing.

Figures 3(c) and 3(d)
show examples for the extrapolation to the thermodynamic limit for j max = 2. Even for this case with the largest color-electric flux, for which we expect the error due to the finite bond dimension to be most pronounced, the error bars resulting from the extrapolation in D are small and we can obtain precise estimates for the ground-state energy density and vector mass gap in the thermodynamic limit.
In the final step, we extrapolate the data obtained in the previous step to the continuum limit ag = 1/ √ x → 0 by fitting a polynomial in ag. To estimate our systematic error, we compare different fits up to quadratic order using different ranges of ag (details about the error estimation procedure are given in Appendix E). For the ground state, we observe that, in general, lattice effects are well under control, independently of the truncation, and we can reliably extrapolate to the continuum limit, as can be seen in Fig. 4. The values obtained for the groundstate energy density, shown in Fig. 5, deviate at most at the percentage level from the result for the continuum theory, even for the simplest nontrivial truncation j max = 1/2. For larger j max our data are closer to the analytic solution for the untruncated lattice Hamiltonian in the limit ag = 1 √ x → 0, especially for smaller masses. In particular, the data obtained for the largest two values of j max show hardly any difference. The dip around m/g = 0.35 for j max = 3/2, 2 is due to the fact that for m/g < 0.4 our estimates for the ground-state energy are lower than the exact results, whereas for m/g ≥ 0.4 we obtain values slightly above the analytical prediction. For j max = 1/2, 1 the values are consistently larger than the exact continuum solution; hence, in these cases there is no dip. Looking at the continuum extrapolation for the vector mass gaps, we observe a noticeably different behavior. While the ground-state energy densities do not show any significantly different behavior for small j max , the vector mass gaps do, as Fig. 6 reveals. In particular, for j max = 1/2 our data suggest that higher than quadratic order corrections in ag are still relevant, which results in large χ 2 d.o.f. in our fits. With the range of lattice spacings available, we do not seem to control lattice effects well enough to give a reliable error estimate for that case. For j max = 1, our data are reasonably well described with a quadratic function in ag (χ 2 d.o.f. around 1); nevertheless, the error of the continuum estimate is still large. In contrast, for j max = 3/2 and 2, quadratic corrections in ag are irrelevant and our data are well described by a linear fit in the range of lattice spacings we study (for details, see Appendix E). The final results obtained after the extrapolation for various masses and truncations are shown in Fig. 7. As the figure reveals, there is a considerable difference between the values obtained for j max = 1/2, 1, 3/2. Only for values of j max ≥ 3/2 our data agree well with the numerical results obtained by a strong coupling expansion in Ref. [55]. In particular, for the largest mass, m/g = 1.6, the data for j max = 3/2, 2 and from Ref. [55] are already close to the nonrelativistic limit, m/g → ∞, for which the vector mass gap is given by ∆ vec = 2m/g. On the contrary, the values obtained for j max = 1/2, 1 severely differ from the nonrelativistic prediction.
The fact that we do not recover the continuum limit for the full theory for j max = 1/2, 1 might be due to several reasons. On the one hand, the truncation to a small value of j max might lead to enhanced lattice effects. While the extrapolations to the bulk limit are in general unproblematic for all truncations we study, for the continuum extrapolations we observe that higher than linear order corrections are relevant for small j max , whereas this is not the case for j max = 3/2, 2. This could indicate that one would need smaller lattice spacings for truncations to a small color-electric flux, to control lattice effects properly. On the other hand, this might be a hint that the continuum limit for small j max does not exist, similar to quantum link models [53,54]. These different types of truncated gauge models, in which the gauge links are replaced by spins, are known to approach the continuum limit by dimensional reduction of an extra dimension.
Contrary to the Abelian Schwinger model, in the SU(2) case the chiral symmetry is restored in the limit m/g → 0 and, hence, the vector mass gap goes to zero and the theory becomes critical. From our data we can extract the critical exponent for the vector mass gap, too. In order to obtain the critical exponent, we fit our data For completeness, we also show the data for jmax = 1/2 (light gray circles), although in this case our lattice spacings do not allow for a reliable estimate.
to a power law γ(m/g) ν in the region of small masses m/g ≤ 0.4, for which the model is still far away from the nonrelativistic limit. The final results obtained for the critical exponents are shown in Table I  our estimates for the critical exponents essentially agree, within error bars, with 2/3, obtained for the large N c limit of the model [75]. The central value for j max = 1 is not so close to 2/3; nevertheless, within the relatively large error bars it is still compatible. In the case of the simplest nontrivial truncation, j max = 1/2, a fit to our data yields 0.639 for the critical exponent. However, due to the large lattice effects in the vector mass gaps, the value is not trustworthy and we cannot reliably estimate the uncertainty, which we expect to be large, too. These observations also have important implications for a potential future quantum simulator. Our data show that the ground-state energy densities obtained from our family of truncated models already give a good estimate for the one of the full model in the continuum, even with the simplest nontrivial truncation. In contrast, the vector state is much more sensitive to truncation effects.
Although with our data it is not possible to fully rule out that for small j max one suffers from enhanced lattice effects, they might indicate that the model does not have a proper continuum limit in those cases.

C. Entanglement entropy
There is renewed interest in understanding the structure of entanglement in the gauge-invariant scenario, motivated in part by a deep connection between entanglement and space-time geometry that has been suggested in the context of the gauge-gravity duality [56,57]. Our reduced gauge-invariant formulation, together with TN techniques, allows us not only to determine the mass spectrum of the theory, but also to compute the entanglement entropy of the corresponding ground state and to analyze the behavior of this quantity as we approach the continuum limit, as well as to tell to what extent these features are sensitive to a truncation in the gauge d.o.f.
Recently it has been shown that for gauge theories the reduced density matrix (RDM) for a subsystem can be written as a direct sum of terms supported on sectors corresponding to different flux configuration of the boundary links [59][60][61]. Specifically for the (1+1)-dimensional case, we can decompose the RDM for the leftmost L sites and links as ρ = ⊕ jρj , where j labels the flux on the Lth link. Hence, the von Neumann entropy can be written as Here, ρ j is the (normalized) RDM corresponding to sector j and p j = tr(ρ j ). For non-Abelian theories, the second term can be further simplified. For a given sector j, the Gauss law fixes the sum of the charge and the incoming flux in the last vertex. As a consequence, ρ j has a block diagonal structure ρ j =ρ j ⊗ 1 j , where 1 j is the identity on the subspace corresponding to j for the combined incoming flux plus vertex charge. Specifically for SU(2) the identity for the sector j is (2j +1) dimensional, which finally yields Looking at Eq. (8), we can identify three contributions to the entanglement entropy. The last part, S dist := j p j S(ρ j ), represents the physical entropy which can be distilled from the system by means of local operations and classical communication (LOCC). The first two terms, S class := − j p j log 2 p j and S rep := j p j log 2 (2j + 1), which we, respectively, call the classical and the representation part, have their origin in the Gauss law, implying that the physical subspace is not a direct product of the Hilbert spaces for the links and the sites. They cannot be extracted with LOCC and appear due to the embedding of the physical subspace in the larger space spanned by the basis states discussed in Sec. II supporting a tensor product structure [59][60][61].
With the MPS approach, it is straightforward to access the different contributions to the entropy. Looking at a state in the full basis, we can compute the RDM for the L leftmost sites in the sector j by simply applying a local projector Π (L) j on link L projecting it onto flux j. In our reduced formulation the gauge d.o.f. are integrated out, but the value of j can be readout from the fermionic content, and the RDMρ for the corresponding state in our basis is still block diagonal. Thus, we can write wherep j = tr(ρ j ). Similarly to the full basis, we can obtainρ j by applying the corresponding projectorΠ j . In our formulation the model is effectively Abelianized; hence, compared to the full basis, the last term cannot be further simplified and does not give rise to a representation contribution. Since the quantum circuit which takes a state from our basis to the full one does not change j, the weights for the different sectors of the reduced density matrices are equal in both bases, p j =p j . Moreover, since S(ρ j ) = S(ρ j ), the distillable entropy is also equal in both bases (see Appendix F for the formal argument). Thus, we can directly compute the different contributions S dist , S class , and S rep and therefore the total entropy in the full basis from our formulation. Notice, however, that the calculation in the reduced basis is much more efficient, because of the smaller physical dimensions we need to manipulate.
In (1+1) dimensions, a massive relativistic quantum field theory corresponds to a spin model off criticality in the scaling limit, for which the correlation length in lattice unitsξ is large. For such a system, the entanglement entropy for the RDM describing half of the system is given by S ∝ (c/6) log 2 (ξ) [76], where the parameter c is the central charge of the underlying conformal field theory describing the system at the critical point. Taking the continuum limit of the lattice formulation, ag = 1/ √ x → 0, corresponds to approaching the limit of diverging correlation length in lattice units [77]. Consequently, for the full theory without truncation, we expect the entropy for the RDM for half of the system to be logarithmically UV divergent as where c 2 , c 3 are constants and we take into account finite lattice corrections as in Ref. [32]. With our numerical data, we can check if the entropy in the ground state for our family of truncated models diverges, too. To this end, we look at the different contributions to the entanglement entropy for a cut along the center of the system in the same range of values for (j max , m/g, x, N, D) as in the previous sections and study the scaling of S for ag = 1/ √ x → 0. In general, we observe that none of the different contributions to the entropy show strong finite-size effects for bipartitions that are far away from the boundaries [see Fig. 8(a) for an example]. Nevertheless, we may expect an oscillating contribution to the entropy that becomes smaller as the system size increases [78,79]. To minimize these effects, we average over the values obtained for 4 bipartitions around the center to estimate the different entropy contributions for the half-chain. As Figs. 8(b) and 8(c) indicate, these averaged values are essentially converged in bond dimension and their dependence on the system size is negligible. Hence, we simply take the values obtained for D = 200 as the central value for every combination of (j max , m/g, x, N ) and estimate our error as the difference with respect to the value obtained for D = 150. Additionally, we take into account a systematic error due to the finite precision in our simulations (see Appendix E for details). To compensate for residual finite-size effects, we take the weighted average for every (j max , m/g, x) for the largest two system sizes available. In a final step, we extrapolate the total entropies, obtained from the sum of the different contributions, to the limit ag = 1/ √ x → 0 by fitting our data to Eq. (9). indicating that the logarithmic term contributes and the entropy is indeed UV divergent. The final results for c for different truncations as a function of the bare fermion mass are shown in Fig. 9. Notice that for the full the- 1 in all our fits), and the values for c are close to the one for the full theory (see Fig.  9). Again, there is hardly any difference between the data for j max = 3/2 and j max = 2. For the smallest two truncations j max = 1/2, 1, instead, the picture is significantly different. In these cases our data are not very well compatible with a logarithmic divergence for m/g ≥ 0.2 (resulting in χ 2 d.o.f. 1 in our fits). In the region of small m/g, for which our data follow Eq. (9) reasonably well, the central charges we obtain differ noticeably from 2. Hence, although the ground-state energy densities are rather insensitive to the truncation, the entanglement entropies of the same states show another sign that we do not recover the proper continuum limit for a small j max truncation.

V. CONCLUSION
We introduce an efficient color-neutral basis for a (1+1)-dimensional SU(2) lattice gauge theory on a finite lattice with OBC. Building on a color-neutral basis used for the strong coupling expansion of the model [55,64], we show how to remove the gauge d.o.f. Moreover, our formulation allows us to truncate the maximum colorelectric flux at a finite value j max in a gauge-invariant manner, yielding a family of SU(2) gauge models with finite-dimensional Hilbert spaces, that coincide with a SU(2) lattice gauge theory in the limit j max → ∞. While general methods exist to truncate arbitrary gauge models with discrete finite or continuous compact Lie groups to a finite dimension [9], the truncation achieved for this particular case is a lot more efficient.
The basis we develop is completely general and can in principle be used with any numerical technique. Here, we combine the use of MPS with an efficient truncation for the color-electric flux to explore different limits. Because of the reduced number of d.o.f., we are able to reach values for j max far beyond the ones reached in previous numerical work with TNs [25,28]. To systematically study truncation effects, we compute the ground-state energy density, the entanglement entropy in the ground state, the vector mass gap and its critical exponent for a family of truncated SU(2) models with a maximum color-electric flux of j max = 1/2, 1, 3/2, 2.
We observe that the continuum estimates for the ground-state energy density are rather insensitive to the truncation. Even for the simplest nontrivial truncation, the deviation between the values obtained from our family of truncated model with respect to the continuum result of the full theory is only at the percentage level. Moreover, the results converge quickly with increasing j max such that between results for j max = 3/2 and 2 we observe hardly any difference.
In contrast, the vector mass gap is a lot more sensitive to a truncation of the maximum color-electric flux. For the simplest nontrivial truncation, j max = 1/2, we cannot control lattice effects in the extrapolations well and reliably estimate the errors. The final value obtained for the mass gap in this case differs significantly from previous numerical results. For j max = 1, lattice effects are becoming smaller, thus allowing for a reliable error estimate. Nevertheless, they are still pronounced, and again the continuum estimate for the vector mass gap is not compliant with previous numerical results within error bars. On the contrary, for j max = 3/2, 2 the continuum extrapolations are unproblematic and we obtain precise values for the vector mass gap which agree with the ones from Ref. [55]. Although our data for j max = 1/2, 1 do not allow us to rule out with certainty that for finer lattices the results would approach the continuum result of the full model, the pronounced lattice effects in those cases might indicate that the continuum limit for these truncated models does not exist, as it is the case for quantum link models [54]. Our findings for the critical exponents for the vector mass gap are essentially in agreement with a calculation in the large N c limit [75].
Looking at the scaling of the bipartite entanglement entropy in the ground state towards the continuum limit, we observe similar effects as for the vector mass gap. The central charges for the two simplest nontrivial truncations differ noticeably from the expected value of 2 for two Dirac fermions and, in particular, for large bare fermion masses our data do not follow the expected logarithmic UV divergence well. On the contrary, for j max = 3/2, 2 our numerical values show a clear indica-tion of a logarithmic divergence and we find values close to 2 throughout our entire regime of bare fermion masses we study. Thus, although the ground-state energy densities extracted for j max = 1/2, 1 are close to the values for the continuum model, this is giving a further indication that for these truncations we do not recover the continuum theory in the limit of vanishing lattice spacing.
In general, our findings for the ground state in the SU(2) case are consistent with those recently reported for the Schwinger model with truncated gauge links [80]. There it was also observed that truncating the maximum electric field to a modest value yields a ground state close to the one of the full model in a wide range of bare fermion masses and lattice spacings.
In our calculations we target only the vector state besides the ground state. Other masses in the theory, such as scalar mass gap or baryon masses, can be computed in a similar fashion with the basis we develop. Moreover, our formulation is not restricted to static problems and can be used to compute time evolution, thus also giving access to dynamical properties.
Additionally, the formulation is also potentially suitable for designing future quantum simulators [41,42,81]. As the number of basis states is drastically reduced with respect to the full basis, and a truncation at a maximum value of j max is straightforward, this could allow for simpler experimental realizations compared to previous proposals. Our results also show that in such a simulator one would be able to obtain good estimates for the ground-state energy for the full theory in the continuum, even with the simplest nontrivial truncation for the color-electric flux. However, other quantities as, for example, low-lying excitations or the scaling of entanglement entropy, seem to be more delicate and we only recover the values for the full theory in the continuum, if j max is chosen large enough.
In this work, we focus on the (1+1)-dimensional case, for which it is possible to remove the gauge d.o.f. completely due the the absence of transversal directions. Although in higher dimensions a complete elimination is not possible, it is feasible to construct formulations that reduce the number of d.o.f. as much as possible [82]. In principle, the physical basis we present here could be modified to realize such a maximal gauge fixing. This would be useful for both the numerical simulations with TN in higher dimensions and the potential quantum simulation of the models, since in both cases the (computational or physical) resources needed would be largely reduced.
Appendix A: Dimension of the physical subspace Here, we compute the dimension of the physical subspace contained in the basis developed in the main text. For all the following, we focus on a system with N sites and the case of vanishing background field, j 0 = 0.
As shown in the main text, an arbitrary basis state can be expressed as |φ = ⊗ N k=1 |α k , |α k ∈ {|0 , |1 − , |1 + , |2 }. To calculate the dimension of the physical subspace, it is convenient to represent the 4 N basis states as directed paths from the root r to one of the leaves in a perfect quaternary tree of depth N . The vertices at level k are labeled with the color-electric flux j k at link k, implied by the fermionic states sitting at the edges along the path from the root to the vertex due to Eq. (6) (cf. Fig. 10). Unphysical states now correspond to directed paths from the root to one of the leaves that contain at least one vertex labeled by a negative number. Looking at a path starting from the root along vertices with non-negative labels to the vertex v at level k − 1, a vertex v with negative label in level k can arise if and only if (i) v is labeled by a 0 and (ii) from v the path is continued along the edge corresponding to |1 − , thus ending up in a vertex v labeled by −1/2. Hence, the first vertex v with negative label along a path corresponding to an unphysical state is always one carrying a −1/2. Moreover, all paths containing v necessarily correspond to unphysical states. The number of paths containing v is simply the number of paths through the perfect quaternary subtree of depths N − k rooted by v , 4 N −k . Consequently, the number of unphysical states is given by Here,t k is the number of vertices carrying label −1/2 at level k, for which the path starting from r does not yet pass any other vertex with a negative label. Because of observations (i) and (ii)t k is equivalent to the number of paths t k−1 starting from r to a vertex at level k − 1 with label 0 that did not yet pass any vertex with a negative label. As we show in the following, t k−1 is exactly given by the Catalan number C k = (2k)!/(k + 1)!k!. As a result, the number of unphysical states is To compute the number of paths starting from the root to a vertex at level k with label 0 that did not yet pass any vertex with negative label, we use to following observations.
1. As explained in the main text, the number of edges l k,− passed with |1 − at any level k < k must not exceed the ones with |1 + , l k,+ , to avoid encountering any negative vertices along the path.
2. Looking at a path from the root with j 0 = 0, to any vertex labeled by 0 at level k, we immediately see that the condition l k,+ = l k,− has to be fulfilled in order to compensate for the flux changes induced by |1 − and |1 + . In particular, this implies that 2l k,+ ≤ k or, equivalently, l k,+ ≤ k/2 . The other k − 2l k,+ edges along the path have to carry |0 or |2 , as those states do not lead to a flux change while going from one layer to the other.
The number of paths of length 2l k,+ which contain at any point at least as many |1 − as |1 + is exactly the number of Dyck paths and given by the Catalan number C l k,+ [83]. Hence, the number of paths fulfilling conditions 1 and 2 at level k is given by where the factor k 2l k,+ takes into account the number of ways that the 2l k,+ symbols |1 + and |1 − can be distributed among the k levels and 2 k−2l k,+ the possible ways of filling the remaining edges with |0 and |2 . In the last step, we use an identity for the Catalan numbers. Thus, the dimension of the physical subspace is given by In Fig. 11, we show a comparison between the scaling of our basis and the basis from Ref. [55], which has dimension d N,Hamer = 3 N (2j max + 1) N −1 . As the figure reveals, for systems with OBC our basis offers a vast improvement over the color-neutral basis from Ref. [55] already in the case of the simplest nontrivial truncation j max = 1/2. Even though the fraction of physical states in our basis, d N,phys /4 N , quickly decreases with system size, the total number of states is still exponentially smaller.

Appendix B: Hamiltonian
In this appendix we show how the terms of the Hamiltonian given in Eqs. (3) and (4) can be formulated in the basis presented in Sec. III. As explained in the main text, the mass term is straightforward as even in the original formulation it depends only on the fermionic occupation number. The color-electric energy term can also readily be formulated in the new basis using Eq. (6). Hence, W 0 is given by where the single-site operators arê From Eq. (B1) one can see explicitly that integrating out the gauge field leads to nonlocal interactions in the color-electric energy term. The hopping term can be obtained by translating the possible hopping processes shown in Fig. 1(c) in the new basis. The possible hopping processes in the new basis are listed in Table II. As the table reveals, the matrix elements for certain transitions depend on the color-electric flux, thus also leading to long-range interactions in the hopping term. The hopping term V can then be expressed in the new basis by defining the operators O i,k and translating the possible transitions from Table II in  Fig. 1(c) expressed in the new basis. The value on the right-hand side shows the corresponding matrix elements for the hopping operator.
operator form: where H.c. refers to the Hermitian conjugates of all terms appearing in the formula above. The color-fluxdependent constants c ± k are given by and are nothing but the matrix elements shown in Table  II. In order to compute these constants, the value for j k has to be reconstructed from the fermionic occupation number via Eq. (6).

Appendix C: Penalty terms
As mentioned in Sec. II and further discussed in Appedix A, the basis would in principle allow for unphysical states implying negative values for j. For our numerical calculations with MPS, we choose to remove those states with appropriate penalty terms shifting unphysical states high enough in the spectrum such that they do not interfere with the low-lying spectrum we are targeting. More specifically, we add the following penalty to the Hamiltonian terms from Eqs. (B1) and (B4) which has a nonvanishing contribution at sites with negative value of j k : where Θ(x) is the Heaviside step function. The constant λ 1 has to be chosen large enough to shift states with a negative value for j high enough in the spectrum such that they do not mix into the low-lying spectrum we are interested in.
Moreover, in our calculations we focus on the vector meson states. To avoid any baryon states, we restrict the total fermion number to the number of sites in the system. This can be easily achieved by adding another penalty term, to the Hamiltonian, where λ 2 again has to be chosen large enough to shift the states high enough in the spectrum to prevent them from mixing with the low-lying ones we are interested in.
In addition to that, we are interested in the subspace with j 0 = 0 = j N . Another energy penalty of the form makes it possible to restrict the calculations to that sector, where λ 3 is again a constant that has to be chosen large enough to penalize unwanted states sufficiently.
For our calculations presented in the main text, we have checked the expectation values for all three penalties and found that they are negligible for λ i = 1000, i = 1, 2, 3.

Appendix D: Distinguishing vector and scalar states
Because of the fact that we are working with finite lattices with OBC, the symmetries which allow for distinguishing between the different meson states are no longer preserved. Nevertheless, following the ideas from Ref. [15], the remnants make it possible to separate the different type of states. However, in the basis formulation presented in the main text, it is not straightforward to write down a pseudomomentum operator as was done in Ref. [15]. Thus, to identify the zero momentum excitations of the ground state, we use a simpler approach. On a lattice with periodic boundary conditions, the zero momentum states correspond exactly to translational invariant states. For our finite lattice this should still be approximately fulfilled as long as the system size is large enough. Because of the staggered formulation a translational-invariant state should be invariant under a cyclic shift by two lattice sites. To assign a pseudomomentum to our states, we compute the expectation value of the operator C (2) , where C (k) describes a cyclic shift by k lattice sites to the right. Moreover, to probe for the charge conjugation number, we proceed again similar to Ref. [15], and apply a cyclic shift followed by exchanging the two states |0 ↔ |2 , N k=1 (|0 2| k + |2 0| k ) C (1) . While this lattice analog of charge conjugation is not a good quantum number in the case of OBC, the phase of this operator allows for distinguishing between vector candidates (charge conjugation number −1) and scalar candidates (charge conjugation number +1). For states with charge conjugation number +1, we observe phases close to 0, whereas for states with charge conjugation number −1, the observed phase is close to π. Together with the dispersion obtained from the pseudomomentum operator this allows us to identify the different states as shown in Fig. 12. Above the ground state we observe a vector candidate with C (2) ≈ 1. Subsequently we discover the momentum excitations of the vector state which are characterized by decreasing C (2) before we finally obtain a scalar candidate with C (2) ≈ 1 again.

Appendix E: Extrapolation procedure and error estimation
Here, we give more details on how we extrapolate our data in bond dimension, system size, and lattice spacing to reach the continuum limit, as well as for the extraction of the critical exponents presented in Table I and the  central charges. In a first step, we estimate the exact values for the ground-state energy and the vector mass gap by extrapolating for each combination of (j max , m/g, x, N ) to the limit D → ∞. As the examples in Figs. 3(a) and 3(b) show, we plot our data versus 1/D and fit a line through the two data points with the largest bond dimension. As an estimate for the exact energy value (vector mass gap), we take the mean value between our value computed with maximum bond dimension, E 0,Dmax (N, x), [∆ vec,Dmax (N, x)], and the extrapolated value with infinite bond dimension, E 0,D∞ (N, x) [∆ vec,D∞ (N, x)]. We estimate the error of the central value in a standard manner by taking half of the difference between these two values δE fit = (E 0,Dmax (N, x) − E 0,D∞ (N, x))/2. Additionally to the extrapolation error, the central value obtained also has a systematic error due to the convergence tolerance of ε = 10 −6 set in the simulations, δE syst = εE 0,D∞ (N, x). Taking into account this error, too, the final error for the ground-state energy (vector mass gap) is given by δE tot = (δE fit ) 2 + (δE syst ) 2 .
With the estimates for the exact ground-state energies and the vector mass gaps, we extrapolate to the thermodynamic limit for each combination of (j max , m/g, x), where we use the asymptotic behavior up to linear order [3], and propagate our errors from the previous extrapolation in D. As an estimate for the error we take the fitting error with a 1σ confidence interval. In a final step, we extrapolate to the continuum limit ag = 1/ √ x → 0 by fitting a polynomial in 1/ √ x up to second order. In general, we take the value obtained by the lowest order statistically significant fit, which achieves χ 2 d.o.f. < 1 as the central value. Statistically significant means that the errors for the fit coefficients are smaller than the actual value of the coefficient. In addition to the statistical error of the fit, we estimate the systematic error as the difference between our central value and the next statistically significant fit of the next highest order and/or omitting the largest lattice spacing. For the ground-state energy density this is for all cases a quadratic fit. As an estimate for the systematic error, we take the difference with respect to the value obtained by a quadratic fit omitting the largest lattice spacing, meaning in the region x ∈ [70,150].
For the vector mass gap we observe largely enhanced lattice effects for small values of j max . In particular, for j max = 1/2 quadratic fits have high values for χ 2 d.o.f. . However, as we have only five different lattice spacings, we cannot take higher-order corrections into account. Thus, we determine the central value with a quadratic fit taking into account all lattice spacings and again estimate the systematic error as the difference with respect to the value obtained by a quadratic fit omitting the largest lattice spacing. Consequently, in this case the error might be heavily underestimated as we are neglecting higher order corrections. For j max = 1, we find that for all m/g ≥ 0.3 quadratic corrections are sufficient and proceed the same way for estimating the central value and its systematic error as for j max = 1/2. For m/g = 0.25, we estimate the central value via a linear fit taking into account lattice spacings corresponding to x ∈ [90, 150]. The systematic error in this case is estimated as the difference with respect to a quadratic fit in the region [70,150]. In the region of smaller masses m/g ≤ 0.2, we find that both linear and quadratic fits are statistically significant. Hence, we estimate our central value with a linear fit through all available lattice spacings and the systematic error as the difference with respect to a quadratic fit in the same region. For the largest two truncations, j max = 3/2, 2, the quadratic correction loses significance and thus we estimate our central value in those cases with a linear fit. The systematic error is then determined as the difference with respect to a linear fit discarding the largest lattice spacing, corresponding to x = 50. The final results obtained for the ground-state energy densities and the vector mass gaps following the procedure described above are listed in Table III.
The values for the critical exponents for the vector mass gaps are estimated in a similar fashion. Our data for j max = 3/2, 2, which are close to the ones from Ref. [55], reveal that for our largest fermion mass m/g = 1.6 we are already relatively close to the nonrelativistic limit. Hence, we restrict ourselves to data for small fermion mass to estimate the critical exponent. For each value of j max , we fit our data to a power law, γ(m/g) ν , for every interval [0.1, (m/g) max ] with 0.25 ≤ (m/g) max ≤ 0.4. As the central value, we take the fit with the smallest χ 2 d.o.f. . To estimate our systematic error we take the difference between our central value and the fit giving the most outlying value. The statistical error is again given by 1σ error bar for the fitting error.
To study the scaling of the entropy towards the continuum limit and obtain an estimate for the central charges, we proceed as described in the main text. First, for D = 200 and every combination of (j max , m/g, x, N ), we average over the values obtained for 4 bipartitions close to the center for each of the different entropy contributions. To estimate our systematic error, we take the difference with respect to the values obtained with D = 150. Additionally, our data have another systematic error due to the finite precision in our simulations, which has to be added on top. For the entropies we cannot give the same precise estimates for this systematic error as for the energies. To get, nevertheless, a rough idea of the order of magnitude, we compare results with convergence tolerance ε = 10 −6 and ε = 10 −10 . Figure  13 reveals that even for the largest value of j max , where we expect the largest differences between these results, it is around 10 −5 . Hence, we simply assume a systematic error of 10 −5 due to the finite precision of our simulations in all cases. Appendix F: Entanglement entropy in the full basis and our formulation In this appendix, we discuss the relation between the different contributions to the von Neumann entropy in our reduced basis formulation and the full basis. In a first step, we show that one can recover the full basis state from the reduced one by means of an isometry which can be written as a sequence of N local isometries, and thus corresponds to a quantum circuit of depth equal to the system size, N . Afterwards, we formally argue why the weights of different sectors, p j , are the same in both bases, and show the explicit relation between the entanglement entropies computed in each basis.

Mapping to the full basis
Here we show how the full basis state can be recovered from our reduced basis formulation by sequentially applying a local, isometric map. The map is given by where C j1 j2 J 1 2 M = J, M |j 1 , 1 ; j 2 , 2 are the usual Clebsch-Gordan coefficients for coupling two angular momenta j 1 , j 2 with z components 1 , 2 to a total angular momentum J with z component M . The symbol α ∈ {0, 1 − , 1 + , 2} labels the decorated fermionic occupation. The state |α is an eigenstate of the operatorŝ n from Eq. (B2) andq from Eq. (B3), with respective eigenvalues n α and q α . The states |n α s correspond to a relabeling of the full basis |n 1 , n 2 . Different from the main text, we label them with the total occupation number and the z component of the related angular momentum, n α = n 1 + n 2 , s = (n 1 − n 2 )/2, to make the dependence explicit. The effect of the map is to introduce extra d.o.f., for the incoming link and j, on the outgoing link, by simultaneously respecting the proper SU(2) composition rules, which is ensured by the Clebsch-Gordan coefficients [see Fig. 14(a)]. Note, that for empty or doubly occupied sites, q α = 0 and the Clebsch-Gordan coefficients are trivial. In the case of singly occupied sites, |q α | = 1/2 and the spin 1/2 of the single fermion couples to the angular momentum of the previous link to ensure a color-neutral superposition. The prefactors 1/ 2(j + q α ) + 1 ensure proper normalization of the resulting state and have to be chosen such that M † loc M loc is the projector on the physical subspace, as we show in the following. A straightforward calculation yields The sum inside the bracket over the squares of the Clebsch-Gordan coefficients can be simplified as follows:  Fig. 14(b), the sequential application of the map corresponds to a quantum circuit of depth N . In all our calculations, we work in the sector j 0 = 0; hence the left electric field necessarily has to vanish and the input left link is |00 (thus, not explicitly shown in the figure).

Classical part of the entropy
Here, we show that the weights of sectors with a particular value of j on a certain link p j are identical in both the full basis and our formulation. For all the following, we assume a system of N sites in a physical state |Ψ in the reduced basis, corresponding, in the full one, to |Ψ = M|Ψ . We consider the bipartition of the system obtained by cutting at the Lth gauge link.
In the full basis, the RDM for the leftmost L sites has block diagonal structure thanks to the gauge constraints, and we can write ρ j = Π j tr L+1,...,N (|Ψ Ψ|)Π j = tr L+1,...,N (Π j |Ψ Ψ|Π j ) = tr L+1,...,N Π j M|Ψ Ψ |M † Π j , where Π j is the projector on total flux j for the Lth link. In the full basis, this projector acts locally on the link and thus can be written Π j = 1 in ⊗ Π (L) j ⊗ 1, where the left factor is the identity on the inner part, i.e., the part where ρ is defined [see Fig. 15(a)].
The corresponding projection in our basis formulation is given bỹ where qα 1 +···+qα L =j takes into account all basis states, for which the sum of the eigenvalues q α k , k = 1, . . . , L, for the single-site operators from Eq. (B3) is equal to j. The corresponding RDM thus reads ρ j =Π j tr L+1,...,N (|Ψ Ψ |)Π j = tr L+1,...,N (Π j |Ψ Ψ |Π j ) It turns out, as we show next, that the action of the projector Π j on a certain value of the flux link commutes with the isometry that changes the basis, namely, Π j M|Ψ = MΠ j |Ψ . This implies the following: To prove the statement Π j M|Ψ = MΠ j |Ψ , we proceed as sketched in Fig. 15. The individual steps are justified as follows.  On the left-hand side, the only effect of the projector for a sectorj which acts on M loc ⊗ 1 (L) α is to restrict the sum over j in Eq. (F1) to values such that j + q α + q α =j, which is precisely the effect of the projector on the right-hand side.
• (d)=(e): Iterating the step above, we can pull the projector through very M loc block, until the edge of the chain, where the input j is fixed to 0 and can be ignored.

Distillable part of the entropy
We assume that for a physical state, in the reduced basis, the unnormalized RDM for the L leftmost sites for the sector with (outgoing) flux j has spectrum λ σ , i.e., ρ j = σ λ σ |σ σ|,p j = tr(ρ j ) = σ λ σ , with the eigenvectors |σ supported on the reduced basis for the L leftmost sites. We can use the relations discussed above to relate the density operators in the full and reduced basis. The RDM in the full basis, given by (F3), can be computed from the RDM for L + 1 sites in the reduced basis, as sketched in Fig. 16(a), since the isometries acting to the right of the projector cancel out in the trace. Thus, the RDM in the full basis is obtained by first applying the isometry M L+1 = M Thus ρ j consists of 2j + 1 blocks, with identical spectrum {λ σ /(2j + 1)}. Consequently, the entropy is given by S(ρ j ) = − j r=−j σ λ σ 2j + 1 log 2 λ σ 2j + 1 = − σ λ σ log 2 λ σ 2j + 1 =p j log 2 (2j + 1) − σ λ σ log 2 λ σ =p j log 2 (2j + 1) + S(ρ j ).
Putting this result together with the fact that each j sector has the same weight in the reduced and the full representations, p j =p j , we find that the relation between the entropies is S(ρ) = − j p j log 2 (p j ) + j p j S(ρ j ) = − j p j log 2 (p j ) + j p j (log 2 (2j + 1) + S(ρ j )) = j p j log 2 (2j + 1) + S(ρ). (F5)