Search for Efficient Formulations for Hamiltonian Simulation of non-Abelian Lattice Gauge Theories

Hamiltonian formulation of lattice gauge theories (LGTs) is the most natural framework for the purpose of quantum simulation, an area of research that is growing with advances in quantum-computing algorithms and hardware. It, therefore, remains an important task to identify the most accurate, while computationally economic, Hamiltonian formulation(s) in such theories, considering the necessary truncation imposed on the Hilbert space of gauge bosons with any finite computing resources. This paper is a first step toward addressing this question in the case of non-Abelian LGTs, which further require the imposition of non-Abelian Gauss's laws on the Hilbert space, introducing additional computational complexity. Focusing on the case of SU(2) LGT in 1+1 D coupled to matter, a number of different formulations of the original Kogut-Susskind framework are analyzed with regard to the dependence of the dimension of the physical Hilbert space on boundary conditions, system's size, and the cutoff on the excitations of gauge bosons. The impact of such dependencies on the accuracy of the spectrum and dynamics is examined, and the (classical) computational-resource requirements given these considerations are studied. Besides the well-known angular-momentum formulation of the theory, the cases of purely fermionic and purely bosonic formulations (with open boundary conditions), and the Loop-String-Hadron formulation are analyzed, along with a brief discussion of a Quantum Link Model of the same theory. Clear advantages are found in working with the Loop-String-Hadron framework which implements non-Abelian Gauss's laws a priori using a complete set of gauge-invariant operators. Although small lattices are studied in the numerical analysis of this work, and only the simplest algorithms are considered, a range of conclusions will be applicable to larger systems and potentially to higher dimensions.

2 such as the quantum theory of strong interactions or quantum chromodynamics (QCD), are notoriously hard to simulate, and often demand applications of nonperturbative numerical strategies. Lattice-based methods, in which quantum fields are placed on a finite spacetime lattice, provide a natural regulation of ultraviolet modes and along with Renormalization Group (RG) methods, recover the continuum limit of the theory. 1 Such a discretized theory also provides the framework for nonperturbative numerical studies, such as those based in path integral (Lagrangian) and Hamiltonian (canonical) formulations. While these two approaches are intrinsically equivalent, symmetries and constraints are manifested differently in each case [5]. Path integral vs. Hamiltonian formulation of lattice gauge theories. From a computational standpoint, the path-integral formulation has emerged as the primary tool in the LGT program given its parallels with the quantum statistical physics [6,7], and its reliance on efficient state-of-the-art quantum Monte Carlo sampling techniques given this connection [8]. However, in order for such an analogy with statistical mechanics to be established, a Wick rotation to Euclidean spacetime is performed such that only imaginary-time correlation functions are directly accessed in this numerical program. This feature introduces two limitations: firstly the connection to real-time quantities is lost, and except in limited cases (e.g., Refs. [9][10][11][12][13]), practical proposals are lacking for mapping a generic Euclidean correlation function that is obtained numerically on a spacetime lattice to dynamical amplitudes as measured in experiments. Second, a non-zero fermionic chemical potential introduces a sign problem by making the sampling weight in the Euclidean path integral imaginary [14], along with an inherently related signal-to-noise problem observed in correlation functions at zero chemical potential but with non-zero baryonic number [15]. On the other hand, the Hamiltonian formulation lacks the manifestly Lorentz covariance of the path-integral formulation, and further requires gauge fixing. In particular, in the most common gauge in which the temporal component of the gauge field is set to zero, the information about 'constants of motion' is lost and the related constraint must be imposed a posteriori on the Hilbert space. From a computational perspective, a Hamiltonian formulation enables both real-time and imaginary-time simulations. However, the most efficient numerical approach in Hamiltonian-based studies is no longer stochastic as in the path-integral formulation, and the cost of a typical numerical simulation scales with powers of the dimension of the Hilbert space (which itself grows exponentially with the size of the system). Extremely efficient Hamiltonian-simulation algorithms have been developed and implemented in recent years for lowdimensional LGTs using tensor-network methods [16], 1 In asymptotically free theories like QCD. but they rely on strict assumptions on the rate of entanglement growth in the physical system [17], assumptions that generally break down as system evolves arbitrarily in time.
Despite the drawbacks encountered in a Hamiltonian approach to LGTs, there has been revived interest in the Hamiltonian-simulation program given the improving prospects of quantum simulation and quantum computation. A plethora of ideas, proposals, and implementations have emerged for simulations of quantum many-body systems in general [18][19][20][21][22], and quantum field theories and LGTs in particular (see e.g., Refs. ), in recent years, in light of advances in existing and upcoming digital and analog quantum-simulation technologies [63][64][65][66][67][68][69][70][71][72][73][74]. Generally speaking, a quantumsimulating/computing hardware will reduce the exponential cost of encoding the Hilbert space of a LGT onto the classical bits down to a polynomial cost, by storing information onto the quantum-mechanical wavefunctions of qubits. Nonetheless, quantum hardware will continue to exhibit small capacity and noise-limited capability for the foreseeable future. As a result, the search for an ultimate efficient formulation of LGTs for the quantum-simulation program is a crucial first step toward harnessing the power of quantum-simulating/computing platforms. In the meantime, as the classical Hamiltonian-simulation algorithms advance, such efficient formulation(s) can facilitate classical studies as well. In what follows, we elaborate on the meaning of an efficient formulation 2 in the case of non-Abelian LGTs, and will make such efficiency considerations explicit by analyzing in depth the case of the SU (2) LGT in 1+1 Dimensions (D) coupled to matter. The focus is exclusively on the cost analysis of exact Hamiltonian-simulation algorithms using classicalcomputing hardware. Nonetheless, this study lays the groundwork for an analysis of resource requirements for simulating the same theory on quantum hardware, to be presented in future work. A summary of the pros and cons of various representations of Kogut-Susskind theory. In 1980s, Kogut and Susskind formulated a lattice Hamiltonian for Yang-Mills gauge-field theories coupled to matter [83] that recovered the continuum limit, and was shown to be equivalent to Wilson's path-integral formulation [6] of the same lattice theory. The Kogut-Susskind (KS) Hamiltonian further made it possible to perform quantum-mechanical perturbation theory around the 3 strong-coupling vacuum, i.e., the ground state of the theory in the limit where the mass term and the electric-field term in the Hamiltonian dominate, see Sec. II. Spectra and dynamics of both U(1) and SU (2) LGTs using the KS Hamiltonian in low dimensions were later studied toward the continuum limit, using high-order strong-coupling expansions, as well as non-perturbative numerical methods, see e.g., Refs. [84][85][86][87][88][89]. Such investigations have continued to date using state-of-the-art algorithms such as those based on Matrix Product States [16,[90][91][92]. Given the requirement of generating and storing the Hamiltonian matrix whose size is determined by the size of the Hilbert space, and given the infinite dimensionality of the Hilbert space of a gauge theory, the truncation on excitations of the gauge degrees of freedom (DOF) is a necessity. 3 Furthermore, only a small but still exponentially large (in the system's size) portion of the space of all possible basis states are those satisfying the gauge constraints, i.e., the local Gauss's laws. As will be demonstrated, exact and inexact Hamiltonian methods are not capable of simulating the full Hilbert space of a gauge theory even for small systems and small cutoff values on the gauge-field excitations. As a result, as a first step in the computation, a mechanism to construct the physical Hilbert space and its corresponding Hamiltonian must be implemented.
The representation chosen for fermionic and gauge DOE in a given LGT dictates the way the Hilbert-space truncation is performed and whether gauge invariance remains intact. It also determines the complexity of the construction of the physical Hilbert space and its associated Hamiltonian matrix for the sake of computation. For non-Abelian LGTs, in particular, the Gauss's law constraints are more complex to impose. In the case of the SU (2) LGT coupled to matter in the original proposal of Kogut and Susskind, the Hilbert space of the gauge DOF are expressed in terms of local angular-momentum basis states that mimic those associated with rigid-body rotations in fixed and body frames [83]. To construct the physical Hilbert space, these on-site basis states are combined with the fermionic DOF expressed in the fundamental representation of SU (2) in such a way that the net angular momentum is zero at each site. Furthermore, an Abelian constraint is satisfied such that the total angular momenta on the left and right of a link connecting two lattice sites are equal. As it will be demonstrated, even in 1+1 D, the computational complexity of imposing such constraints grows quickly with the size of the system, and with the cutoff imposed on the total angular momentum. The physical states, in general, become linear combinations of a large and growing number of terms in the basis chosen, adding to the complexity of Hamiltonian generation.
With open boundary conditions (OBC) in 1+1 D, it is possible to eliminate any dependence on the gauge DOE with the use of a gauge transformation, along with the application of Gauss's laws at the level of the Hamiltonian operator itself, as already proposed as a viable efficient basis for Hamiltonian simulation of the SU (2) LGT in Refs. [90,91]. Such a trick significantly reduces the computational cost of imposing Gauss's laws on the Hilbert space and eliminates the need for a cutoff on the gauge DOF. However, as will be shown, there remain redundancies in this representation compared with the physical Hilbert space of the KS theory in the angular-momentum basis that grows slowly with the system's size. Furthermore, this formulation can not be extended to higher dimensions since there are not sufficient Gauss's laws present to eliminate all gauge DOF. So it will be beneficial to examine a formulation with better generalizability prospects. A recent bosonization proposal [93,94] for LGTs is explored as well, in which a subset of Gauss's laws corresponding to the Cartan sub-algebra of the SU(N) group coupled to one flavor of fundamental fermions are augmented with an additional U(1) Abelian Gauss's law to allow the elimination of fermionic DOF. While this procedure works in any dimension, it trades the finite-dimensional Hilbert space of the fermions with intrinsically infinitedimensional Hilbert space of the gauge DOF, including an additional U(1) field. The Hilbert space of these bosonic fields must be cut off for practical purposes, which could lead to systematic uncertainties in computations with finite resources. Furthermore, constructing the Hamiltonian matrix in the physical Hilbert space remains computationally involved in the bosonic formulation.
The complexity associated with non-Abelian LGTs in the KS theory in its original formulation is the motivation behind the development of a recent framework called the Loop-String-Hadron (LSH) formulation for the SU (2) LGT coupled to matter, which is valid in any dimensions [40,95]. It is founded upon the prepotential formalism of pure LGTs [96][97][98][99][100], which is fundamentally a representation that re-expresses the angularmomentum basis in terms of the harmonic-oscillator basis of Schwinger bosons [101]. 4 As a result, the SU(2) gaugelink and electric-field operators are expressed in terms of harmonic-oscillator creation and annihilation operators and allow gauge-invariant operators to be formed out of gauge and fermionic DOF at each site. These operators, therefore, excite only the states in the physical sector of the Hilbert space, as long as an Abelian Gauss's law is satisfied, which requires the number of oscillators at the left and right of the link to be equal. 4 The LSH formulation constructs a complete set of properly normalized gauge-invariant operators and expresses the Hamiltonian in terms of this complete basis [95]. As will be shown, the Hilbert space of the KS theory in the angular-momentum basis after imposing the Abelian and non-Abelian Gauss's laws, and for a given cutoff on the gauge DOF, is identical to that of the LSH Hamiltonian. Nonetheless, the computational cost of generating the associated Hamiltonian is far less in the LSH framework given its already gauge-invariant physical basis states. Consequently, the simplicity of the fermionic representation (with OBC) is enjoyed by the LSH formulation as well but without associated redundancies, and with the prospects of straightforward applications to higher dimensions. Outline of the paper. While all the different formulations studied here have been introduced, and to some extent implemented, in literature, the conclusions briefly stated above and those that will follow, are new and have resulted from a thorough comparative analysis that is conducted in this work. In particular, an analysis of the size of the full and physical Hilbert spaces as a function of the system's size and, when applicable, the cutoff on the gauge DOF is presented in Sec. III for all the formulations of SU (2) LGT in 1+1 D enumerated above (and reviewed in Sec. II). Here, empirical relations are obtained from a numerical study with small lattice sizes. These results lead to a discussion of the time complexity of exact classical Hamiltonian-simulation algorithms in Sec. V. Section VI contains an analysis of the impact of the cutoff on the spectrum and dynamics of the theory. A detailed discussion of the global symmetries of SU (2) LGT in 1+1 D is presented in Sec. IV, which allows the decomposition of the physical Hilbert space of the theory to even smaller decoupled sectors, hence simplifying the computation. While not a focus of this work, a brief comparative study of the KS SU(2) theory in 1+1 D with a QLM formulation [82] is presented in Appendix A. Given the extent of discussions, and the spread of observations and conclusions made throughout this paper, Sec. VII will summarize the main points of the study more crisply, along with presenting an outlook of this work.
In summary, the results presented here should offer a clear path to the practitioner of Hamiltonian-simulation techniques to evaluate the pros and cons of a given formulation of the SU (2) LGT in 1+1 D in connection to the simulation algorithm used. A similar study for the 2+1-dimensional theory can shed light on the validity of the conclusions made for higher-dimensional cases. 5 Moreover, the conclusions of this work will guide future studies of non-Abelian LGTs in the context of quantum simulation.
Within the Hamiltonian formulation of LGTs introduced by Kogut and Susskind, the temporal direction is continuous while the spatial direction is discretized. Each site along the spatial direction is split into two staggered sites, as shown in Fig. 1, such that matter and anti-matter fields occupy even and odd sites, respectively. The number of sites along this direction is denoted by N and is called the lattice size throughout. The spacing between adjacent sites after staggering is denoted as a. For the SU (2) LGT in 1+1 D, the KS Hamiltonian can be written as: Here, H where N 1 = N − 1 for PBC and N 1 = N − 2 for OBC. The fermion field ψ is in the fundamental representation of SU(2) and consists of two components, i.e., ψ = ψ 1 ψ 2 . The gauge linkÛ (x) is a 2 × 2 unitary matrix operator which emanates from site x along the spatial direction and ends at point x + 1, as shown in Fig. 1 Here, N 2 = N − 1 for PBC, N 2 = N − 2 for OBC, and g is a coupling. Further, Ê L andÊ R are the left and the right electricfield operators, respectively, as shown in Fig. 1. These satisfy the SU(2) Lie algebra at each site, where abc is the Levi-Civita tensor and the spatial dependence of the fields is suppressed in these relations. Further, the electric fields on different sites commute. The electric fields and the gauge link satisfy the canonical commutation relations at each site, where T a = 1 2 τ a , and τ a is the a th Pauli matrix. The corresponding commutation relations for fields with different site indices vanish.
Finally, H Here, N 3 = N −1 for both PBC and OBC, and m denotes the mass of each component of the fermions. In this theory, a fermion SU(2)-charge-density operator defined at each site, which satisfies the SU(2) Lie algebra. It further satisfies the following commutation relation at each site, Such a commutation relation vanishes for fields at different sites. This SU(2)-charge-density operator also commutes with the electric fields and the gauge link. With these commutation relations, and those given in Eqs. (4) and (5), one can show that the Hamiltonian in Eq. (1) commutes with the following operator, As a result, the Hilbert space of the theory is classified into sectors corresponding to each of the eigenvalues of the Gauss's law operators G a , and these eigenvalues are the 'constants of motion'. The physical sector of this Hilbert space is that corresponding to the zero eigenvalue of this operator.
A. Angular-momentum formulation The first step in forming the Hilbert space of a LGT for the sake of computation is to map the vacuum and the excitations of the fields to a state basis. In the absence of the magnetic Hamiltonian, which is the case in 1+1 D LGTs, the most efficient basis is formed out of eigenstates of the electric-field operator. The direct product of the fermionic eigenstates and the electric-field eigenstates forms the full Hilbert space. This is called the electricfield basis, or the strong-coupling basis, i.e., in the g → ∞ limit, the interaction terms in Eq. (2) that involves transitions between different eigenvalues of the electric-field operator becomes insignificant compared with the electricfield term, Eq. (3), and the Hamiltonian becomes diagonal in the electric-field basis. Since the electric fields satisfy the SU(2) algebra, a familiar representation is the angular-momentum representation. In fact, as pointed out by Kogut and Susskind, the left and right electric field can be mapped to the body-frame (Ĵ b ) and spaceframe (Ĵ s ) angular momenta of a rigid body. Explicitly, R on each link. Given this correspondence, one may write the electricfield basis states for the KS formulation as for each site x. Here, f 1 and f 2 quantum numbers refer to the occupation number of the two components of the (anti)matter field, ψ 1 and ψ 2 , each taking values 0 and 1, corresponding to the absence and presence of (anti)matter, respectively: at each site, and similarly for the other component of ψ.
Here, δ denotes the Kronecker-delta symbol. Furthermore, the angular-momentum basis states satisfŷ at each site x, where for brevity the site indices are suppressed. Here, J L , J R = 0, 1 2 , 1, 3 2 , · · · , and m R and m L quantum numbers satisfy −J R ≤ m R ≤ J R and −J L ≤ m L ≤ J L , as dictated by the angular-momentum group algebra. The action of the gauge-link operator on this basis can be written as: where α, β = ± 1 2 . 7 Note that J L and J R on each link are equal, and as such we have defined J ≡ J in this relation. Ellipses denote states that precede and follow those shown at site x and x + 1, respectively.
The physical states |φ (KS) can be formed by identifying proper linear combinations of the basis states in Eq. (10) such that the Gauss's laws are satisfied at each site, and by constructing the direct product of these combinations for adjacent sites along the lattice, following additional gauge and boundary conditions as is detailed below. First, given the Gauss's law operators defined in Eq. (9), the physical states |φ (KS) are required to satisfy G a (x) |φ (KS) = 0. Explicitly, for a = 1, 2, 3, and for every x where x = 0, 1, · · · , N − 1 along the one-dimensional lattice. So the Gauss's laws can be simply interpreted as the angular momentaĴ L (corresponding toÊ L ),Ĵ f withĴ 2 f = 3 4 (corresponding to the presence of one and only one fermion), andĴ R (corresponding to −Ê L ) should add to zero at each site. When there is no fermion present or two fermions are present,Ĵ f = 0 and the left and right angular momenta are the same. Moreover, as mentioned before, J L and J R quantum numbers need to be equal on each link. These two requirements, in addition to the boundary conditions imposed on the J R value at site x = 0 and the J L value at site x = N − 1, constrain the Hilbert space to a physical gauge-invariant one, as analyzed in Sec. III A 7 Note that:

B. Purely fermionic formulation
The KS Hamiltonian in Eq. (1) combined with the Gauss's law constraints on the Hilbert space, in essence, leaves no dynamical gauge DOF in 1+1 D beyond possible boundary modes. In particular, with OBC where the incoming flux of the (right) electric field is set to a fixed value, the value of electric-field excitations throughout the lattice is fixed. This, in fact, is a general feature of LGTs in 1+1 D, as is evident from the proof outlined below. As a result, the KS Hamiltonian acting on the physical Hilbert space can be brought to a purely fermionic form, in which the identification of (anti)fermion configurations is sufficient to construct the Hilbert space. This eliminates the need for adopting a state basis for the gauge DOF, and for solving the complex (non-diagonal) Gauss's laws locally which is the case in an angularmomentum basis. Such an elimination of gauge DOF in LGTs in 1+1 D was first discussed in Ref. [85] and is used in recent tensor-network simulations of the SU (2) LGT in Ref. [91]. Here, we present a generic derivation of such a purely fermionic representation, before analyzing its Hilbert space in the following section.
Consider the following gauge transformation on the fermion fields in the KS Hamiltonian: Note that the products of gauge links are defined as y<x U (y) = U (0)U (1) · · · U (x − 1). Consequently, the gauge links must transform as such that the gauge-matter interaction Hamiltonian, H (KS) I , remains invariant: Now considering the unitarity condition on the gauge links, i.e., U † U = I, reveals that Eq. (21) simplifies to where I is the unity matrix whose dimensionality is equal to that of the fundamental representation of the gauge group, e.g., two in the case of SU (2). Therefore, the interaction Hamiltonian becomes where N 1 = N − 2 as noted after Eq. (6).

7
Now given the relation among the gauge link and the left and right electric fields belonging to the same link [106], one obtains the following relation in the new gauge: This relation, combined with the OBC set to E a R (−1) = a 0 for a = 1, 2, 3, and the Gauss's laws G a |φ (KS) = 0 with G a defined in Eq. (9), fully fixes the values of E L and E R at all sites on the one-dimensional lattice in the physical Hilbert space: with ρ a defined in Eq. (7). Consequently, the electricfield Hamiltonian H where N 2 = N − 2 for OBC as noted after Eq. (3). The consequence of applying Gauss's laws to arrive at Eq. (28) is that the local electric-field Hamiltonian in the original formulation is replaced with arbitrary-range fermionfermion interactions in the fermionic Hamiltonian. Finally, the mass term in the new gauge remains the same, as is expected from gauge invariance: where N 3 = N − 1 as noted after Eq. (6). Note that upon expanding Eq. (28), terms with a fermionic structure similar to the mass term arise, effectively modifying the mass in the new representation. The procedure outlined above establishes that any explicit dependence on the gauge link and electric fields are eliminated in the KS Hamiltonian with OBC, giving rise to a purely fermionic Hamiltonian whose terms are specified in Eqs. (24), (28), and (29), and which is identical to the original KS Hamiltonian only in the physical Hilbert space. As a result, any state in this formulation can be written in terms of a complete fermionic occupationnumber basis, where as before, f 1 and f 2 refer to the occupation number of the two components of the (anti)matter field, ψ 1 and ψ 2 , respectively, each taking values 0 or 1.
Gauge transformation, along with the imposition of the local Gauss's laws with OBC, led to the elimination of the gauge DOF in the previous section. Unfortunately, this procedure can obtain a purely fermionic theory only in 1+1 D, as in higher dimensions the number of constraints at each lattice site is not sufficient to eliminate the gauge DOF in all spatial directions. One could reversely consider eliminating the fermionic DOF with the use of the Gauss's laws, as proposed in Ref. [93], to obtain a fully bosonic theory. This protocol works in all dimensions, but in the case of SU (2N ) theories, requires enlarging the gauge group to U (2N ) to accommodate a sufficient number of constraints needed to eliminate the fermions. 9 One further needs to keep track of the fermionic statistics by encoding in the purely bosonic interactions, the nontrivial signs associated with the anti-commuting nature of the fermions [94]. The extended theory can be shown to be equivalent to the original theory for all physical purposes, as long as the cutoff on the new gauge DOF of the extended symmetry is set sufficiently high, see Sec. III C. In the following, the bosonized form of the SU (2) LGT in 1+1 D is derived, following the procedure outlined in Ref. [93] for general dimensions. Consider the Gauss's laws in the KS formulation of the SU (2) LGT in 1+1 D, given in Eq. (18). Although there exist three Gauss's laws at each site, only the Gauss's law corresponding to the a = 3 component of Gauss's law operator in Eq. (9) provides a diagonal relation in the angular momentum/fermionic basis. In other words, two of the Gauss's laws mix basis states with different quantum numbers, and only one of the Gauss's laws leads to an algebraic relation among the gauge and fermionic DOF. Explicitly, for a basis state However, in order to fully express the {f 1 , f 2 } quantum numbers at each site in terms of the {J R , m R , J L , m L } quantum numbers surrounding the site, at least one more independent relation is needed. Such a relation can be obtained by adding an extra U(1) symmetry to enlarge the gauge group, effectively modifying each link on the lattice by a U(1) link U 0 , i.e., where U is the SU(2) link. This introduces a staggered 8 U(1) charge, 10 along with a corresponding U(1) electric field E(x) on each link emanating from site x.
The physical Hilbert space of the extended U(2) theory is the direct product of the physical Hilbert spaces of the KS SU(2) theory and the U(1) theory, i.e., |φ U (2) = |φ (KS) ⊗ |φ U (1) , where |φ (KS) was introduced in Sec. II A, and withÊ for all values of the electric field that satisfy the U(1) Gauss's lawĜ 0 |φ U (2) = 0 witĥ Explicitly, the U(1) gauss's law acting on the new basis Here, OBC is considered with E(−1) = 0 , where 0 is a constant integer. From Eqs. (32) and (36), the {f 1 , f 2 } quantum numbers at each site become redundant, as they can be written as for i = 1, 2. Here, s 1 = −1, and s 2 = 1. As a consequence of Eq. (37), the action of the mass Hamiltonian ∝ ψ † (x)ψ(x) on |f 1 , f 2 (x) can be written as the action of the corresponding operators in the set ,Ê(x)} on the bosonic basis states in the physical Hilbert space, effectively rendering a purely bosonic term. However, the action of the mattergauge interaction Hamiltonian ∝ ψ † (x)U (x)ψ(x+1)+h.c. will be non-trivial due to the non-commuting nature of the fermions, and this feature must be built in the purely bosonic formulation explicitly. In other words, the gaugematter interaction Hamiltonian must carry the information regarding fermionic signs in its purely bosonic form.
While there are a number of protocols for transforming the fermions to hardcore bosons (spins) such that the 10 The staggered term in the U(1) charge ensures that filled even and odd sites have opposite charges, corresponding to the presence of matter and antimatter at the site, respectively. fermionic anti-commutation relations are preserved, such as the familiar Jordan-Wigner transformation [107], an alternative protocol that preserves the locality of interactions is presented in Ref. [94]. For a generic SU (N ) theory, this protocol amounts to augmenting the theory with an additional Z 2 gauge symmetry, where the new local Gauss's law corresponding to the auxiliary gauge keeps track of the fermionic signs, see also Refs. [108,109]. Back to the case of bosonized SU(2) theory, the extended U(2) theory (necessitated by the need for an extra diagonal Gauss's law) includes the Z 2 symmetry as a subgroup. Therefore, the protocol of Ref. [94] does not require introducing an additional Z 2 gauge group. In other words, the transformation from fermions to hardcore bosons can already proceed by exploiting the local U(1) electric fields introduced above. Nonetheless, since the Hamiltonian terms are either local or nearest neighbor, and that only a 1+1 D theory is considered in the current work, all such transformation, local or non-local variants, are of comparable (low) complexity. Distinctions among different transformations become more relevant in the context of quantum simulation, in which the fermions need to be mapped to qubit DOF. Such considerations will be studied in future work.
Besides their coupling to the fermions in the modified matter-gauge interaction Hamiltonian, which guarantees the U(1) Gauss's law constraint, no further dynamics is introduced for the U(1) gauge DOF. As a result, apart from the issue of the fermionic statistics that needs to be dealt with via a separate transformation as discussed above, the Hamiltonian of the extended U(2) theory is a straightforward extension of the KS Hamiltonian presented in Sec. II A, as shown in Ref. [93]. The extended Hamiltonian involves nearest-link interactions, as a result of the replacement in Eq. (37) and the transformation to hardcore bosons, but is otherwise local. The physical Hilbert spaces of the SU(2) theory and the U(2) theory are isomorphic, meaning that in the limit where the U(1) gauge link approaches unity, the Hamiltonian matrix elements in the original theory is recovered from those in the extended theory. This is established straightforwardly for OBC, while for PBC, the isomorphism holds only in a given topological sector [93]. In the remainder of this paper, we analyze the dimensionality of, and the resource requirement for constructing, the physical Hilbert space of the bosonized theory compared with the KS theory, along with the effect of the U(1) cutoff on the spectrum.

D. Loop-String-Hadron formulation
An alternate reformulation of Kogut-Susskind Hamiltonian formalism in terms of Schwinger bosons, known as the prepotential formalism, has been developed over the past decade [96][97][98][99][100][102][103][104][105]. In a recent work [95], the prepotential formalism of the SU (2) LGT has been made complete to include staggered fermions, explicit Hamiltonian, and the associated Hilbert space. In this sec- tion, the LSH formalism in 1+1 D will be introduced. Later on, we demonstrate the advantage of this formulation compared with the original KS theory within the angular-momentum basis in the physical sector, and with the purely fermionic and bosonic formulations.
Within the prepotential framework, the original canonical conjugate variables of the theory, i.e, electric-field and link operators are replaced by a set of harmonicoscillator doublets, defined at each end of a link as shown in Fig. 2. Both the electric-field as well as the link operators can be re-expressed in terms of Schwinger bosons to satisfy all properties of these variables spelled out in section II A. However, the most important feature of the prepotential formalism is that the link operators U , originally defined over a link connecting neighboring sites (x, x + 1), are now split into a product of two parts: where U L(R) is left (right) part of the link attached to site x (x + 1). As a result of this decomposition, the gauge group is now totally confined to each lattice site, which allows one to define gauge-invariant operators and states locally. For the pure gauge theory, these local gauge-invariant operators and states can be interpreted as local snapshots of Wilson-loop operators of the original gauge theory. One can now construct a local loop Hilbert space by the action of local loop operators on the strong-coupling vacuum defined locally at each site. At this point, it must be emphasized that mapping the local loop picture to the original loop description of the gauge theory requires one extra constraint on each link. This constraint demands that the states must satisfy where N L(R) counts the total number of Schwinger bosons residing at the left (right) end of a link connecting sites x and x + 1. 11 This constraint is a consequence of the relation E 2 L = E 2 R on the link and is equivalent to the constraint J L = J R in the angular-momentum basis. 11 In the notation of Ref. [95], these are indicated as N L(R) .
The inclusion of the staggered fermionic matter in the SU (2) LGT is straightforward, and combines smoothly with the local loop description obtained in the prepotential framework. The reason is that both the prepotential Schwinger bosons and the matter fields associated with a given site transform in the fundamental representation of the local SU(2). One can now combine matter and prepotential to construct local string operators, besides local loop operators. Acting on the strong-coupling vacuum, these build a larger local gauge-invariant Hilbert space, including string and 'hadron' states. This complete description is named the LSH formalism in Ref. [95]. The LSH formalism is briefly described in the following, focusing on necessary steps for working with this formalism in one spatial dimension.
Within the LSH framework, a gauge-invariant and orthonormal basis is chosen, that is defined locally at each site and is characterized by a set of three integers: for all x. The three quantum numbers signify loop, incoming string, and outgoing string at each site, respectively. 12 The allowed values of these integers are given by It is clear from the range of the quantum numbers that n l is bosonic, whereas n i and n o are fermionic in nature.
In terms of the LSH formalism, the operator building the local string Hilbert space consists of SU(2)-invariant The LSH states are illustrated on a 1D lattice with N = 6, where N denotes the number of staggered sites. Here, a gauge-invariant state is characterized by loop quantum numbers, nl, and string quantum numbers: ni for an incoming string and no for an outgoing string. Any state may consist any number of open or closed strings (depending on the boundary conditions as well as the cuto↵ on the boson excitations), and any number of hadrons up to N . A sample of LSH states are shown on this lattice, denoted by a set of strings and/or hadrons. The number of solid thick lines passing through a site denotes the loop quantum number nl. The number of solid disks (red or green) at a site from (to) which a thick line starts (ends) denotes outgoing (incoming) string quantum number no(ni). If two strings start and end at the same site, as shown e.g., in the second panel at x = 3, it is equivalent to a hadron sitting on top of a loop, and is denoted as |nl = 1, ni = 1, no = 1i. A single hadron present at a site is denoted by |nl = 0, ni = 1, no = 1i (e.g., the last panel at x = 0).
Within the LSH framework, a gauge-invariant and orthonormal basis is chosen, that is defined locally at each site and is characterized by a set of three integers: for all x. The three quantum numbers signify loop, incoming string, and outgoing string at each site. The al-lowed values of these integers are given by It is clear from the range of the quantum numbers that n l is bosonic, whereas n i and n o are fermionic in nature.
In terms of the LSH formalism, the operator building the local string Hilbert space consists of SU(2)-invariant | n l = 0,n i = 1,n o = 1 1 | n l = 0,n i = 0,n o = 1 2 | n l = 0,n i = 1,n o = 0 3 | n l = 0,n i = 1,n o = 1 4 | n l = 0,n i = 1,n o = 1 5 | n l = 0,n i = 1, The LSH states are illustrated on a 1D lattice with N = 6, where N denotes the number of staggered sites. Here, a gauge-invariant state is characterized by loop quantum numbers, nl, and string quantum numbers: ni for an incoming string and no for an outgoing string. Any state may consist any number of open or closed strings (depending on the boundary conditions as well as the cuto↵ on the boson excitations), and any number of hadrons up to N . A sample of LSH states are shown on this lattice, denoted by a set of strings and/or hadrons. The number of solid thick lines passing through a site denotes the loop quantum number nl. The number of solid disks (red or green) at a site from (to) which a thick line starts (ends) denotes outgoing (incoming) string quantum number no(ni). If two strings start and end at the same site, as shown e.g., in the second panel at x = 3, it is equivalent to a hadron sitting on top of a loop, and is denoted as |nl = 1, ni = 1, no = 1i. A single hadron present at a site is denoted by |nl = 0, ni = 1, no = 1i (e.g., the last panel at x = 0).
Within the LSH framework, a gauge-invariant and orthonormal basis is chosen, that is defined locally at each site and is characterized by a set of three integers: for all x. The three quantum numbers signify loop, incoming string, and outgoing string at each site. The al-lowed values of these integers are given by It is clear from the range of the quantum numbers that n l is bosonic, whereas n i and n o are fermionic in nature.
In terms of the LSH formalism, the operator building the local string Hilbert space consists of SU(2)-invariant | n l = 0,n i = 0,n o = 1 1 | n l = 1,n i = 0,n o = 1 2 | n l = 2,n i = 1,n o = 1 3 | n l = 1,n i = 1,n o = 0 4 | n l = 1,n i = 0,n o = 0 5 | n l = 0,n i = 1,n o = 0 6 | n l = 0,n i = 1,n o = 1 1 | n l = 0,n i = 0,n o = 1 2 | n l = 0,n i = 1,n o = 0 3 | n l = 0,n i = 1,n o = 1 4 | n l = 0,n i = 1,n o = 1 5 | n l = 0,n i = 1,n o = 1 6 The LSH states are illustrated on a 1D lattice with N = 6, where N denotes the number of staggered sites. Here, a gauge-invariant state is characterized by loop quantum numbers, nl, and string quantum numbers: ni for an incoming string and no for an outgoing string. Any state may consist any number of open or closed strings (depending on the boundary conditions as well as the cuto↵ on the boson excitations), and any number of hadrons up to N . A sample of LSH states are shown on this lattice, denoted by a set of strings and/or hadrons. The number of solid thick lines passing through a site denotes the loop quantum number nl. The number of solid disks (red or green) at a site from (to) which a thick line starts (ends) denotes outgoing (incoming) string quantum number no(ni). If two strings start and end at the same site, as shown e.g., in the second panel at x = 3, it is equivalent to a hadron sitting on top of a loop, and is denoted as |nl = 1, ni = 1, no = 1i. A single hadron present at a site is denoted by |nl = 0, ni = 1, no = 1i (e.g., the last panel at x = 0).
Within the LSH framework, a gauge-invariant and orthonormal basis is chosen, that is defined locally at each site and is characterized by a set of three integers: for all x. The three quantum numbers signify loop, incoming string, and outgoing string at each site. The al-lowed values of these integers are given by It is clear from the range of the quantum numbers that n l is bosonic, whereas n i and n o are fermionic in nature.
In terms of the LSH formalism, the operator building the local string Hilbert space consists of SU(2)-invariant | n l = 0,n i = 0,n o = 1 1 | n l = 1,n i = 0,n o = 1 2 | n l = 2,n i = 1,n o = 1 3 | n l = 1,n i = 1,n o = 0 4 | n l = 1,n i = 0,n o = 0 5 | n l = 0,n i = 1,n o = 0 6 | n l = 0,n i = 1,n o = 1 1 | n l = 0,n i = 0,n o = 1 2 | n l = 0,n i = 1,n o = 0 3 | n l = 0,n i = 1,n o = 1 4 | n l = 0,n i = 1,n o = 1 5 | n l = 0,n i = 1,n o = 1 6 The LSH states are illustrated on a 1D lattice with N = 6, where N denotes the number of staggered sites. Here, a gauge-invariant state is characterized by loop quantum numbers, nl, and string quantum numbers: ni for an incoming string and no for an outgoing string. Any state may consist any number of open or closed strings (depending on the boundary conditions as well as the cuto↵ on the boson excitations), and any number of hadrons up to N . A sample of LSH states are shown on this lattice, denoted by a set of strings and/or hadrons. The number of solid thick lines passing through a site denotes the loop quantum number nl. The number of solid disks (red or green) at a site from (to) which a thick line starts (ends) denotes outgoing (incoming) string quantum number no(ni). If two strings start and end at the same site, as shown e.g., in the second panel at x = 3, it is equivalent to a hadron sitting on top of a loop, and is denoted as |nl = 1, ni = 1, no = 1i. A single hadron present at a site is denoted by |nl = 0, ni = 1, no = 1i (e.g., the last panel at x = 0).
Within the LSH framework, a gauge-invariant and orthonormal basis is chosen, that is defined locally at each site and is characterized by a set of three integers: for all x. The three quantum numbers signify loop, incoming string, and outgoing string at each site. The al-lowed values of these integers are given by It is clear from the range of the quantum numbers that n l is bosonic, whereas n i and n o are fermionic in nature.
In terms of the LSH formalism, the operator building the local string Hilbert space consists of SU(2)-invariant bilinears of one bosonic prepotential operator and one fermionic matter field, yielding overall fermionic statistics, whereas the local loop Hilbert space is constructed by the action of SU(2)-invariant bilinears of two bosonic prepotential operators. Such operators will not be introduced in detail, instead the Hamiltonian will be written shortly in this operator basis. Characterization of gaugeinvariant states on a one-dimensional lattice consisting of six staggered sites in terms of the three quantum numbers is illustrated in Fig. 3 via a few examples covering any situation that can occur within this theory.
Let us define a set of LSH operators consisting of diagonal and ladder operators locally at each site as following: 13n 13 In the notation of Ref. [95],n l ,n i , andno are indicated asN l ,N i , andNo, respectively. Further, in that referenceλ ± is indicated Here, the site index x is made implicit for brevity. Being an SU(2) gauge-invariant basis, one no longer needs to satisfy the SU(2) Gauss's laws at each site. However, the neighboring sites still need to be glued together by the Abelian Gauss's law, i.e., Eq. (39). In terms of the LSH operators, the Abelian Gauss's law reads as: Upon acting on the LSH basis states and comparing with Eq. (39), one obtains: (54) where, N L (x) and N R (x) count bosonic occupation numbers at each end of the link connecting site x and x + 1. Pictorially, the left and right sides of Eq. (39) are represented in Fig. 3 by the number of solid lines on the left and right side of each link, respectively. Another important relation iŝ which establishes the relation between the {f 1 , f 2 } quantum numbers of the original formulation and the string quantum number of the LSH formulation. The Hamiltonian of the SU (2) LGT coupled to matter in the LSH formulation can be written in terms of the LSH operators and is given by: Here, H is the electric-energy term, and H (LSH) M is the mass term. Explicitly, in terms of the LSH operators defined in Eqs. (44)-(51), each part of the Hamiltonian can be written as [95]: where (57) contains the LSH ladder operators in the following combinations (suppressing the site indices): andŜ The strong-coupling vacuum of the LSH Hamiltonian is given by It is easy to verify that Eq. (68) satisfies the Abelian Gauss's law, Eq. (52).
This completes the introduction of the LSH formulation for the SU (2) LGT in 1+1 D. In later sections, the finite-dimensional Hilbert space of the theory will be constructed by imposing a cutoff on the N L and N R quantum numbers, and the associated cost of the classical simulation within this framework will be analyzed.

III. PHYSICAL HILBERT-SPACE ANALYSIS
As introduced in the previous section, the naive basis states in the KS LGTs spans a Hilbert space that is predominantly unphysical. The physical sector corresponds to the zero eigenvalue of the Gauss's law operator in Eq. (9). As mentioned before, in contrast to the U (1) LGT, in SU (2) LGT the Gauss's law is not a single algebraic constraint on the eigenvalues of the electricfield operator, but instead, it mixes states with different electric-field quantum numbers, and is therefore a nondiagonal constraint when expressed in the electric-field basis. 14 A major complexity in the Hamiltonian formulation of non-Abelian LGTs is to diagonalize the Gauss's law operator locally to form the physical Hilbert space, as otherwise the computation is prohibitively costly even in small systems.
A question worth addressing is how beneficial it is, from a computational perspective, to work with a formulation that solves the Gauss's law at the level of operators as opposed to states (e.g., the LSH formulation) compared with a formulation that sustains a simple mapping of the Hilbert space to operators in the Hamiltonian but requires solving the Gauss's laws for basis states subsequently (e.g., the KS formulation in the angular momentum representation). Such a cost analysis is presented in Sec. V, but it requires understanding and analyzing in more detail the steps involved in forming the physical Hilbert space in each formulation and the dimensionality of the Hilbert spaces involved. Another interesting question is how fast the dimensionality of the Hilbert space and its physical subsector grows as a function of the lattice size and the cutoff on the electric-field excitations in each of the formulations considered. Such questions are studied in various depth in this section for all the formulations introduced in Sec. II, and briefly in Appendix A for the QLM. 15

A. Gauge-invariant angular-momentum basis
Despite significant state reduction after imposing physical constraints on the full basis states, without a finite cutoff on the electric-field quantum numbers, the physical Hilbert space will still be infinite-dimensional. In the remainder of this paper, we impose: J L , J R ≤ J max , and denote the cutoff as Λ = 2J max . Examining the dependence of the dimension of the Hilbert space, as well as that of observables, on this cutoff is one of the objectives of this work. After the imposition of this cutoff, the states in the physical Hilbert space can be obtained following the procedure outlined in Sec. II A given the boundary conditions specified. For PBC, the J R value at site x = 0 and the J L value at site N − 1 are set equal. For OBC, the J R value at site x = 0 is set to a constant The dimension of the physical Hilbert space, called M throughout, for lattices of the size up to N = 10 and cutoffs up to Λ = 10 is provided in Tables III and IV of Appendix B for PBC and OBC, respectively. There are a few interesting features to observe: For PBC, asymptotically the dimension of the physical Hilbert space grows linearly as a function of the cutoff for all x. This feature is evident from the plot of M as a function of Λ for various N as shown in Fig. 4(a). This is a consequence of the observation that as Λ increases, the number of new allowed states quickly saturates, i.e., introducing an additional possibility for the J R,L quantum numbers amounts to only adding 2N N possible states. For N = 2, 4, 6, 8, 10, the cutoff Λ at which the growth of states become linear afterwards is 0, 1, 2, 3, 4, respectively. The best empirical fits to this linear dependence are shown in the plot. Second, as expected, the dimension of the physical Hilbert space grows exponentially with the system's size at a fixed cutoff, as plotted in Fig. 5(a). The growth, up to constant factors and higher-order terms in the exponent, can be approximated by M ∼ e pN . The coefficient of N in the exponent approaches a constant value as a function of cutoff, as shown in Fig. 6(a). This value can be obtained from a fit to points shown in the plot, as depicted in the figure. For moderate N values such that the higher-order terms in the exponent are negligible, this p value can be used to approximate the dimension of the physical Hilbert space with PBC as Λ → ∞.  For OBC, the dimension of the physical Hilbert space grows as a function of Λ until it becomes a constant for Λ ≥ N (Λ ≥ N + 2 0 for an arbitrary 0 ), as depicted in Fig. 4(b). The reason for this behavior is that the J quantum number only changes (by 1 2 ) from the left to the right side of site x if the site's total fermionic occupation number is equal to one. If the J R value at site x = 0 is set to 0 , it can become at most J L = 0 + N/2 at the last site. Increasing the cutoff beyond this value will not change the states present in the physical Hilbert space. This growth of the dimension of the physical Hilbert space to this saturation value at a fixed N can be approximated by an exponential form, M ∼ e qΛ . The coefficient of Λ in the exponent for various values of N is plotted in Fig. 7, and is seen to asymptote to a constant value at large N . The fit to this asymptotic value is shown in the plot. This value can be used to approximate the dimension of the physical Hilbert space for an arbitrary large N and any Λ. Similarly, the dependence of the dimension of the physical Hilbert space on the lattice size can be approximated by an exponential form, M ∼ e pN , for a fixed cutoff, and up to constant factors and higher order terms in the exponent. The coefficient of N in the exponent asymptotes to a constant value at large Λ, as shown in Fig. 6(b).
The size of the full Hilbert space before implementing physical constraints can be approximated by 16 with PBC, where j = {0, 1 2 , 1, · · · , Λ 2 }. To compare this with the dimension of the physical Hilbert space with PBC, one can again write the lattice-size dependence of the M as e pN . The coefficient of N in this exponent as a function of Λ can be plotted for both the full and physical Hilbert space, as is shown in Fig. 8. As is evident, even for small values of the cutoff, the full Hilbert space grows much faster with the system's size than the physical Hilbert space. For example, with Λ = 5, the p values differ by ≈ 7. This means that for a lattice size N = 10, for example, the dimension of the full Hilbert space is ≈ 30 orders of magnitude larger than that of the physical Hilbert space. As a result, it is not plausible to perform a classical Hamiltonian simulation with a manageable cost if the physical constraints are not imposed a priori. Implementing the physical constraints, nonetheless, introduces further complexity at the onset of the calculation and amounts to an additional preprocessing cost. We will come back to this point when comparing the simulation cost between the KS and LSH formulations in Sec. V.
A physical basis state is generally a superposition of the original angular-momentum basis states, see e.g., the example in Eqs. (70) below. The number of terms in each superposition can become exponentially large in system's size. This creates significant complexity when generating the Hamiltonian matrix, due to the need to keep track of the Hamiltonian action on each constituent basis state. The maximum number of terms in a physical state is plotted in Fig. 9 as a function of Λ for PBC, demonstrating this exponential growth. We will come back to this feature in Sec. V when analyzing the computational cost of the Hamiltonian simulation.

B. Purely fermionic formulation
As discussed in Sec. II B, the basis states that represent the Hilbert space of the purely fermionic representation of the KS formulation with OBC consist of the direct product of on-site fermionic states, see Eq. (30), giving rise to M = 4 N basis states, where N denotes the size of the lattice in 1+1 D as before. The dimension of the Hilbert space of the fermionic theory is larger than the dimension of the physical Hilbert space of the KS formulation in the angular momentum (and LSH) basis for cutoff values that allow the full physical Hilbert space to be constructed with OBC (i.e., Λ ≥ N + 2 0 ). The ratio of the former to the latter is shown in Fig. 10  tion of the lattice size. This form shows that the ratio of the dimensions of the two Hilbert spaces asymptotes slowly to a fixed number.
To understand this mismatch between the number of (physical) states in both formulations, despite the fact the fermionic formulation is constructed to fully represent the physical Hilbert space, inspecting the following example will be illuminating. Consider the N = 2 theory in the ν = 1 sector, where ν denotes the normalized fermion occupation number on the lattice defined in Eq. (76) below. In the KS formulation in the angular-momentum basis, the only four physical basis states are: 17 , , where each triplet in the square brackets denotes at the corresponding site x, and the direct product symbol is suppressed in such triplets for brevity. On the other hand, in the purely fermionic representation of the same theory, the six basis states are , As is seen, while all the six possible fermionic configurations in the ν = 1 sector are present in the physical basis states of the KS formulation in the angular-momentum basis, only two proper linear combinations of states 2)-5) in the fermionic formulation appear in the KS formulation in the angular-momentum basis. Further inspection of the two representations reveals that the spectrum of both theories matches exactly for all values of the couplings, but with degeneracies present in the fermionic case. To conclude, the fermionic representation of the SU (2) LGT in 1+1 D with OBC has redundancies in the representation compared with the KS formulation in the angular-momentum basis, however it avoids complex linear combinations of basis states that arise in the latter due to the imposition of the Gauss's laws. As will be discussed in Sec. III D, the LSH formulation of the SU (2) LGT is free from the redundancies of the fermionic formulation, while at the same time it does not involve a cumbersome physical Hilbert-space construction.

C. Purely bosonic formulation
The physical basis states of the bosonized SU(2) theory with OBC are, at the first sight, the direct product of the physical basis states of the KS theory discussed in Sec. III A and the electric-field basis states satisfying the extra U(1) Gauss's law. Recall that the U(1) symmetry was introduced in the bosonized form to allow the elimination of fermionic DOF in favor of bosonic DOE in the SU(2) theory. The statement above is only true if the cutoff on the U(1) electric field is set sufficiently high such that all fermionic configurations allowed in the physical Hilbert space of the SU(2) theory can be realized. To make this statement more explicit, consider the example studied in Sec. III B, where N = 2 and ν = 1 in the KS theory with OBC, and the incoming fluxes of the SU(2) and U(1) electric fields are set to zero. As was shown, while there are six allowed fermionic configurations in the purely fermionic representation, these reduce to four linear combinations of basis states for Λ ≥ 2 in the physical Hilbert space of the KS formulation in the angular-momentum/fermionic basis (still encompassing all six possible fermionic configurations with ν = 1). Note that with Λ ≥ 2, the physical Hilbert space of the SU(2) theory is complete. Now consider the purely bosonic formulation, with Λ 0 denoting the cutoff on the U(1) electric-field excitations. Obviously for Λ 0 = 0, the only state allowed is the strong-coupling vacuum states, i.e., state 1) in Eq. (70), and the physical Hilbert space of the bosonized U(2) theory has dimension one in the specified sector. For Λ 0 = 1, there are three states contributing, corresponding to states 1), 2), and 3) in Eq. (70), times the U(1) electric-field states (1) for state 1) and (1) for states 2) and 3). Finally, for Λ 0 ≥ 2, all four states in Eq. (70) are allowed, and the corresponding U(1) electric-field states are those given above for states 1), 2), and 3), and |2 (0) ⊗ |0 (1) for state 4). In general, the dimension of the physical Hilbert space of the extended U(2) theory approaches that of the original SU(2) theory with Λ ≥ N , and reaches a saturation value at Λ 0 = N . This trend has been depicted in Fig. 11 for N = 2, 4, 6, 8-site theories. Such an extra cutoff dependence can, in particular, be important in encoding the purely bosonic Hamiltonian onto qubits in a Hamiltonian simulation, in which one trades the few-dimensional Hilbert space of the fermions with the Hilbert space of U(1) bosons that grows with the lattice size when OBC are imposed. Nonetheless, in higher-dimensional gauge theories coupled to fermions, the bosonization may present some benefit by avoiding the non-local fermionic encodings, although the need for a sufficiently large U(1) cutoff may be more significant in higher dimensions. The pros and cons of such a reformulation of the original gauge theory in the context of quantum simulation requires further investigation.

D. Loop-String-Hadron formulation
As described in Sec. II D, the Hilbert space of the LSH formulation is spanned by basis states for x = 0, 1, 2, ..., N − 1, subject to the Abelian Gauss's law constraint N L (x) = N R (x) along each link connecting sites x and x + 1. N L (x) and N R (x) quantum numbers are expressed in terms of the LSH quantum numbers according to Eqs. (53) and (54). In the following, an efficient procedure for generating the physical Hilbert space of the LSH formulation will be presented for both OBC and PBC. It can be shown that this Hilbert space is in one-to-one correspondence with the physical Hilbert space of the KS theory in the angular-momentum representation, and that LSH formulation is a more economical encoding of such a Hilbert space given its reliance on fully gauge-invariant DOF.
OBC fixes the incoming electric flux into the lattice. In the language of LSH quantum numbers, this condition reads: N R (−1) = 0 . Given this, and using Eqs. (39), (53), and (54) consecutively, it is straightforward to show that Eq. (73) implies that for the one-dimensional lattice with OBC, n l quantum number at any site is completely fixed by the boundary condition and string quantum numbers n i and n o at all the sites to its left. In other words, once the string quantum numbers are specified throughout the lattice, all LSH quantum numbers are known, and hence a particular gauge-invariant state is specified. Note that for an N -site lattice, fixing the fermionic quantum numbers give rise to 4 N basis states, that is the same as the number of basis states one obtains with the purely fermionic formulation of the KS theory as described in Sec. III B. However, at this point, it should be noted that there exist certain string configurations (for a fixed value of 0 ) that make the right-hand side of Eq. (73) negative on one or more sites on the lattice. Such spurious string configurations must be discarded, and thus the Hilbert space of the LSH formulation is smaller in size than the purely fermionic formulation. In fact, the dimension of the LSH Hilbert space comes out to be of exactly the same as that of the physical Hilbert space of the KS theory in the angular-momentum basis. However, the cost of generating the physical Hilbert space is much less than that of the KS formulation in the angular momentum basis, as the states no longer need to satisfy SU(2) Gauss's laws at each site (since this has already been taken care of by the LSH construction), and the only re-maining Gauss's law that is Abelian in nature is solved analytically. Moreover, while working with a cutoff Λ, the LSH Hilbert space is constrained to only contain those string configurations that yield 0 ≤ N L/R (x) ≤ Λ for all x. An example of the LSH Hilbert space with OBC on a lattice of size N = 2 and with Λ = 2 and ν = 1 is given in Fig. 12. These basis states are in one-to-one correspondence with the physical basis states in the angular-momentum basis, i.e., the states enumerated in Eq. (70).
yielding (Λ + 1)4 N states to start with. Identifying 0 ≡ N R (−1) and following the same prescription outlined above for OBC, all possible states in the LSH Hilbert space can be constructed subject to cutoff Λ. Note that with PBC, one obtains many copies of the same {n i , n o } configurations corresponding to different winding numbers, i.e., the number of closed loops that go around the lattice. A detailed account of the global symmetries of the KS theory will be presented in the next section.
The LSH Hilbert space, both for OBC and PBC, is identical to the physical Hilbert space of the KS theory, in the sense that each state in the LSH basis corresponds to one and only one state in the KS physical Hilbert space and vice versa. Therefore, all discussions regarding the scaling of the physical Hilbert space presented in Sec. III A are valid for the LSH formulation as well. The only distinction is that the LSH Hilbert space and the associated Hamiltonian can be generated with far less computational complexity, as will be further discussed in Sec. V.

IV. REALIZATION OF GLOBAL SYMMETRIES
The physical Hilbert space, projected out by imposing the Gauss's laws as studied in the previous sections, can be further characterized by global symmetry sectors as well as topological configurations. Identification of these symmetries can further simplify the Hamiltonian simulation, as the Hilbert space that needs to be studied can be further divided to smaller sectors. These symmetries are manifested differently in the case of PBC and OBC, and are discussed in the following.
The total fermion occupation number is conserved, as the operatorQ commutes with the Hamiltonian with both boundary conditions. Note that in the LSH framework, this quantum number is simply Q = x [n i (x) + n o (x)]. With  OBC, Q can be any integer in the interval [0, 2N ]. With PBC, Q can only be an even integer in the same interval, as an odd total fermionic occupation number creates an imbalance between the net flux of electric field into site x = 0 and out of site x = N − 1, which contradicts PBC. For convenience, the global charge associated with the total number of fermions can be normalized by the lattice size as with ν ∈ [0, 2]. The strong-coupling vacuum lies in the ν = 1 sector of the Hilbert space.
Besides the conservation of the total number of fermions, there is an additional quantum number that divides the Hilbert space of each Q sector to multiple disjoint sectors in general. This quantum number in the LSH language can be written as In addition to the total fermionic number and the net flux, the Hamiltonian and the associated Hilbert space are symmetric under charge conjugation. While there is no U (1) charge associated with the fields in the SU (2) LGT, the system is invariant, up to µ → −µ, if the presence of fermions on the lattice is exchanged with their absence, i.,e, corresponding to a particle-hole exchange 1) 2
symmetry. One can verify that the associated charge conjugation operator,Ĉ, commutes with the Hamiltonian of both PBC and OBC, and thatĈ whereQ ≡Q − NÎ, withÎ being an identity operator. Therefore, for any given state in the Hilbert space with charge Q:QĈ This implies that the charge-conjugated state exhibits a charge 2N − Q. For an even Q-charge sector, the resulting spectrum is invariant under µ → −µ, and hence the pair of charge-conjugated Hilbert spaces for {Q, 2N − Q} are physically identical. For odd Q-charge sectors, the charge-conjugated Hilbert space is only physically equivalent to the original Hilbert space once µ is replaced with −µ. Finally, for PBC it is useful to characterize the states by a winding number variable l, such that for any state |ψ in the Hilbert space with a given cutoff Λ, [Tr(U (0)U (1)U (2) . . . U (N − 1))] l |ψ is also a valid state of the physical Hilbert space with cutoffs up to Λ + l, where 0 ≤ l < ∞. Since the operator [Tr(U (0)U (1)U (2) . . . U (N − 1))] l does not commute with the Hamiltonian, it is not a symmetry of the theory, nonetheless, it provides a useful characterization of the states. As was explicitly realized in the previous section, the dimension of the Hilbert space for OBC is finite for arbitrary Λ. However with PBC, the Hilbertspace dimension grows linearly with the cutoff (once a saturating value of the cutoff is reached to accommodate all the fermionic configurations). In the linear region, the slope is obtained from the number of all fermionic configurations with charge Q = 0, which is 2N N . As a result, the complete physical Hilbert space with PBC contains copies of the particular gauge-invariant Hilbert space with winding numbers varying from 0 to Λ. Such a winding-number characterization of PBC Hilbert space is evident in the example shown in Fig. 13, where basis states 1), 2), 3), and 6) are repeated for different values of the winding numbers, l = 0, 1, 2, 3. The breakdown of the PBC Hilbert space in terms of the fermionic occupation quantum number Q and the winding number l is worked out for a related example in Table IV. Finally, it should be noted that with PBC, the theory exhibits a discrete translational symmetry, and the eigenstates of a discrete momentum operator can be formed as well, see e.g., Ref. [36] for such a classification of momentum eigenstates in the case of the lattice Schwinger model.

V. COMPARATIVE (CLASSICAL) COST ANALYSIS
A classical algorithm for Hamiltonian simulation, in general, involves three steps: I ) Hilbert-space construction, II ) Hamiltonian generation, and III ) observable computation. The Hilbert space can be constructed by identifying the theory's DOF, symmetries, and a convenient basis to express the states. In the case of LGTs, where a major portion of the Hilbert space is irrelevant, to reduce the computational cost, one needs to project to the physical Hilbert space. This entails one of the following. One may impose the (non-trivial) Gauss's laws by reformulating DOF as is the case in the LHS formulation. Alternatively, the Gauss's law constraints can be imposed a posteriori on convenient basis states. Another option is to generate the physical states by the consecutive action of the Hamiltonian on a trivial physical state such as the strong-coupling vacuum [89,90]. The associated computational cost of this step, therefore, depends largely on the Hamiltonian formulation used, as well as the algorithm itself. In the following, we analyze the first two approaches, noting the third approach is of compa-rable cost as it requires the Hamiltonian matrix to be acted on states by a number of times that grows exponentially with the system's size. After the basis states in the (physical) Hilbert space are identified, the next step of the simulation is to generate the Hamiltonian matrix. This step, obviously, depends on the formulation used as well. For example, some formulations may provide simpler operator structures, which could affect the sparsity of the matrix generated. Finally, the Hamiltonian matrix can be used to compute observables, such as spectrum, and static or dynamical expectation values of operators. This step often entails matrix diagonalization and matrix exponentiation, which can be sped up by efficient sparsematrix algorithms especially when acted on a sparse state vector.
Having introduced various formulations of the SU (2) LGT in 1+1 D and analyzed their physical Hilbert-space dimensionality with regard to the system's size, electricfield cutoff, and boundary conditions, one can now analyze the classical-simulation cost within each formulation. For this purpose, we focus on the KS formulation in the physical Hilbert space, the LSH formulation, and the purely fermionic formulation, all with OBC, and will briefly comment on the case of PBC and the purely bosonic formulation in the end.

Purely fermionic formulation
Within the purely fermionic formulation, redundant gauge symmetries are removed algebraically by an appropriate gauge transformation and after applying the Gauss's law repeatedly, as explained in Sec. II B. As such, the projection to the physical Hilbert space is essentially free, and the time complexity is: 18 Here, in principle, there is an additional cost associated with generating 4 N fermionic configurations N i=1 |f 1 , f 2 (i) with f 1,2 ∈ {0, 1}. Nonetheless, with an efficient Kronecker-product algorithm introduced in the next subsection, the Hamiltonian can be generated without the need to generate and store these fermionic configurations.

Loop-String-Hadron formulation
An efficient algorithm and its associated cost for generating physical Hilbert space of the LSH formulation with OBC goes as follows. One first generates 4 N string configurations. The cost of generating each configuration 20 can be realized as converting an integer label k associated with each configuration to a binary number, which goes as log(k). This step can therefore be conducted with the time complexity O(N ) for a lattice of N sites. The binary digits in each generated configuration are then labeled by n i and n o string numbers, by e.g. assigning them to the even and odd digits, respectively. This step is essentially free. Next comes the generation of the n l quantum numbers. As mentioned before, with OBC, the Gauss's law is used to fix this number for any given string configuration. Since n l must be fixed at all links, O(N ) number of operations is needed. There are additional O(N ) operations required to pick each generated n l and check it against the requirement of not exceeding the cutoff as well as being a non-negative integer, as explained in Sec. III D. However, this step can be simultaneously performed as generating n l quantum numbers consecutively, to reduce the cost. The total cost of generating the physical Hilbert space with the LSH formulation is therefore: Note that in order to reduce the dimensionality of the Hilbert space, one could additionally restrict the states to a given global-symmetry sector. For example, if only interested in the charge ν = 1 sector, there are O(N ) operations involved to check the ν-number of each string configuration, reducing the number of configurations needed for generation of n l quantum numbers from 4 N to 2N N ≈ 4 N √ N . Since this is not an exponential speed up as a function of the size of the system, such finer decompositions of the physical Hilbert space will not be considered in the rough estimate of the computational cost in the remainder of this section. Such symmetry considerations, however, will be advantageous in practical implementations.

Angular-momentum representation
An efficient algorithm for the generation of the physical Hilbert space of the KS in the angular-momentum basis with OBC starts by making a gauge-invariant state at site x, i.e., one that satisfies the non-Abelian Gauss's laws in Eq. (18). If there is one and only one fermion at site x, then a gauge-invariant state is obtained from the relation with the notation defined in Sec. II A. While m f only takes values ± 1 2 in these sums, each sum over m R , m Rf , and m L involves of the order of J R operations. This is because the value of J R fixes the value of J Rf to be equal to J R ± 1 2 , and that in order for the final total angular momentum to be zero, the value of J L needs to be equal to J R ± 1 2 . Now to generate a set of all possible physical states, such construction at the site should be repeated for all possible values of J R , i.e., 0 ≤ J R ≤ Λ. As a result, the number of operations required to generate a complete set of physical states at a given site is O(Λ 4 ). If, however, there is either no fermion or there are two fermions at site x, the above relation is modified by setting 1 2 → 0, the expression simplifies to only two summations, and the final number of operations required to generate a gaugeinvariant set of states is O(Λ 3 ), which is subdominant compared with the first case and can be ignored in the limit Λ 1. Note that in both cases, there is an additional cost involved amounting to checking and removing the generated J L values that violate 0 ≤ J L ≤ Λ, but this step can be checked simultaneously in the sum above, and the total asymptotic cost in the limit Λ 1 remains the same. The next step involves an N -fold Kronecker product of states in each set to connect states at adjacent sites throughout the lattice. This adds a cost with the time complexity O Λ 4N . Finally, the boundary condition on J R at site x = −1 must be imposed, along with the Abelian Gauss's law to ensure that J R and J L belonging to the same link are set equal. This can be achieved by a search and elimination algorithm, and involves an additional cost that scales as O Λ 4N , which is the conservative scaling of the number of basis states formed in the previous step. As a result, the total time complexity of generating the physical Hilbert space is: T where ≈ sign in this section is meant the approximate scaling in the limit: N 1. The time complexity if Λ is fixed to a constant much smaller than N is: T (J) An alternative algorithm can be realized by first generating 4 N fermionic configurations throughout the lattice. Each configuration generation involves a time complexity that scales as O(N ). Next, given the value of J R at the boundary x = −1, all J L and J R values can be produced throughout the lattice, given the known fermionic occupation at each site and the Abelian Gauss's law. This involves a maximum number of operations that goes as O(N + N 2 N ), since at each site and given a J R value, there may be two possibilities for the J L value as the fermion occupation may be equal to one. Now given the set of configurations for fermions and angular momenta generated, the non-Abelian Gauss's law can be implemented following the relation in Eq. (83), introducing an additional cost O(Λ 3 ). The Kronecker-product cost is the same as before but there will be no need to impose the boundary condition and Abelian Gauss's law anymore as these are already implemented. In summary, this algorithm involves a time complexity that scales as T (J) For Λ ∼ N and N 1, this later algorithm therefore is asymptotically faster than what was described earlier.
The cost, however, remains super-exponential in N this limit.

Purely fermionic formulation
In the fermionic formulation of the KS theory with OBC, the Hamiltonian becomes non-local with O(N 2 ) terms in the Hamiltonian, see e.g., the electric Hamiltonian in Eq. (28). However, a considerable advantage is that the operator structures encountered are only of the type ψ † (x)ψ(y), which make the Hamiltonian generation amenable to Kronecker-product algorithms, eliminating the need to generate and store the fermionic Hilbert space a priori. Explicitly, each ψ † (x)ψ(y) operator can be written as at the Kronecker product of I 4×4 at all sites but x and y, A 4×4 at site x and B 4×4 at site y, in an ordered manner, where I is the identity matrix, and A and B are matrices formed by the action of ψ † (x) and ψ(y), respectively, on the four allowed fermionic configurations at the respective sites. Since the on-site matrices are sparse and involve O(d i ) elements (with d i = 4 being the dimensionality of the matrices), performing the Kronecker product along a lattice of length N comes with the time complexity O(4 N ). The total time complexity of this algorithm for generating the Hamiltonian is therefore: Note that there are additional O(N ) operations involved in finding the position of x and y along the chain, but that is a subdominant cost compared with the subsequent Kronecker-product operation.

Loop-String-Hadron formulation
In order to generate the LSH Hamiltonian, first note that the dimensionality of the physical Hilbert space can be approximated by 4 N −1 at large N , estimated by comparing the asymptotic ratio of the physical Hilbert space dimension in the fermionic formulation to that of the LSH formulation, as shown in Fig. 10. On the other hand, in the LSH formulation, the Hamiltonian remains local, with the total number of operators scaling as O(N ). Generating the Hamiltonian matrix elements amounts to picking one state out of the Hilbert space at a time, acting by the Hamiltonian operator on the state to arrive at another state (which requires O(N ) operations), and find the position of the state's assigned index in a previously produced look-up table of states indices, with a time complexity that scales as O(4 N −1 ), and is therefore dominant compared with Hamiltonian operation cost. This step identifies the position of the element in the Hamiltonian matrix and its value. The total time complexity of Hamiltonian generation for the LSH formulation, therefore, scales as Note that the more efficient Kronecker-product algorithm that was applied in the fermionic case could not be taken advantage of in the LSH formulation as the one-to-one mapping between the Kronecker product of on-site physical states and the global physical state is lost given the imposition of the boundary condition and the Abelian Gauss's laws. Such a convenient feature is lost in the angular momentum representation of the physical states as well.

Angular-momentum representation
In the limit Λ ∼ N , this step introduces another superexponential cost to the Hamiltonian simulation in this basis, in addition to the cost of generating the physical Hilbert space as derived in the previous subsection.

General scaling relations
The cost of Hamiltonian matrix manipulation required to evaluate observables depends upon the dimensionality of the Hamiltonian, its sparsity, and the sparsity of the state vector in computing expectation values. For a square matrix with dimensions M × M , the matrix density is defined as the ratio of the number of non-zero elements in the matrix to M 2 . Computation of spectrum in the Hamiltonian formulation amounts to evaluating the (first m) eigenvalues, while time-dependent expectation values require the action of the exponential of a matrix on a state vector. A rigorous analysis of the computational cost of these matrix manipulations given a desired accuracy in the estimated outcome is beyond the scope of this work. Given that the goal of this section is to make a comparison among different Hamiltonian formulations of the SU (2) LGT, with a focus on select algorithms, we point out a generic time complexity measure for the operations without attempting to find rigorous theoretical bounds on the algorithms considered. The extremal eigenvalues of a matrix can be obtained using an iterative algorithm, which generates convergence to the desired eigenvalues. To demonstrate the time complexity of such algorithms, consider a Hermitian matrix with the following eigenvalue equation: A |v a = λ a |v a . The power series algorithm for extremal-eigenvalue estimation requires choosing a random initial unit vector |b = a c a |v a , |c a | 2 = 1. One repeats applying A to |b , A k |b = n c n λ k a |v a , and normalizes the output state after each iteration. This procedure suppresses the contribution of smaller eigenvalues and quickly converges to the maximal eigenvalue. The Arnoldi algorithm [110] enhances this procedure by generating a Krylov subspace with the set of output vectors |b , A |b , A 2 |b , ..., A m−1 |b that are orthogonalized against each other. To generate the lowest eigenvalues, A is replaced with A −1 in the procedure outlined. The limiting feature of this approach is the cost of matrix multiplication. As the procedure is repeated, it becomes more costly to expand the Krylov subspace, due to the growing density of the output vectors. Applying the Lanczos algorithm [111][112][113][114], a generalization of the Arnoldi algorithm to Hermitian matrices, requires O(ξρ H M 2 ) operation for an M × M Hamiltonian with matrix density ρ H . Empirical observations have shown that ξ ∼ 3 2 m, where m is the number of desired extremal eigenvalues [112].
Among methods to compute matrix exponentiation is the truncated Taylor series method [115]. For the unitary time-evolution operator, the k-th order truncated series is defined as: e −iHt k = k−1 l=0 (−iHt) l /l!. The number of terms, k, needed in the Taylor expansion to reach the desired accuracy depends on time, Hamiltonian's dimension, and on some form of a Hamiltonian norm, such as the size of its largest eigenvalue. With identical eigenvalues and similar Hilbert space sizes in all formulations considered, the number of terms to be computed will be of the same order in all cases and will not be constrained any further. The time complexity of multiplying an M × M Hamiltonian with density ρ H and an M -dimensional vector with density ρ v goes as O(ρ H ρ v M 2 ) [116], hence a kth order truncation of e −iHt on state vector |v scales as: . Although the initial state can be chosen to be sparse, the time evolution will eventually saturate the sparsity bound determined by the action of the k-th order Taylor series on the state, which results in mixing the initial state with states in the same symmetry sector.
To simplify the cost analysis, in the following the time complexity will be used as a general measure of the spectral and dynamical observables and the companying factors mentioned above will be dropped. Once again, it should be noted that due to the increase in ρ v as the system evolves, the time-evolution cost will ultimately approach the time complexity O(M 2 ).

Computing cost in different formulations
Since the Hamiltonian in the angular-momentum basis in the physical Hilbert space is the same as that in the LSH basis, the cost analysis of the observable computation is the same in both formulations. What makes a difference in computational costs of one versus the other is clearly the Hilbert-space and Hamiltonian-generation cost, as will be studied more closely in the next subsection. As a result, only a cost comparison is presented here for the fermionic and LSH formulations.  (2) LGT in 1+1 D with OBC is shown in Fig. 10. As the plots demonstrate, while the dimension of the LSH Hilbert space remains smaller than that of the fermionic representation for all values of N , the density of the Hamiltonian matrix in the fermionic formulation decreases compared with that of the LSH formulation as N increases. The determining factor in the time complexity of the eigenvalue computation and matrix exponentiation is, nonetheless, ρ H M 2 , and the ratio of this quantity between the two formulations is plotted in the lower panel of Fig. 10. As is observed, the ratio grows slowly as O(N ) as N → ∞, demonstrating the slightly higher cost of observable computation in the fermionic formulation. In the comparative analysis of the next subsection, the asymptotic cost of observable computation for the LSH formulation (hence the angular-momentum formulation in the physical sector) will be estimated simply as  (2) LGT in 1+1 D with the fermionic formulation (F), LSH formulation, and the angular-momentum formulation in the physical sector (J), as a function of lattice size N for large N .

D. Total cost and comparisons
Given the time complexity of various steps of a Hamiltonian simulation outlined above for the three formulations of SU (2) LGT in 1+1 D studied, the cumulative cost of the simulation can be estimated by adding T I , T II , and T III for each formulation from the previous subsection. Figure 14 plots this asymptotic scaling cost, i.e., in the limit where N 1. For the angular-momentum formulation, two scenarios are considered: The cutoff on the gauge excitations is set to a constant value, and the cutoff is set to N , corresponding to the saturating value of the cutoff with OBC. As already evident from the analytical scaling relations, the least time complexity is offered by the fermionic representation, while the angularmomentum formulation provides the worst scaling and is unfit for even small-scale simulations. As an explicit comparison, for N = 20, angular-momentum formulation (with Λ = N ) requires 160 orders of magnitude larger computing resources than the LSH formulation, while the fermionic formulation requires 20 orders of magnitude lesser resources than the LSH formulation.
It is also interesting to examine the cost of each step of the simulation, i.e., Hilbert-space constructions (I), Hamiltonian generation (II), and observable computation (III), for the algorithms outlined. While the conclusions are already evident from the analytic scaling forms provided in the previous subsection, they can be observed more crisply from the plots shown in Fig. 15. For both the fermionic and LSH formulations, the least costly step is constructing the physical Hilbert space, while interestingly for the angular-momentum formulation this step is the most costly step involved. Hamiltonian generation is the most costly step of the simulation for the fermionic and LSH formulations. Once the Hamiltonian is generated, as already explained, the observable computation costs roughly equal computing resources in all cases. In other words, in the angular-momentum formulation  (2) LGT in 1+1 D with the fermionic formulation (F), LSH formulation, and the angular-momentum basis in the physical sector (J), as a function of the lattice size N for large N .
Step (I) refers to Hilbert-space construction, step (II) refers to Hamiltonian generation, and step (III) denotes observable computation assuming a generic scaling for sparse matrix manipulations, see the text. The step (I) for the fermionic formulation is of O(1) with the chosen algorithm and is not shown. much of the computing effort is put into generating a Hamiltonian that is of much smaller dimensionality than the naive construction and is the same as that of LSH. That is the ultimate reason behind why the LSH formulation is developed and promoted, as it simplifies a priori the construction of the Hamiltonian and its associated Hilbert space without requiring the implementation of these steps on the inefficient angular-momentum basis states.
Before completing this section, a few comments are in order. First, it should be noted that the algorithms described can be made far more efficient by parallelizing operations when possible. The degree of parallelization will depend on the formulation considered, as well as the algorithm devised for each step of the computation. Second, the benefit of fermionic representation over the LSH formulation is weakened by the fact that only the LSH formulation can be generalized to other boundary conditions and to higher dimensions, hence its great benefit over the other formulations when the gauge DOF are present is a significant result. In particular, if the fully bosonic formulation of Sec. II C is considered within the angularmomentum basis, its computational cost will be greater than that of the angular-momentum representation with both bosons and fermions present, given the introduction of an extra U(1) cutoff in constructing the Hilbert space and the associated Hamiltonian, with a saturation value that scales as N . Last but not least, it must be now obvious that the extreme advantage offered by the LSH formulation over the angular-momentum formulation will persist with PBC as the complexity of Hilbert-space construction and the Hamiltonian generation is even more significant with PBC. While in the LSH formulation with PBC, the Hamiltonian operator still projects one basis state to another, in the angular-momentum formulation, not only the maximum number of basis states in a given physical state can scale as O(Λ N ), as shown in Fig. 9, but also the Hamiltonian operator, in general, transforms a basis state to a linear combination of many distinct basis states.

VI. SPECTRUM AND DYNAMICS IN TRUNCATED THEORY
The analysis of the previous section revealed that the most suitable formulation to construct the Hilbert space, generate the Hamiltonian, and analyze system's properties is the LSH formulation. It reproduces the same Hamiltonian as the KS case in the angular-momentum basis in the physical sector, and it does so by eliminating the need to diagonalize the Gauss's law operator in the angular-momentum basis at the level of states. As was discussed in Sec. III A, this latter step generally leads to a proliferation of terms in a given physical state when expressed in the angular-momentum basis, see Fig. 9, increasing the cost of the Hamiltonian-matrix generation. The fermionized form leads to the same Hilbert space with OBC as the KS formulation in the physical Hilbert space and the LSH formulation but involves redundancies in the spectrum as mentioned above, making it less suitable for computation than the LSH framework. It further introduces non-local fermionic interactions which make it less efficient to simulate with certain classical and quantum-simulation algorithms. The bosonic formulation, once augmented by the U(1) gauge fields with sufficiently large truncation on their excitations, reproduces the KS theory exactly, as discussed in Sec. III C. As a result, it is less favored compared with the LSH framework when it comes to computational cost. Finally, the QLM does not correspond to the KS Hamiltonian unless certain limits are implemented. As a result, the comparison between the KS and QLM spectrum and dynamics will not be meaningful away from those limits, as is the case in this study. Given these considerations, in this section, the LSH formulation (or identically KS formulation in the physical Hilbert space) will be mostly considered for the analysis of the dependence of the spectrum and dynamics on the cutoff on the gauge-field excitations, but an example of the effect of the extra gauge DOF on the spectrum of the purely bosonic formulation will be analyzed as well.
For convenience, the KS Hamiltonian in Eq. (1) can be multiplied by 2 ag 2 to yield a dimensionless scaled Hamil- where x = 1 a 2 g 2 , µ = 2m g 2 a , and N 1 , N 2 , and N 3 are defined in Sec. II. The limit x → 0 corresponds to the strongcoupling limit of the theory, while the limit x → ∞ at a fixed m g provides a trajectory in parameter space along which the continuum limit can be taken. The matrix elements of this Hamiltonian can be formed using the KS angular-momentum or LSH bases, giving rise to identical results in the physical sector, which serves as a strong check of the newly-developed LSH formulation for the 1+1 D case. While efficient classical simulations such as those based on tensor networks have enabled studies of SU(2) lattice gauge theories with a large number (hundreds) of sites [16,[90][91][92], enabling the continuum limit of the results to be taken systematically, such considerations are not the focus of this work. Instead, the aim is to arrive at qualitative conclusions regarding the cutoff dependence of several quantities of interest in small lattices where exact diagonalization is achievable on a typical desktop computer. A comprehensive study of SU (2) physics within the LSH formulation using state-of-the-art classical Hamiltonian algorithms is underway and will be presented elsewhere. Note that such methods still encounter the issues stemming from inefficient formulations of DOF and implementation of non-Abelian symmetries especially toward higher dimensions, and hence the analysis of this work is relevant in optimizing such studies.

A. Spectrum analysis
Since the gauge-matter interaction term, i.e., the term proportional to x in Eq. (89), mixes the electric-field eigenstates, it is anticipated that for large x (small ga), the truncation of the Hamiltonian matrix, imposed by the cutoff on the electric-field excitations, will introduce a more significant deviation from the exact spectrum compared with the small-x values. The ratio of the scaled energy, E , obtained from diagonalizing the scaled Hamiltonian in Eq. (89) to the corresponding x value, nonetheless, will approach a constant value as x → ∞, which represents the continuum limit once N → ∞. The same is true for the asymptotic large-x behavior of the ratios of energies at a given value of x. As a result, these ratios provide a more suitable means to analyze the scaling behavior of spectral quantities. Figure 16 plots the dimensionless quantity E 2N x for the ν = 1 sector of the KS Hamiltonian with N = 6 and with PBC, for several values of x and Λ. As is visually evident, Λ = 2 appears sufficient to obtain almost cutoff-independent energy eigenvalues for a large fraction of the spectrum associated with a given Λ when x < 1. The higher-energy states exhibit a more significant cutoff dependence, as is expected, and such cutoff dependence becomes more prominent as x increases. Nonetheless, as x → ∞, the shown quantity asymptotes to fixed values for all states in the spectrum, hence stabilizing the effect of the cutoff.
To make these features more explicit, three eigenstates of this systems are selected, corresponding to the 1st, 21st, and 283rd lowest eigenenergies. The 1st eigenenergy exists for all Λ values considered. The 21st eigenenergy is the first eigenenergy that does not appear in the spectrum of the Λ = 0 theory, as the Λ = 0 theory has a 20-dimensional physical Hilbert space given the quantum numbers specified. Finally, the 283rd eigenenergy does not exist in the spectrum of the Λ = 1 theory, as the Λ = 1 theory has a 282-dimensional physical Hilbert space given the quantum numbers considered. Figure 17 plots the quantity ∆E as a function of Λ for various x values and for the three states identified. 20 20 For x 1 as noted above, the lowest-lying energy levels are min-As is seen for all x, as Λ increases, the quantity plotted exponentially approaches the exact value, which is taken to be the value at Λ = 8. This scaling behavior may only exhibit itself for sufficiently large Λ, as is seen clearly in the case of the 21st state at x = 1. This indicates that a sufficiently large Λ is required to be able to approximate the Λ → ∞ values using these asymptotic relations. This is similar, in spirit, to the numerical procedure outlined in Ref. [117], where an iterative process is outlined to include higher-energy 'shells' in the Hilbert space in order to drive the quantities to the scaling regime. In this region, analytic relations based on RG can then be applied to arrive at cutoff-independent observables. In the same plots, the values of Λ at which the scaled energy deviates by 10% or less from the E (Λ = 8) value are denoted as dashed vertical lines. As is seen, the required Λ values to reach this moderate accuracy increases with increasing x but stabilizes at fixed values toward the continuum limit.
For higher-energy states in the theory, as is evident from the behavior of the 283rd eigenenergy, the energy corresponding to the Λ = 8 value is not a good approximation for the exact value and higher cutoffs must be included in the analysis. The scaling of the shown quantity, nonetheless, is expected to follow the same exponential form for sufficiently large Λ adjusted to the large eigenenergy considered. Finally, to demonstrate the convergence of the quantity ∆E E to a constant value at large x, the same eigenenergies picked in Fig. 17 are taken for the smallest Λ values for which those eigenenergies exist, and are plotted as a function of x in Fig. 19. After a threshold x value (that is larger for higher-energy states), the dependence of the ratio on x can be approximated by r+ae −ax , where r, a, b are constants to be fit. The values of the constant r obtained from the fits shown correspond to the asymptotic large-x value of the ratios, and are provided in the plots for each of the eigenenergies selected.
For OBC, similar qualitative features are seen for the cutoff dependence of the spectrum, as are evident from Figs. 26-28 in Appendix C. Such cutoff dependence, however, is less severe than the PBC case in general, consistent with the fact that the gauge fields are not truly dynamical DOF in gauge theories in 1+1 D with OBC. In particular, once the value Λ = N is reached, the spectrum of the theory is the same as that corresponding to Λ → ∞. Nonetheless, for large lattice sizes, it is not computationally viable to implement this saturation value, since the Hilbert space increases as e qΛ with the q value obtained from a fit to results at small lattice sizes. As a result, the scaling behavior presented for select spectral quantities in Appendix C will still serve to guide approaching cutoff-independent results from computations performed at smaller values of Λ.
imally affected by the truncation cutoff and the scale dependence for Λ > 1 is almost vanishing, so this case has not been included in the figure.  19. The same quantity as in Fig. 19 corresponding to the 283rd state in the spectrum for Λ = 2. With coarser resolution, the asymptotic (x → ∞) value of the function can be obtained from a fit to the exponential form shown in the left panel. By zooming into the large-x region of the plot, a finer structure can be observed in the data as a function of √ x, revealing an exponential asymptote to the continuum value that starts at a much larger value of √ x as shown in the right panel, with an asymptotic value that is within sub-percent of that obtained in the left plot. The numerical values associated with these plots are provided in Supplemental Material.
As discussed in Sec. III C, the bosonized KS theory exhibits a physical Hilbert space that is identical to that of the KS theory with fermions in the limit where the cutoff on the U(1) electric-field excitations of the extended U(2) theory is equal or larger than the saturation value Λ 0 = N . To examine the effect of the extra U(1) cutoff, consider, as an example, the bosonized theory with OBC and with N = 8 and ν = 1. In this sector, the full physical Hilbert space is reached once the cutoff on the SU(2) electric-field excitations is set to Λ ≥ 6, giving rise to a total of 784 states. However, such a value is only achieved if Λ 0 ≥ 6 as well. In particular, the number of allowed states in the physical Hilbert space is 0, 0, 208, 576, 759, 783 for Λ 0 = 0, 1, 2, 3, 4, 5, respectively. The values of the scaled energies obtained from H = 2 ag 2 H for a small and a large value of the coupling x (divided by 2N x) are plotted in Fig. 20 for Λ 0 = 2, 3, 4, 5, 6, demonstrating a strong Λ 0 dependence in the spectrum. For comparison, the total number of states whose scaled energies in the shown unit are less than or equal to 10% of the exact value (corresponding to Λ 0 = 6) are shown in the figure.

B. Dynamics analysis
Dynamical quantities, i.e., those obtained from timedependent expectation values, exhibit cutoff dependence as well. In fact, intuitively one expects that when the system evolves, the unitary evolution operator e −iHt ultimately mixes the initial state with all states within the Hilbert space with the same (conserved) quantum numbers [118], where H denotes the Hamiltonian of the system and t is time. In other words, various parts of the system entangle quickly as the system evolves. This feature is at the heart of the inefficiency of tensor-network methods in dynamical studies as the evolution time increases. As a result, quantities evolved to large times compared with 1/||H|| suffer from Hilbert-space truncation errors more significantly. Here, ||H|| is some form of a Hamiltonian norm, such as the absolute value of the largest energy eigenvalue.
To demonstrate this point quantitatively in the case of the KS SU (2) LGT in 1+1 D, consider a select quantity obtained from where |vac(x = 0) denotes the strong-coupling vacuum, t is a dimensionless time: t ≡ t 2x , and all other quantities are defined after Eq. (89). Here t is time in units of a. The quantity in Eq. (90), therefore, counts the occupation number of particles (as opposed to antiparticles) on the lattice as a function of time, and is identically zero for the strong-coupling vacuum. Figure 21 plots N (t) as a function of Λ for various values of x and t in the KS Hamiltonian within the physical Hilbert space with N = 6 and ν = 1 with PBC. Note that t is the absolute time, and hence N (t) at different x values can be directly compared at a fixed t. At small x values, the mixing among states expressed in the electric-field basis is small, and the dependence of the time-evolved quantity on the cutoff can be approximated by Ae −BΛ + C even for rather small Λ. As x becomes larger, larger values of Λ are needed to enter this asymptotic scaling region, at which point the Λ → ∞ value of the quantity can be estimated based on the exponential scaling. This feature stabilizes as the function of x for large x to ensure convergence to a continuum limit (once N → ∞ is taken). Not surprisingly though, at larger values of t, it takes a large value of Λ to enter the asymptotic region. For example, at t = 100 no convergence is seen in the quantity considered for x = 100 and x = 400 in the range of cutoffs considered. This example demonstrates that studying dynamics in this theory with high precision will be computationally demanding, given the extreme sensitivity to the cutoff considered at long evolution times.
Another feature of dynamical quantities is that their continuum limit seems to be achieved more slowly, requiring more computational resources compared with static quantities. Although such a limit should be considered as an ordered double limit in which N → ∞ and x → ∞, this point can be demonstrated at a fixed N and for various x values. Examples of this feature are shown in Fig. 22, in which the quantity N (t) is plotted as a function of x for the KS theory with PBC and Λ = 4 at a select time, t = 100. This observation is consistent with the slow convergence of high-lying energies in the spectrum as a function of x (see Fig. 19), considering that time evolution leads to contributions to observables from all states in the same symmetry sector. Nonetheless, if percent-level precision is needed in this quantity, one would not need computations to be performed at extreme values of x.  • t = 20 • t = 20 It is worth noting that the point at which the quantities can be described by exponentially converging functions as Λ → ∞ or as x → ∞, depends on the size of the system, the choice of m/g, and the quantity considered. To formally understand the underlying mechanism for such an scaling and its breakdown, beyond the empirical observations of this work in small systems, is the subject of a future investigation and potentially has its root in the Nyquist-Shannon sampling theorem [77]. A more systematic understanding of Hamiltonian renormalization can shed light on this question, as we have alluded   to in Sec. VII.

VII. CONCLUSIONS AND OUTLOOK
The present paper provides a stepping stone to establishing the most efficient Hamiltonian formulation(s) for the case of a non-Abelian gauge field theories in 1+1 Dimensions, with implications for other non-Abelian gauge theories and for higher dimensions. The theory studied in this paper is the SU (2)  The dimension of the Hilbert space grows quickly with the size of the lattice, N , as ∼ e pN . In the angular-momentum basis, a naive counting of the full Hilbert space reveals that the exponent grows logarithmically with the cutoff on the electric-field excitations. With open boundary conditions, a saturation value of Λ = N + 2 0 ensures that the full untruncated theory is achieved, but at the cost of expediting the growth of the full Hilbert space by ∼ e N log N . The physical Hilbert space in the same representation still grows exponentially with the system's size, but with a p value that approaches a constant of O(1) as a function of the cutoff. The empirical values are obtained in Sec. III. Additionally, as a function of the cutoff on the electric-field excitations, the dimension of the physical Hilbert space approaches a constant with OBC but grows linearly with PBC. The considerable reduction in the Hilbert-space dimension by restricting computations to the physical sector comes at a significant cost, since imposing the non-Abelian Gauss's laws locally amounts to a super-exponential cost in generating the Hilbert space. There is a similar cost associated with generating the Hamiltonian matrix, simply due to the fact that each physical state in the original angular-momentum basis is a linear combination of many basis states, with the number of terms in given physical states growing exponentially with the cutoff. Furthermore, the action of Hamiltonian on such states generates other populated linear combination of states. Given these features, the angular-momentum basis of the KS theory appears the least appealing framework for Hamiltonian simulation.
In 1+1 D, the gauge DOF are non-dynamical with OBC and can be fully fixed by the boundary value of the electric field. This means that restricting to the physical Hilbert space comes as no computational cost. Nonetheless, the fermionic basis states lead to redundancies, and while the spectra of both the original theory in the physical sector and the fermionic theory are identical, there are degeneracies in the spectrum in the latter case. The difference in the dimension of the physical Hilbert spaces in these formulations is not significant though, and the ratio remains constant as the system's size grows. The larger Hilbert space is compensated by a slightly faster decline in the Hamiltonianmatrix density of the fermionic formulation, and the cost of computing observables, i.e., eigenvalue evaluation or matrix exponentiation, is only O(N )  (2) LGT in 1+1 D studied in this work, and the connection among them.
higher in the fermionic formulation compared with the angular-momentum formulation in the physical sector. Nonetheless, the purely fermionic formulation is the least costly among all the formulations considered in the limit of large lattice sizes, as expected. The drawback is that this formulation has no analog in higher dimensions as the boundary conditions and Gauss's law constraints are not sufficient to fully remove the gauge DOF. As a result, other formulations with better generalizability perspective and similar computational-resource requirements are favorable.
Alternatively, the fermionic DOF can be removed at the expense of enlarging the gauge group to U(2) and introducing additional complexity due to the need for sufficiently large cutoffs on the additional U(1) electric field excitations, and for appropriate encoding of the fermionic statistics. This theory, once the saturating value of the U(1) cutoff is imposed with OBC, and the Gauss's laws associated with both the U(1) and SU(2) symmetries are solved, has an identical physical Hilbert space to that in the original KS theory. All the complexities arising from imposing the non-Abelian constraints in the angular-momentum basis remain relevant to this formulation, and hence no particular advantage is gained by the purely bosonic formulation with classical Hamiltonian-simulation algorithms. Such a formulation, nonetheless, is argued to be beneficial in digital quantum simulations as it avoids nonlocal fermion-to-hard-core-boson transformations.
The complexities with the angular-momentum representation of the KS theory and the increasing mixing of the local basis states in the physical Hilbert space are avoided altogether in the LSH formulation. Here, the building blocks of the Hamiltonian are operators that act directly on local gaugeinvariant basis states, namely strings and loops, and only an Abelian Gauss's law on the link between the lattice sites is left to be imposed a posteriori. Furthermore, the action of each operator in the Hamiltonian maps a given state to one and only one state, a feature that is absent in the angular-momentum basis in the physical Hilbert space. While the string quantum numbers can only take values 0 or 1, the loop quantum numbers are non-negative integers, and their spectrum must be truncated in practical applications. With the same truncation cutoff, the Hamiltonian matrices in the physical Hilbert space with the angularmomentum and LSH bases are identical, nonetheless the cost of generating such a matrix is exponentially suppressed in the LSH formulation. Similarly, while in the original angular-momentum basis, the construction of Hilbert space is the most costly step in the simulation algorithm, this step is the least costly with the LSH formulation, and offers up to hundreds of orders of magnitude reduction in the cost of the simulation compared with the original formulation for lattices with tens of sites, as demonstrated in Sec. V. Although the purely fermionic formulation remains slightly less complex and computationally less costly to simulate (de-spite its all-to-all fermionic interactions), the LSH framework is already generalizable to higher dimensions, and is shown to retain the same crucial features as the 1+1 D case. Namely, a complete basis of local gauge-invariant states exists to eliminate the need for the imposition of non-Abelian Gauss's laws (that are more complicated to implement in higher dimensions) at the level of basis states, and that the Abelian Gauss's law is still implemented in only one dimension, thanks to a point-splitting procedure developed in Ref. [100,119]. As a result, this framework appears to be the most appealing among all the formulations considered for the purpose of Hamiltonian simulation of the SU (2) LGT. The generalization of such an efficient formulation to the SU (3) LGT coupled to matter is, therefore, an important next step in approaching the ultimate goal of Hamiltonian simulation of QCD.
Hamiltonian-truncation methods, common in a variety of problems from condensed-matter physics to nuclear physics, have been well developed over decades and have a direct connection to the concept of Renormalization Group. In particular, Wegner [120] and Glazek and Wilson [121] independently studied the concept of renormalization within a Hamiltonian approach. They established powerful methods such as similarity renormalization group to continuously flow the transformed Hamiltonian to a banded-diagonal form such that the high-energy and low-energy modes in a given theory are maximally decoupled. This procedure leads to smaller truncation errors resulting from considering only a partial sector of the Hilbert space in generating and processing the Hamiltonian matrix. More concretely, such approaches put the system in the asymptotic scaling regime such that the limit of Λ → ∞ can be taken using asymptotic formulae, see e.g., Ref. [117]. In other words, a Hamiltonian truncation approach can be used to derive the RG equations for the parameters of the theory, promoting the truncated bare Hamiltonian to an effective renormalized Hamiltonian. A dedicated study of renormalized Hamiltonian gauge theories which extends upon previous work is a necessity as one moves toward using the Hamiltonian approach in the simulations of gauge theories. In this first study, however, we have only presented empirical observations for the scale dependence of the spectrum and dynamics in the SU(2) gauge theory in 1+1 D in small lattices to illuminate several seemingly general features. In particular, the exponential dependence of observables on the excitation of gauge fields is established for energies and given dynamical expectation values, but a suf-ficiently large cutoff is needed to enter this scaling regime for high-energy states, as well as the longtime-evolved states. Approaches in constructing a renormalized theory in the Hamiltonian formulation will be studied in future work to shed light on the interesting empirical features observed. In this context, the formal analysis of digitization effects in bose-fermi systems [75,76] and in scalar field theories [77] will potentially be relevant.
Despite the steady progress in advanced simulation algorithms such as tensor networks, the enormous cost of Hamiltonian simulation, even with the most efficient formulations identified, suggests that the quantum simulation of LGTs of interest may be the ultimate path to reaching accurate and largescale dynamical simulations. It is, therefore, important to extend the analysis of this work to verify the efficiency of the LSH formulation compared to other formulations in the context of quantum simulation, and to find the most efficient quantum algorithms for simulating it, similar to studies conducted to date for various quantum field theories, such as scalar field theories [24,77] and the lattice Schwinger model [53]. The insights obtained from the thorough comparative study of this work are the key to embarking on such investigations in a quantum-computing setting. Analog quantumsimulation protocols within the LSH formulation are being developed, see e.g., Ref. [122], and similar developments are planned for digital protocols. In particular, explicitly solving Gauss's laws in a digital quantum computation and projecting out the physical Hilbert space of the non-Abelian SU (2) LGTs in any dimension has only been achieved using the LSH formalism [40]. Presumably, a mixture of various formulations, such as purely bosonic formulations or partially fermionic formulation with the LSH basis, could provide more computational benefits in this context, but further investigations need to be conducted to reach accurate conclusions. E R (x) (A22) E L (x + 1) (A23) 24. Rishon constituents at two adjacent staggered sites x and x+1 in a one-dimensional lattice within the QLM formulation of SU (2) LGT.
This construction re-expresses the different parts of the SU (2) LGT Hamiltonian in terms of rishon operators as s (x + 1) ↑ (x − 1)n ↑ (x) +n wheren [τ ] s (x) ≡ĉ implying that the local symmetry is instead SU (2) ⊗ U (1). Hence, in order to recover the SU(2) gauge symmetry of interest, one needs to add the following U(1) symmetry-breaking term to the Hamiltonian, It should be noted that the U(1) symmetry-breaking term in Eq. (A10) is only non-vanishing if the total number of fermionic rishons on each link is equal to two [123]: The QLM is only equivalent to the KS SU(2) gauge theory in the continuum limit and through a dimensional reduction [81,82]. As a result, the spectrum and dynamics of the theory outlined here are not the same as those of the KS LGT obtained in this paper, even when the same lattice size and couplings are considered, a feature that can be verified by numerical computations. Nonetheless, it is interesting to ask whether the QLM is computationally efficient in the context of Hamiltonian simulation. While a thorough analysis of the computational cost of constructing the physical Hilbert space, generating the Hamiltonian, and computing observables with the QLM are not analyzed in this paper, one important observation can be made regarding the dimensionality of the physical Hilbert space of the QLM. Figure 25 plots the ratio of the dimension of the physical Hilbert space in the QLM to that in the KS formulation (or equivalently the LSH formulation) when the cutoff is set to its saturating value with OBC. As is seen from the empirical fit to the values for the first lowest N values, this ratio grows exponentially with N , and so from a computational standpoint, such a finite-dimensional representation of the SU (2) LGT is still costly. Note that here, one 21 To compare with Ref. [92], note that E L ≡ −J L and E R ≡ J R . should also account for a comparable computational cost to the LSH formulation in generating the Hilbert space: there are six types of fermions present locally, requiring 2 6 configurations to be generated at each site. Additionally, the Gauss's law constraints and the fixed-rishon number per link must be imposed when constructing the Hilbert space. It is worth noting that while in the LSH formulation, the infinite-dimensional bosons are present (and are cut off at some finite value), with OBC their value can be fixed given the string configurations. In higher dimensions, there may still be an advantage in working with the finite-dimensional QLM since with other formulations, the Hilbert space grows with the cutoff on the link quantum numbers as these can no longer be fully fixed with boundary conditions. Nonetheless, the LSH formulation appears to be a competitive formulation of the SU (2) LGT, that not only is equivalent to the original KS theory, but also its economical resource requirements can bring it to the same footing as the QLM when it comes to quantum-simulation proposals on analog and digital simulators, an avenue that will be explored in the upcoming studies.         Note that for Λ = 6, the full Hilbert space with this quantum number is achieved, nonetheless the quantity shown reaches the exact values with Λ < 6 in all cases. As a result, no fit to the Λ dependence of the data is performed in contrast to the PBC case. The numerical values associated with these plots are provided in Supplemental Material.