Randomized measurement protocols for lattice gauge theories

Randomized measurement protocols, including classical shadows, entanglement tomography, and randomized benchmarking are powerful techniques to estimate observables, perform state tomography, or extract the entanglement properties of quantum states. While unraveling the intricate structure of quantum states is generally difficult and resource-intensive, quantum systems in nature are often tightly constrained by symmetries. This can be leveraged by the symmetry-conscious randomized measurement schemes we propose, yielding clear advantages over symmetry-blind randomization such as reducing measurement costs, enabling symmetry-based error mitigation in experiments, allowing differentiated measurement of (lattice) gauge theory entanglement structure, and, potentially, the verification of topologically ordered states in existing and near-term experiments. Crucially, unlike symmetry-blind randomized measurement protocols, these latter tasks can be performed without relearning symmetries via full reconstruction of the density matrix.

Randomized measurement protocols are based on changing the basis via unitaries U drawn from an appropriate ensemble E, i.e. ρ → U ρU † , followed by measurement in this basis, and classical postprocessing or quantum variational techniques to determine quantities of interest. 1The effectiveness and cost of a scheme depends on the choice of E and the particular quantities one wants to compute.For instance, estimating expectation values, Tr[O m ρ], of M q-local operators {O m } M m=1 can be done with qubitlocal random rotations, O(3 q log(M )) measurements, and efficient classical postprocessing [9].For nonlinear quantities, i.e. those that depend on ρ k (k ≥ 2), one approach is to choose an ensemble E that forms an (approximate) unitary k-design [21,22].
A fundamental problem is that standard randomized measurement protocols do not account for the symmetry structure of states, e.g., by randomizing unnecessarily over unphysical bases.In particular, a (classically known) symmetry Ŝ of ρ = s ρ s ([ Ŝ, ρ] = 0), where s labels Ŝ-eigensectors, is lost after randomization.Therefore, it is advantageous to perform symmetry-conscious randomized measurement by using randomizing unitaries such that [U, Ŝ] = 0 for all U ∈ E. This information is compactly summarized in Fig. 1, comparing symmetry-ignorant versus the symmetry-conscious random circuits U proposed in this manuscript.Our goal is to systematically study the construction of symmetry-conscious randomizing circuits and their use cases for exploring physical phenomena.The primary results and findings of our work are as follows (a) In Section 2, we discuss a comprehensive approach to symmetry-conscious randomization for qubit-based models with inherent symmetries.We discuss one application, symmetry-conscious unitary k-designs, as an approach to global randomization.This offers several advantages, including a significant reduction in measurement complexity and allowing for symmetry-based error mitigation [92].
(b) In Section 3, we focus on a concrete LGT example to demonstrate the application of symmetryconscious unitary k-design based randomization.
Our main finding is that symmetry-conscious designs enable the simultaneous measurement of both the distillable and symmetry components of entanglement, a capability that symmetryignorant designs lack.
(c) In Section 4, we present our key result, a protocol designed to detect topological order (TO) experimentally by assessing the gap of the entanglement spectrum (ES).The ES is a presentation of a state in terms of an entanglement Hamiltonian in accordance with Li and Haldane's conjecture for TO states [63].This has attracted significant attention in recent years, but has so far remained within the realm of theoretical exploration.To address this, we leverage symmetry-conscious random measurements in combination with entanglement Hamiltonian tomography techniques.This approach represents a promising step towards experimental identification of topological phases which have already started to materialize in experiments [134,135].
The realization of an approximate unitary k-design is not necessary for our protocol, but it is a convenient approach to illustrate the scheme given that Haar-random measurement channels can be easily inverted.Any tomographically complete symmetry-respecting scheme could serve the same purpose.
Finally, in Section 5, we summarize our results and further discuss applications and extensions.The manuscript is supplemented by several appendices where we discuss various details of the employed numerical techniques.

Symmetry-conscious randomization
We begin with a general description of symmetryconscious unitary circuit construction-i.e., circuits that preserve the symmetry structure of an input state ρ = s ρ s , U = s U s .Our approach, partly inspired by Refs.[21,22], uses the fact that symmetryconscious randomization can be performed using local generators that are present in the Hamiltonian of a given system.The first step is to identify local q-qubit generators (1 ≤ q ≤ m) for symmetryconscious Haar random unitaries from amongst the m-body terms of the relevant physical Hamiltonian H.For pragmatic reasons, such terms are good candidates for experimentally-implementable symmetryconscious randomization: if an experimentalist can engineer the local Hamiltonian terms for the purposes of (Trotterized) time evolution, with adequate control of coupling strengths, they can also perform randomization generated by such terms.The next step is to ensure that such terms can generate (Haar-)random unitaries over an m qubit Hilbert space within every symmetry sector.The goal is to use this mlocal symmetry-conscious randomization as a building block for randomization over the full Hilbert space.General U (2) rotations are in either a 'ZXZ' or 'ZYZ' generator decomposition with angles α, β and γ (and a global phase).In an m-qubit circuit one must identify the corresponding operators that act on q qubits and which 'embed' 1-qubit rotations in every symmetry sector.
Different possible families of circuits can be generated by different choices of arranging these blocks within a larger circuit; the m qubits interacting via a m-local Hamiltonian term, need not match the same (typically, geometrically local) m qubits in the corresponding Hamiltonian.While we focus primarily on generating families of unitary circuits that form symmetry-conscious unitary k-designs, one could also consider shallower circuits, e.g., for measuring local observables.
To leverage global Haar randomness, we determine if the selected local, symmetry-respecting terms are sufficient to generate global, symmetry-respecting Haar random unitaries within each symmetry sector.This is not generically the case [17].Verifying that sufficient randomization is, indeed, possible with the selected set of terms must be done on a case-bycase basis and accounts for the primary challenge in extending this approach to new systems.Here, we will demonstrate that such local terms are sufficient to generate approximate unitary k-designs for a few different examples of interest (with an emphasis on LGTs).

A symmetry-conscious unitary k-design example: particle number symmetry
We illustrate our approach first for spin systems with particle number symmetry, ŜN ≡ j (σ z j + 1)/2, where j labels a lattice site and σ a j (a = x, y, z) are Pauli matrices acting on the site.We emphasize that symmetry-conscious randomization schemes for particle number symmetry have been considered before [12,17,136], but it serves as a simplified setting in which to develop the particular approach we will take to symmetry-conscious randomized measurement protocols.The lessons we learn here can then be extended to LGTs.
Symmetry-conscious randomization is achieved with components consisting of (m=2)-qubit gates of the form where the rows and columns label (from top down) the ↓↓ (0-particle), ↑↓, ↓↑ (1-particle) and ↑↑ (2-particle) sectors; If α, β, and γ are selected such that [u 1 ] is drawn from a circular unitary ensemble (CUE) [24], θ and ϕ are drawn evenly from [0, 2π), then Eq. ( 1) is a block-structured unitary acting Haar-randomly on the blocks of equal particle number. 2 It has the following circuit realization: In line with our general approach, the intuition behind Eq. ( 3) comes from inspecting a particle number conserving Hamiltonian: U j1,j2 (β) is essentially a hopping term between j 1 and j 2 with amplitude β/2, the z-rotations before and after are densitydensity interactions with amplitudes α/2 and γ/2, respectively.Together, these generate SU (2) in the 1-particle block.The relative phases θ and ϕ are because we wish to embed U (2), not SU (2), and because the 0-and 2-particle blocks should be U (1)randomized. (A global phase is irrelevant.) As a first application, we now focus on global randomization, arranging these components into approximate unitary k-designs.To do so, we construct an ensemble E of N E n-qubit symmetry-conscious random circuits, as shown in Fig. 2, consisting of ℓ layers where n(n − 1)/2 pairs of qubits ((n − 1)(n − 2)/2 if n is odd) are randomly connected by Eq. ( 1).For sufficiently large ℓ and N E , this realizes an approximate unitary k-design in every symmetry block s.
We verify that such circuits yield an approximate k-design numerically for k = 2 by considering the following moments, separately for every s, ⟨U s ij U s * kl ⟩.For a 2-design it holds that while for a 1-design [137,138].We numerically simulate an ensemble of N E circuits and compute their deviation from a 2-design ϵ as the absolute difference between Eq. (5) and Eq.(4), averaged over indices and multiply by d s (d 2 s − 1) to make ϵ dimension-independent.Up to a rescaling of the approximation ratio, this choice of quantifying the error between our random ensemble and a 2-design is equivalent to other standard choices, such as the diamond norm or the the frame potential.This is demonstrated in Appendix A.
Fig. 3 A primary concern for constructing symmetryconscious random circuits is circuit complexity-in particular, the required number of layers ℓ needed to reach a given precision as a function of the system size.To investigate this, we conducted numerical simulations on systems involving up to n=10 qubits.In Appendix B, we provide a detailed account of the number of layers necessary to represent a k-design (specifically a 2-design), comparing our findings with exact Haar random sampling results.For the system sizes we can consider numerically this convergence is very rapid.For instance, for n = 10 qubits with N E = 8192 random unitaries the difference between sampling from our circuits versus sampling directly from the CUE is already very small within ℓ ≈ 15 layers for all symmetry sectors.The exact value of this error floor is set by N E , falling off with the expected ∼ 1/ √ N E scaling, see Appendix B. As our numerical methods are restricted to relatively small systems and considering that k-designs serve as one application rather than the core focus of our work, we refrain from asserting an exact analytical form of the error scaling that holds true for large n.Nonetheless, it seems reasonable to conclude from our data that the behavior scales as a low-degree polynomial, and that symmetry-conscious circuits require comparable depths when compared to symmetry-agnostic circuits.For convenience, we continue to work in the large layer limit, in practice ℓ = 128, thus massively overdoing the actual number of layers needed.We will return to analytic estimates in future work.

Estimating k-purities from k-designs
Next, we explore measuring k-purities using a prototypical model with particle number symmetry; we consider the following spin Hamiltonian in (1 + 1)d, on a lattice with N sites labelled by j ∈ [0, N − 1] and with periodic boundary conditions (PBCs); σ ± j ≡ 1 2 (σ x j ± iσ y j ), a is the lattice spacing and m a mass term.Eq. ( 7) conserves particle number, [H, j 1 2 (σ z j +1)] = 0.As inputs for the particle number symmetry-conscious scheme in Fig. 2, we numerically determine ground states of H via exact diagonalization [139].Because these are simple, unentangled computational basis states in the limit 1/2a → 0, we work in the opposite limit with m • a = 0.05 (the ground states are half-filled).We focus on a subsystem A with N A = N/2 sites.Symmetries of ρ A are particle number s ≡ n A ∈ [0, N A ] in the subsystem, noting that if particle number is fixed globally, We extract k-purities, Tr(ρ k A,s ) by measuring the probabilities P U (b, s) of bitstring b (and symmetry sector s) with N M shots in the basis defined by the random unitary U ; k-purities for k ≥ 2 are directly related to the k-Rényi entropies S (k) Following the approach taken in Refs.[21,22,140], stochastic moments ⟨P U (b, s) k ⟩ ≡ (1/N E ) U ∈E P U (b, s) k are related to k-purities via [22], Tr ρ j A,s aj (8) where A,s refer to the k-moments of the probabilities and k-purity per symmetry sector s, respectively, with b∈s P U (s, b) = Tr[ρ A,s ] ≡ p s ≤ 1 and s p s = 1.We assume an ideal quantum machine, and the total measurement cost is In Fig. 3(b), we show the required shot number N M to measure S (2) A,s with precision better than ϵ = 0.05, for fixed N E = 1,428 and ℓ = 128 layers, comparing different system sizes and symmetry sectors.The fit reveals an approximately linear dependence on d s , error bars indicate standard error of the mean for seven independent trials.(In Fig. 11 of Appendix B we investigate the N E dependence in the infinite shot limit.)Together, our results show that the cost of the symmetry-conscious approach, N E • N M , is proportional to d s instead of d A (the size of the subsystem).In many cases d s ≪ d A , yielding a significant advantage.We shall see that this advantage is exponential for lattice gauge theories.3 Lattice gauge theory entanglement Lattice gauge theories (LGTs) are systems with an extensive number of local constraints in the form of Gauss laws defining a physical sub-Hilbert space.One of the key applications is the simultaneous measurement of distillable and symmetry components of entanglement.This capability is exclusive to symmetry-conscious designs and is not attainable through symmetry-ignorant schemes; it would otherwise require full state tomography.We consider Z 2 LGT coupled to staggered matter in 1+1 dimensions, with Hamiltonian,
Fig. 4(b) and (c) depict random measurement strategies for this model, we investigate a N A site subsystem with boundary conditions ending in a matter site.The first strategy, (b), is an extension of the circuit in Section 2, Eq. ( 3), made gauge invariant by introducing Wilson lines connecting sites j 1 and j 2 , (11) and with circuit representation, )}; the angles α, β, γ, ϕ and θ are randomly drawn as before.The abbreviated middle qubit bundle in Eq. ( 12) refers to qubits representing gauge links (orange in Fig. 4), j 1 and j 2 start and end on matter sites (blue in Fig. 4).This results in a 2 + |j 1 − j 2 | qubit unitary which is not a feasible strategy near-term.Because of this we will focus on a 3-qubit unitary strategy, Fig. 4(c), at the cost of somewhat deeper circuits to obtain k-designs.We demonstrate numerically that these circuit form approximate 2-designs in Appendix C, and focus here on measuring entanglement.We note that (b) and (c) generate k-designs in every symmetry sector, despite the fact that no explicit 2-qubit entangling operation is performed between gauge sites (orange).This is a consequence of gauge symmetry: For this simple model the gauge link degrees of freedom are not truly independent because they could have been eliminated using Gauss' law. 3 3 A caveat is that with PBCs a gauge zero mode cannot be

Measuring distillable and symmetry components of entanglement
We now put our randomizing circuits to use by demonstrating their utility for determining the entanglement structure of ground states.In the literature, this work included, the term "entanglement structure" is used fairly ambiguously to denote anything characterizing entanglement beyond entanglement entropies.This includes the separation of entanglement entropies into distillable and symmetry entanglement, but, more generally, also includes the structure in terms of an entanglement Hamiltonian, eigenvalue spectrum (the so-called Schmidt spectrum), and symmetries of a reduced density operator.We will consider all of these in this section.
In particular, we consider a bipartition ρ A ≡ Tr Ā(ρ) of the lattice with ) but additionally has symmetries beyond those of the non-gauge spin model.They are illustrated in Fig. 4(d): Because of Gauss' law, on the right (left) boundary j = j R ≡ N A (j = j L ≡ 0), out-going (in-going) electric fields can be written as These operators are symmetries of ρ Gauss' law respecting) and result in a four-block symmetry structure.Together, we label the sectors as s ≡ n s L ,s R A where s L/R = ↑/↓ and n A ∈ [0, N A ] is particle number.
In Fig. 5(a), we show k-purities, Tr[ρ k A,s ] reconstructed by inverting Eq. ( 8), measured for N = 10 (N A = 5), m • a = 0.1, e/m=8 and ℓ = 32.We plot them as a function of N E ; measurements are obtained in the infinite shot limit (N M → ∞).The k = 1 result p s ≡ Tr[ρ A,s ] is recovered exactly by design, higher (k ≥ 2) purities are recovered for sufficiently large N E (data is shown up to N E = 2 14 ).Not shown is the cost in N M to obtain constant error which, as in Sect.2, scales with d s .In Fig. 5(b), we estimate the von-Neumann entropy per symmetry block, where ρA,s ≡ ρ A,s /p s .We make use of a 4 th order finite-difference approximation of the derivative to derive Eq. ( 14) from the measured k-purities; error bars are obtained from comparing a 4 th and 3 rd order derivative.Finally, Fig. integrated out.A k-design acting on the full system, not just a subsystem, would include a modification which depends on how the remaining bosonic mode is digitized; see e.g.[141,142] for bosonic random measurement schemes.Finally, we note that a difference in the number of unitaries required for approximating a 2-design versus the seemingly lower requirements for measuring 2-fidelities evident in Fig. 5.It is important to note that the accuracy of reproducing the latter is inherently dependent on the state itself.In contrast, the analysis we perform in section 2 provides an upper bound, ensuring convergence for any input state.

Classical shadows
Next, we explore classical representations of ρ A , starting with the classical shadow formalism of [9].The basic idea is to randomize ρ A → U ρ A U † and perform a computational basis measurement yielding a bitstring b from which the sector s can be read off.A symmetry-conscious shadow is U † s |b, s⟩⟨b, s|U s , where the subscript s indicates that one works in block s with dimension d s .The ensemble of random rotations yields a CUE-random quantum channel and with many measurements and, consequently, many shadows one obtains classical sectorwise state representations by taking the expectation value where LGT serves as a useful case study, we next consider 2 + 1 dimensions where our ability to extract LGT entanglement structure provides a window to studying topologically ordered systems.We add that symmetry-conscious shadows based on deep scrambling k-designs are not ideal for estimating local observables, a single layer scheme is better suited for this task.

Experimental verification of topological phases
Finally, our scheme can facilitate the experimental identification of topological phases, an inherently difficult task because of the non-detectability of such phases by local measurements [134,135].Our proposed strategy is to combine symmetry-conscious randomization with classical shadow or entanglement Hamiltonian tomography to measure the entangle- ment gap of topologically ordered states, drawing on the foundational work by Li and Haldane [63].
For convenience, we will continue utilizing k-designs in a large layer limit-much larger than necessary.It is important to reiterate that these choices are not integral to our protocol but are adopted for illustration and convenience, because of the simple channel inversion associated with them for classical shadows and for easier comparison across different schemes.Any tomographically complete scheme is suitable for this task, including shallow depth circuits that are more likely to be realistically employed in near-term experiments.
We focus on a model related to the toric code, used in [134], ).This model is graphically illustrated in Fig. 7(a), consisting of spin 1/2 degrees of freedom placed on the links (j, b) of a two-dimensional rectangular N x × N y lattice, where j = (j x , j y ) and b = x, y is the direction of a link.The Hamiltonian is given by and superselection sectors are determined by Gauss laws ([H, Ĝj ] = 0), Ĝj ≡ where the product is over neighboring links to each lattice site j, see Fig. 7(a) for illustration where l runs over links labelled (1,2,3,4); physical states obey Ĝj |ψ⟩ = +1 |ψ⟩.We assume periodic boundary conditions along the x-direction, and fixed or periodic boundary conditions (BCs) in the y-direction; for fixed BCs, Eq. ( 18) involves three links at the yboundary.
Our primary result for this model, and, consequently, for the paper, is to demonstrate that symmetry-respecting randomized measurement schemes enable an appealing and experimentally tractable route to verifying topological order.Ultimately, we shall see that by leveraging the approach known as entanglement Hamiltonian tomography (EHT) [28,30,32] and an ansatz for the Entanglement Hamiltonian inspired by the Bisognano-Wichmann theorem [144,145], we can accurately extract much of the Schmidt spectrum of reduced density matrices ρ A of the ground state in this model, allowing us to detect topological order via the presence of a so-called 'entanglement gap' in the spectrum.

Symmetry-respecting randomized circuits
We are interested in the entanglement properties of ρ A ≡ Tr Ā(ρ), where A is a bipartition obtained by separating the lattice along the x-direction, with entanglement cuts at j x = 0 and j The boundary is such that system A contains the y-direction links at j x = 0 and j x = N A x − 1, see Fig. 7(c) or [56] for more details, with We seek a gauge invariant, symmetry-respecting family of random circuits that form k-designs acting on every symmetry sector of ρ A .Our ansatz is illustrated in Fig. 7(b), made of alternating even-odd half-layers, consisting first of R z (α) ≡ exp{iασ z i } rotations randomly placed at one side i ∈ {a, b, c, d} of a plaquette (orange squares), followed by , and again by a R z (γ) ≡ exp{iγσ z i }, placed at the same random i ∈ {a, b, c, d}.For every plaquette, the angles α, β, γ are drawn according to a ZXZ decomposition of a CUE matrix, where the phase δ is irrelevant.In Appendix D, we verify numerically that these random circuits form an approximate unitary 2-design (k-design).
Before continuing, we discuss the symmetries of ρ A ≡ Tr Ā[ρ], depicted in Fig. 7(c).The Gauss laws at entanglement boundaries allow us to write Here, σ z 4 is the 'electric field' operator just outside A, and Gauss law respecting).We label simultaneous eigensectors of all Eq.(20) as s ∈ {↑, ↓} 2 2Ny ; an example is s = ↑↑ ↓↑ (for N y = 2) where the first N y bits are the eigenvalues of Eq. (20) at j x = 0 and the other N y at j x = N A x − 1.Additionally, a 'ribbon' operator commutes with ρ A , [ρ A , V A x ] = 0, where l ∈ C indicates the links intersected by a contour C through the centers of plaquettes, from one boundary to the other in an arbitrary path, see the green dashed line in Fig. 7(c).For fixed y-BCs V A x is not independent (it is determined by fixing all sectors Eq. ( 20)), but for y-PBC its eigensectors are independent and labelled by an additional ↑/↓, so that s ∈ {↑, ↓} 2 2Ny +1 .

Classical shadows
As a warm-up, we first study again k-purities for the new model, thereby verifying correctness of the randomization, before moving on to reconstructing a classical shadow representation of the state that we will use to infer the presence or absence of TO.In Fig. 8

Detecting topological order through entanglement Hamiltonian tomography
Aiming at higher precision, we explore an alternative approach, Entanglement Hamiltonian Tomography (EHT) [28,30,32].The basic idea is to parameterize the reduced density matrix by an Entanglement Hamiltonian (EH), A (heuristic) parameterization of ground state EHs, inspired by the Bisognano-Wichmann (BW) theorem [144,145], is where H O are the local operators comprising the physical Hamiltonian, Eq. ( 17), and β O a local 'temperature', varying with the distance of O from the entanglement cut(s) (also depending on j y if translation invariance in y is broken, i.e. for fixed y-BCs).The applicability and accuracy of this ansatz was investigated for Z 2+1 2 in [56].To extract the {β O }, we follow [28] and first measure probabilities P U (b, s) in N E random bases and symmetry sector s A horizontal red line is the infinite volume limit, ϵc = 0.33 ± 0.01.Right: Symmetry-resolved Entanglement spectra ξ s,λ for ϵ = 0.075, 03, 0.5.
with N M shots each.We then minimize with respect to {β i } via classical post-processing, where ρA,s ≡ ρA,s is the EH with H s,O restricted to a symmetry sector s.Because the BW optimization is performed sector-wise, matrices of size d s are involved, versus the dimension of A, d A .In the infinite measurement limit, the optimization will yield one universal set {β O } for all s, but in practice we work with normalized ρA,s ≡ ρ A,s /p s and PU (b, s) = P U (b, s)/p s , so that the extracted {β O,s } depend on that normalization and differ from {β O }, see Appendix D for details.
Results of the BW-EHT-optimization for the Schmidt spectrum P s,λ are displayed as diamond symbols in Fig. 8(b) along with shadow results, for N E = 50 and N M = 1024.Despite comparable cost (N E • N M = 51, 200) relative to the classical shadow approach (N S = 2 16 = 65, 536), BW-EHT reproduces the eigenvalue spectrum much more accurately; values as small as 10 −6 are approximately recovered and even eigenvalues as small as 10 −11 are not far off.The inset of Fig. 8(b) shows the entanglement spectrum (ES), i.e. the spectrum of the EH, which is also well reproduced.The apparent advantage of the BW-EHT approach comes at the expense of generality, it is tailored for ground states (it can be extended to nonequilibrium states [28,110]) while classical shadows work regardless of the state.
Enabled by the performance of the BW-EHT optimization, we focus on a practical application: de-tecting topological order (TO) of quantum states.Ground states of Z 2+1 2 are separated (in the infinite volume limit) into topologically ordered, ϵ < ϵ c , and trivial states, ϵ > ϵ c , with a phase transition at a critical coupling ϵ c .Li and Haldane's entanglementboundary conjecture [56,63] asserts that TO states are 'entanglement-gapped', i.e. their ES has separated low energy (large Schmidt values) and a high energy (small Schmidt values) parts.Further, the low lying part is (up to rescaling) identical to the spectrum of a conformal field theory (CFT) describing gapless excitations at the edge of the system.We focus here on measuring the existence of an entanglement gap ∆ξ to detect TO, which has been shown as very a robust order parameter for the TO transition in this model even for very small systems [56].
Without loss of generality, to reduce finite size effects, we focus on periodic boundary conditions in y (a torus) for the BW-EHT analysis.In Fig. 8(c) we show the entanglement gap ∆ξ for N x × N y = (3 + 3) × 2, N E = 50, N M = 1024, and ℓ = 64 (black diamonds), compared to exact results (black solid line).Error bars represent the combined statistical error due to finite N E and N M , see Appendix D for details.A vertical red line indicates the infinite volume extrapolated value ϵ c = 0.33 ± 0.01 [146].We also show result for N x × N y = (3 + 3) × 3 (green dashed line) and N x × N y = (3 + 3) × 4 (blue dashed line), taken from [56] and approaching the infinite volume limit to within less than 10%.Side panels show the BW-EHTreconstructed sector-wise ES for ϵ = 0.075, 0.3, 0.5, demonstrating the closing of the entanglement gap (gray shaded area) at ϵ c where our results approximately reproduce the phase transition.
We could not numerically simulate systems larger than N x ×N y = (3+3)×2.For example, a N x ×N y = (3+3)×4 lattice of 48 qubits (20 qubits in the subsystem) exhausts our classical computational resources. 5 However, the classical (shadow-or BW-EHT-) analysis is simple for such a system if it were prepared in experiment because, while d A = 2 2N A x Ny−Ny = 2 20 = 1, 048, 576, the analysis is restricted to symmetry blocks of only d s = 2 N A x Ny−Ny = 2 8 = 256 states.This is a significant (in fact, exponential) reduction in the space over which the state is randomized, but d s still grows exponentially with the subsystem size, albeit much slower than d A .

Conclusion and Outlook
In this manuscript, we proposed randomized measurement protocols for lattice models that leverage symmetries, focusing primarily on LGT entanglement structure exploration.We devised deep-scrambling circuits that realize symmetry-conscious k-designs and illustrated their use in simple gauge and nongauge model examples.Our approach is intuitive and, therefore, easily generalizable: by examining the physical Hamiltonian one can readily identify basic symmetry-preserving interactions which can be used as the generators of a randomized measurement scheme.Consequently, if a particular physical Hamiltonian can be realized, so can our measurement scheme.
Symmetry-conscious randomized measurement schemes like those considered here have lower sampling costs compared to symmetry-ignorant schemes by avoiding randomizing over non-relevant Hilbert space parts.In particular, one obtains a sampling cost (to realize a 2-(k-)design) that scales with block size d s , instead of Hilbert space dimension d H .This reduction can be exponential, e.g., for particle number conserving systems away from half-filling, or for LGTs due to randomizing only over the physical sector of Hilbert space.In constructing k-designs, while still efficient, they also incur a somewhat larger circuit complexity.
entropies, and separated their symmetry-and distillable components.
We then focused on Z 2 in 2 + 1 dimensions where the intricate structure of gauge symmetric states can lead to topologically ordered (TO) phases.These are currently receiving great attention, including experimental realizations in AMO and solid-state platforms [134,135], motivated by applications such as fractional quantum Hall effect states [160,161] or fault-tolerant quantum computation and storage [124][125][126][127][128][129][130][131].A difficulty is that TO cannot be probed by measuring local operators, a serious impediment for its experimental verification.Our approach to overcome this is based on measuring the entanglement structure of such systems.While the importance of entanglement as a robust indicator of topological order was realized long ago [162][163][164], we developed a concrete random-measurement scheme, following the logic of Li and Haldane [63], that uses a state presentation in terms of Entanglement Hamiltonians (EH) and is based on measuring entanglement gaps of their (symmetry-resolved) spectrum using a tomographic protocol based on the Bisognano-Wichmann theorem [28].Remarkably, performing random measurements on very small subsystems as small as N A x × N y = 3 × 2, we observe a relatively sharp TO-to-trivial phase transition.While our focus was on Z 2+1 2 , the protocol can be easily generalized to other systems.
A benefit of our approach, not explicitly explored in the main text, is that symmetry-conscious randomization allows for a rudimentary, but useful, nearterm error-mitigation strategy similar to that discussed in [92].A feature of symmetry-conscious randomization is that symmetries of the input states are not lost and can be measured.Thus, machine errors that violate those symmetries are detectable after randomization, suggesting that e.g., a postselection of measurement results can improve the computation (at the cost of reduced statistics).Another potential application is approximate Haar random state preparation for thermal state algorithms [111].Finally, we expect symmetry-respecting randomized measurement schemes to be useful to investigate thermal systems, including e.g, systems with non-Abelian conserved charges [165][166][167][168][169][170][171][172][173][174][175][176].
There are many future extensions of our work.For example, while our approach significantly reduces algorithmic costs compared to a symmetry-ignorant scheme, extracting entanglement entropies and structure still relies on classical post-processing which ultimately scales exponentially with system size, an issue which can be addressed with quantum variational [30] and machine learning techniques [16,18].Extending the robust numerical analysis performed here to provide analytical performance guarantees for the circuit depth and sampling complexity of the random circuits in this work is also of clear interest.We also emphasize that realizing an approximate k-design is a sufficient, but not necessary, condition for randomized measurement protocols, and more studies regarding optimal randomization for certain observables are needed. 6Developing a formalism for shadow and entanglement tomography protocols that applies to systems with limited control or is independent of the circuit model would be useful for analog quantum simulation.Finally, while there are encouraging indications [37][38][39], the robustness of our scheme against experimental imperfections and noise should be investigated in detail.
Finally, we point out related work [12,17], following a similar idea for fermionic systems with particle number symmetry, and also demonstrating a significant cost advantage.While not programmable enough to realize k-designs, we think this is a very useful approach for fermionic entanglement tomography.We also point out Ref. [179], proposing randomized measurement schemes that take advantage of a priori knowledge about observables of interest to improve sampling complexity and Refs.[180,181], which, building on Ref. [27], propose randomized measurement schemes to extract symmetry-resolved purities.Their scheme is based on local random unitary transformations, in contrast to the symmetrypreserving unitaries we consider here.
gauge-theory dynamics on near-term quantum hardware.This work was enabled, in part, by the use of advanced computational, storage and networking infrastructure provided by the Hyak supercomputer system at the University of Washington [182].

A Equivalence of error metrics
Our choice of error metric is the average difference between matrix elements of B s , Eq. ( 4), calculated from an ensemble of unitaries generated by random circuits and those obtained from an ensemble of unitaries that form an exact 2-design, as given in Eq. ( 5).This difference is normalized by a factor of d s (d 2 s − 1) where d s is the dimension of the relevant symmetry block.Because of the otherwise prohibitive cost, in practice, for large N we compute the error by only averaging over a sample of all matrix indices of B s , where we increase the number of indices sampled until we see convergence.Mathematically, where S is the set of indices i, j, k, l, i ′ , j ′ , k ′ , l ′ sampled and the subscripts on B s denote the ensemble with which the expectation values are taken with respect to in the definition of B s .This choice of quantifying the error between our random ensemble of circuits and an exact unitary 2design was chosen for two reasons: (1) it enables a computationally tractable approach wherein we sample a subset of matrix elements until we see convergence; (2) once normalized (by a factor of d s (d 2 s − 1)), it is manifestly dimension independent once our random ensemble of circuits has converged to an approximate unitary 2-design.Despite these benefits, it is not the only reasonable choice of error metric, nor is it a standard one.Therefore, in this appendix, we demonstrate analytically that it is equivalent to more typical definitions, up to a rescaling of the approximation ratio for the case where S is the set of all index combinations.We then numerically demonstrate that sampling (the non-zero) matrix elements is also a valid (and computationally accessible) choice of error metric.
A particularly common definition of an approximate unitary k-design states that an ensemble E of unitaries forms a δ-approximate unitary k-design if and only if [183] |F where F (k) is the k-th frame potential for an ensemble of unitaries, defined as For the Haar ensemble F Haar = k! for k ≤ d s [184].Ultimately, we will numerically compare our error metric to the frame potential definition, but for the purposes of analytically showing equivalence of our error metric it is convenient to use a different formulation known to be equivalent to this one up to a rescaling of the approximation factor [185].In particular, an ensemble E of unitaries forms a δ-approximate kdesign if and only if for all balanced monomials M (U ) of degree ≤ k in the matrix elements of U where the subscript k denotes an expectation value with respect to an exact k-design.A balanced monomial of degree k in the matrix elements of U is defined as any product in the matrix elements of the form for some choice of indices.
To show the equivalence between Eq. ( 25) and Eq. ( 28) we return explicitly to the case of k = 2.
An ensemble E that forms a δ-approximate 2-design necessarily also forms a δ-approximate 1-design, so in Eq. ( 28) we can restrict our attention to balanced monomials of degree 2. Consequently, for unitaries U s acting on a symmetry sector s, and a general degree where the matrix elements are the matrix elements for an exact unitary 2-design, as defined in Eq. ( 5) and, as noted in the main text, (A s ) kl ij ≡ ⟨U s ij U s * kl ⟩.The subscript 2 denotes an average taken with respect to an exact unitary 2-design and the subscript E denotes an average taken with respect to the corresponding ensemble of unitaries.
For a δ-approximate unitary k-design the term labeled I is bounded via Eq.( 28).The terms labeled II and II * are equivalent up to the choice of indices.Recalling that for an exact 1-design (and thus also for a 2-design) and applying Eq. ( 28) for an approximate unitary 1-design, we can bound these terms as Plugging back into Eq.( 30) we obtain for any choice of indices.It is straightforward to go through the converse of this argument and show that if and |A E − A 2 | < δ ds (indices suppressed), then the approximate two-design property for ⟨M ⟩ E holds, up to a rescaling of δ.Consequently, up to an inconsequential rescaling of the approxima-tion factor, any norm on the difference of the tensors 2 ) is a valid choice for defining an approximate unitary 2-design consistent with the definition in Eq. (28).
Our error metric, if all indices are sampled, is just such a norm (namely, an element-wise 1-norm, up to dimension-dependent factors).For the sake of computational tractability, however, we randomly sample a collection of the non-zero matrix elements of B s until we see convergence.If the sampled average has converged statistically, i.e. has small error bars, this is a good indication that the sampled average accurately reflects the total average, since the total is simply a sum of sample averages.Since the matrix elements of B s 2 can take only a handful of values and are all either zero or of order d −3 s , we expect there are no "outlier" matrix elements whose errors are systematically larger than others.We have explicitly verified that the sampling scheme converges to the results obtained by evaluating all index combinations, including those where Eq. ( 5) yields zero.In all our numerics, we ensure that the number of samples taken is such that this convergence occurs.
Numerical results demonstrating the equivalence of our error metric when using finite samples of the nonzero matrix elements of B s , as given in Eq. ( 25), to the error metric in terms of the frame potential in Eq. ( 26) are shown in Fig. 9. Here, we consider a system of N = 10 sites with particle number symmetry as presented in Sect. 2 of the main text.For this comparison, we compare the two error metrics as a function of the number of unitaries N E drawn directly from CUE within each symmetry sector.In the next section, we elaborate on how we demonstrated that the random circuits we describe in the main text do, in fact, converge to sampling from the CUE for a sufficient number of layers.However, temporarily leaving this aspect aside allows us to purely compare the two error metrics.We find that the two error metrics are equivalent up to a rescaling of the approximation ratio by a factor ∼ d s for a sufficiently large ensemble of random unitaries.This is consistent with, but tighter, than our analytic results.Also consistent with our analytics, the numerics indicate sub-leading d s -dependent factors in the rescaling of the error.

B Particle number analysis
In this appendix, we provide additional details for the particle number symmetry analysis, presented in Sect. 2 of the main text.In particular, we present the numerical evidence that sufficiently deep random symmetry-respecting circuits of the sort we describe there form an approximate 2-design and provide details on the numerics leading to Fig. 3.
For both of these purposes, we compute the error ϵ as defined in Eq. (25) for ensembles E of N E random unitaries generated both by our random circuits and by drawing directly from the CUE.Note that in the latter case, the error is purely due to the fact that these ensembles have only a finite number of elements N E , as for N E → ∞ the ensemble will be a 2-design by construction.Consequently, as can be seen in Fig. 10a, the error ϵ scales as ∼ 1/ √ N E .For a sufficient number of layers, the ensembles generated by our random circuits have identical error ϵ as the ensembles drawn directly from CUE, indicating that these circuits indeed generate samples that form an approximate 2-design.This is demonstrated in Fig. 10b, where we plot the difference between the Figure 9: Ratio of our error metric ϵ [given in Eq. ( 25)] to the more standard frame potential error metric ϵ f [given in Eq. ( 26)], for an approximate unitary 2-design per symmetry sector (labeled by particle number n f ), normalized by ds, versus number of unitaries NE sampled, indicating that the two error metrics are equivalent up to a rescaling of the approximation ratio for sufficiently large NE .Results are for a system of N = 10 sites with particle-number symmetry as described in Sect. 2. For both error metrics, we sample directly from the CUE within each symmetry sector.Our error metric is computed by averaging over |S| = 2000 non-zero matrix elements.
error ϵ obtained by taking N E = 8192 samples from the random circuits and taking N E = 8192 samples directly from the CUE versus number of layers used in the circuits.Within a short depth ℓ ≈ 15, the difference reaches a floor set by N E (i.e. the standard 1/ √ N E sampling error), indicating that at this depth, drawing samples from our circuits is equivalent to drawing samples directly from CUE.All of the remaining error in approximating a 2-design is purely due to the finite number of samples N E .
While the saturation point (in terms of number of layers required) in Fig. 10(b) shows a slight dependence on the dimension d s and similar numerics for smaller system sizes shows a dependence on the number of sites N , the fact that our numerics are limited to around N ≤ 10 sites prevents us from extracting an asymptotic scaling form of these dependencies.For the pragmatic approach taken in this work-namely, demonstrating that our approach is a viable one for extracting quantities of interest for specific systems with symmetries of interest, such as Z 2 LGT, for small to moderate system sizes-these numerics are sufficient.However, a detailed analytic analysis of the scaling of the error with the number of layers for the particular models considered here (or, perhaps, generally for random symmetry-respecting circuits), remains a compelling prospect for future work.Such analyses have been done for symmetryignorant designs [186] and similar approaches should apply here [178].
B.1 Details for Figure 3(a) Fig. 10 demonstrates that ℓ = 128 is well past the number of layers needed to faithfully sample from a 2-design on each of the symmetry blocks.Thus, in this large layer limit, the scaling with N E is independent of the circuit construction, as every block simply represents a random CUE matrix in this limit.Fitting the curves in the lower right panel of Fig. 10a (ℓ = 128) and similar curves for N = 4, 6, 8 sites and, then, extrapolating to determine the number of samples N E needed to reach an error of ϵ = 0.01 leads to Fig. 3(a).
In the large layer limit, the scaling with Hilbert space dimension d s observed in this figure applies equally well to symmetry-ignorant schemes if d s is replaced with the full Hilbert space d H , i.e. the scaling observed is simply a property of sampling from any set of unitaries that forms a unitary 2-design in the N E → ∞ limit.We use this fact and the fit in Eq. ( 6) to generate the inset of Fig. 3(a) showing the relative gain r s in number of circuit samples for the symmetryconscious over a symmetry-ignorant scheme.

B.2 Details for Figure 3(b)
We now turn to providing additional details for Fig. 3(b), which shows the sample cost scaling for estimating sector-wise k-purities for subsystems of size N A = N/2 for a system in the ground state of the particle-number symmetry-preserving Hamiltonian in Eq. (7).Recall, that we consider periodic boundary conditions and systems of size N = 4, 6, 8, 10, 12, 14, 16, 18 with couplings such that ma = 0.05 (i.e., a parameter regime where the ground state is entangled).k-purities can be extracted via randomized measurement schemes by utilizing the identities in Eq. ( 8), which hold when the expectation values ⟨P U (b, s) k ⟩ taken with respect to the ensemble of random unitaries E converge to the value obtained for an ensemble that is an exact unitary k-design.
A sufficient condition for this convergence is that E forms an approximate unitary k-design; however, note that it is possible that these expectation values converge to a fixed error for a smaller number of samples N E than is needed to converge to the same fixed error in being a unitary k-design.This is indeed what is observed, as can be seen by comparing the number of samples needed for a convergent estimate of the kpurities, shown in Fig. 11(a), to the number of samples needed to reach an approximate unitary k-design, shown in Fig. 10(a).This is because an approximate unitary k-design will reproduce expectation values of all operators of degree k, whereas for this scheme, we must only reproduce the k-purities.It is also important to note that the accuracy is inherently dependent on the state under consideration; the k-design bound provides a worst case scenario.
In Fig. 11(a), we have shown the number of ensem- bles N E required to estimate the 2-and 3-purities to a relative error of ϵ = 0.05 in the infinite shot limit N M → ∞.The necessary number of ensembles peaks and then begins to decrease as a function of block dimension.This trend cannot continue indefinitely (one must always implement at least one random unitary); consequently, we expect this behavior to saturate for large enough block dimension, as the variance of the infinite shot purity estimator approaches a constant in the large Hilbert space dimension limit [33].
To create Fig. 3(b) we fix the number of ensembles N E = 1428, well beyond the the number of ensembles needed to predict the 2-purity and 3-purity to within 5 percent for all cases considered.Since N E is constant, this allows us to consider the scaling of the sample cost to purely depend on the number of shots per random unitary from the ensemble (i.e. the number of measurements N M made in each random basis).Therefore the cost N M plotted in Fig. 3(b) is representative of the full sample cost of estimating the 2-purity of the subsystem states ρ A .
As a sanity check on the N E scaling, in Fig. 11(b) we plot the true k-purities for the states used to create Fig. 3(b) and Fig. 11(a).Note that the states, even for large Hilbert space dimension, have purities of order one.As the cost in N E for purity estimation is expected to be largest in the pure state case [33], this shows that the trend in Fig. 11(a) is not simply because states at large N (large d s ) are less pure.LGT for a subsystem of size NA = 5 (N = 10).Number of indices sampled is 900.

C Details for the (1+1)d Z 2 LGT example
In the main text, we presented circuits forming symmetry-conscious k-designs for Z 2 LGT in (1+1)d with matter.We demonstrated that these circuits allow to measure k-purities and von Neumann entropies within each symmetry sector, as well as separately extracting the symmetry and distillable entanglement, and the symmetry-resolved Schmidt spectrum using classical shadows.A sufficient, but not necessary condition for such randomized measurement schemes to be successful is that the randomizing circuits form approximate unitary k-designs.In this appendix, we explicitly demonstrate that the circuits in question do, in fact, form a sector-wise approximate unitary 2-design.In particular, we show that they reproduce the correct 2-design matrix elements (see Eq. ( 5)) for sufficiently deep circuits.
Representative results are shown in Fig. 12 for a subsystem of size N A = 5 (9 qubits) of a N = 10 (matter) site system (20 qubits).As described in the main text, Gauss laws at the entanglement boundaries lead to symmetries of ρ A in the subsystem.We demonstrate that, within each sector, the random circuits described in Sect.3.1 form a 2-design, by computing the error defined in Eq. ( 25) with respect to N E random circuits with ℓ = 32 layers for all nontrivial symmetry sectors calculated. 7To compute this error, we average over |S| = 900 random non-zero matrix elements.We see good agreement between sampling from our random circuits versus sampling directly from the CUE, indicating that our circuits do indeed form approximate unitary 2-designs.LGT for a subsystem of size 3×2 (with fixed boundary conditions in y).Number of indices sampled is 900.

D Details of the (2+1)d Z 2 LGT example
In this appendix, we provide details of the analysis of Z 2+1 2 , discussed in the main text.

D.1 Approximate unitary k-designs and kpurities
The determination of k-purities follows exactly that in (1+1) spacetime dimensions.We explicitly show that the (2+1)d circuits in Fig. 7(b) explicitly realizes a 2-design by repeating the analysis of Section 2. Representative results for a subsystem of size 3 × 2 are summarized in Fig. 13 showing the error Eq. ( 25), for every symmetry sector s, demonstrating agreement between sampling from our circuits for a sufficient number of layers and sampling directly from the CUE in each sector.Further, we see convergence with increasing samples N E with the standard 1/ √ N E scaling.

D.2 Classical shadow analysis
The classical shadow analysis for Z 2+1   10 −2 −10 −3 is achieved for the largest samples (where p s • N S ⪆ 10 4 ); the BW-EHT ansatz at similar cost typically reaches a precision better than 10 −5 − 10 −6 for the same configuration.Convergence with increasing shadow number N S of the shadow-reconstructed density matrix towards the exact one is evident and shows a power-law behavior consistent with the scaling of the (1+1)d case within error bars.

D.3 Entanglement Hamiltonian tomography analysis
Finally, our Bisognano-Wichmann theorem based entanglement Hamiltonian tomography (BW-EHT) protocol follows Ref. [28], except that we perform the optimization in every symmetry-sector s separately.The approach is based on an alternative representation of a reduced density matrix with Schmidt representation ρ A = λ P λ |λ⟩⟨λ|.In particular, one defines an Entanglement Hamiltonian as Eigenvalues are given by P λ = exp(−ξ λ ), where ξ λ are the eigenvalues of H A .
Because ρ A can be split into its corresponding symmetry sectors as ρ A = s ρ A,s it follows that also H A = s H A,s , where s labels quantum numbers of the spectrum ξ λ,s of H A,s .Inserting Eq. ( 23), for the Entanglement Hamiltonian we obtain the following BW-EHT ansatz for the state within symmetry sector labeled by quantum number s: where H A,s is is a deformation of the physical Hamiltonian, i.e. local couplings β i replace the physical couplings.The ansatz is such that the state is normalized within each symmetry sector, i.e., Tr s [ρ s ] = Tr s [ρ]/p s = 1, where p s is the probability of being Here, E = {U s } is the ensemble of random circuits restricted to the block labeled by s, and P U (b, s) is the probability of measuring outcome bit string b, normalized such that b P U (b, s) = 1 for all s.P U (b, s) is determined by classically simulating an (ideal) circuit for a given number of shots.In practice, the optimization is performed using python's simplicial homology global optimization (scipy.optimize.shgo)[187] with the following parameters scipy.optimze.shgo(chi_squared,bounds,n=32, sampling_method='sobol', options=opt_dict) with sampling method 'sobol', and n = 32 sampling points in the construction of the simplicial complex, and very large bounds i.e. typically β i ∈ [−30.0,30.0]; all other options are set to their default values.
An example of the results of this analysis is shown in Fig. 15, where we show the symmetry-resolved Schmidt P s,λ and entanglement spectrum ξ λ,s of the Z 2+1 2 ground state at ϵ = 0.075, 0.3, 0.5, for N x ×N y = (3 + 3) × 2, periodic boundary conditions in y and ℓ = 64.To estimate the error from applying a finite number of random circuits and estimating the effect of shot noise on obtaining P U (b, s), we additionally perform the following analysis: We compute the exact state using exact diagonalization.We then numerically minimize, within each symmetry sector s, the relative entropy between the exact density operator σs , normalized so that Tr can be easily computed.Minimizing Eq. ( 37) with the exact same numerical optimization as used for Eq. ( 36) provides the BW-EHT result in the infinite measurement bases and infinite shot limit.The deviation of our circuit simulation from this result provides the error that we show in Fig. 8(c).

Figure 1 :
Figure 1: (a) Illustration of our proposed symmetry-conscious randomized measurement scheme, which preserves the symmetry structure of states, compared to symmetry-blind approaches.(b) Applications and benefits include cost reduction when measuring entanglement and finding classical representations of quantum states, allowing a rudimentary symmetry-based error mitigation scheme, and enabling symmetry-resolved measurement of entanglement structure.

Figure 2 :
Figure 2: Random circuit scheme for systems with particle number symmetry.Blue squares connected by black lines represent the unitary u in Eq. (3).

Figure 3 :
Figure 3: (a) Particle number-respecting random circuits NE required to approximate a 2-design with precision better than ϵ = 10 −2 in every symmetry sector (ℓ = 128) for a total number of sites N = 4, 6, 8, 10.Inset: Relative sampling cost reduction compared with a symmetry-ignorant scheme.(b) Measurements NM required to estimate sector-wise 2purities with precision better than ϵ = 0.05, for fixed NE = 1, 428 and ℓ = 128.Here, the states considered are reduced density matrices on subsystems of dimension NA = N/2 (i.e. on a bipartition of the lattice.)
5(b) shows the decomposition of the von-Neumann entropy, S = −Tr A [ρ A log(ρ A )]
Tr s [X] I s ; Tr s and I s are the sector-wise trace and identity, respectively.The bar indicates normalization, i.e. σA,s ≡ σ A,s /p s , Tr s [σ A,s ] = 1, p s ≡ Tr s [σ A,s ] is simply the number of shadows measured in one sector relative to N S .Fig. 6 compactly summarizes the results of this analysis, showing the symmetry-resolved Schmidt spectrum of ρ A in (a), comparing exact results (empty squares) versus shadows (filled triangles) for m • a = 0.05, e/m = 1, ℓ = 32 and N S = 2 17 .The eigenvalue spectrum is well reproduced down to P s,λ ≈ 10 −2 − 10 −3 ; the inset shows the full spectrum and the top of (a) the accurately recovered probability for each block, p s .Fig. 6(b) shows the relative entropy, S(ρ A,s ||σ A,s ) ≡ Tr s [ρ A,s log(ρ A,s ) − ρA,s log(σ A,s )], between exact ρA,s and shadow-reconstructed σA,s , as a function of shadows per sector p s • N S 4 .The fit (black curve) indicates approximate power law ∼ (p s N E ) − 1 2 convergence.We note that automatically fewer shadows are sampled in less important sectors.While Z 1+1 2

Figure 7 :
Figure 7: (a) Illustration of Z2 2+1 LGT, including Hilbert space and Gauss law constraints Gj.(b) Random measurement circuits: We work with an even-odd alternating layers consisting of randomly placed 'electric' rotations Rz(α), Rz(γ) and plaquette rotations U □ (β) (orange squares), approximating k-designs for sufficient circuit depth.(c) The symmetry structure of ρA originates from Gauss laws at entanglement boundaries and a non-local 'ribbon' operator V A x spanning the two entanglement cuts.
(a) we display k-purities Tr[ρ k A,s ] (k = 1, 2, 3, 4) of the Z 2+1 2 ground state at ϵ ≡ g/K = 0.1, following the approach outlined in Sect.3, as a function of N E and in the infinite shot limit (N M → ∞), with fixed y-BCs, N x × N y = (3 + 5) × 2 and ℓ = 64.As before, the k = 1 results, p s = Tr[ρ A,s ], are exact by design while k = 2, 3, 4-entropies are reproduced with increasing N E .Finally, aiming at reconstructing the entanglement spectrum from classical shadows to diagnose TO, in Fig. 8(b), we show the Schmidt spectrum P s,λ for ϵ = 0.2, comparing shadow-reconstructed (N S = 2 16 = 65, 536, triangles) versus exact results (squares).Different colors represent symmetry sectors s with weights p s ≡ Tr s [ρ A,s ] in the top panel.Large Schmidt values are well reproduced down to approximately 10 −2 , beyond which we observe significant deviations that prevent us from reliably determining the entanglement gap.

Figure 10 :Figure 11 :
Figure 10: (a) Error relative to an exact unitary 2-design ϵ (Eq.(5)) per particle number symmetry sector (labeled by n f ) versus number of unitaries NE .For clarity only particle number sectors n f ≤ 5 are shown; as they are of equivalent dimension, particle number sectors for n f > 5 have the same behavior as the particle number sector (10 − n f ).The ensemble E is sampled both from the particle number symmetry-respecting random circuits with different numbers of layers ℓ ∈ {4, 8, 16, 128} (triangles) and sampled directly from the CUE within each symmetry sector (squares).Results shown are for n = 10 sites and the error is averaged over 2000 non-zero matrix elements.Extrapolation of fits to the error as a function of NE at ℓ = 128 layers for the data pictured here and equivalent data for N = 4, 6, 8 sites is used to produce Fig.3aof the main text.(b) The difference between the error ϵ obtained by taking NE = 8192 samples from the random circuits and taking NE = 8192 samples directly from the CUE versus number of layers used in the circuits.By ℓ ≈ 15, the difference reaches a floor set by NE (i.e. the standard sampling error) with a small dependence on the dimension ds of the corresponding sector.

Figure 12 :
Figure 12: Average error ϵ in 2-design matrix elements (Eq.(5)), normalized by ds(d 2 s −1), versus number of random unitaries NE sampled from either the symmetry-respecting random circuits of Section 3.1 of depth ℓ = 32 (triangles) or from direct sampling from the CUE within each symmetry sector (squares) for Z 1+1 2