U(1) Wilson lattice gauge theories in digital quantum simulators

Lattice gauge theories describe fundamental phenomena in nature, but calculating their real-time dynamics on classical computers is notoriously difficult. In a recent publication [Nature 534, 516 (2016)], we proposed and experimentally demonstrated a digital quantum simulation of the paradigmatic Schwinger model, a U(1)-Wilson lattice gauge theory describing the interplay between fermionic matter and gauge bosons. Here, we provide a detailed theoretical analysis of the performance and the potential of this protocol. Our strategy is based on analytically integrating out the gauge bosons, which preserves exact gauge invariance but results in complicated long-range interactions between the matter fields. Trapped-ion platforms are naturally suited to implementing these interactions, allowing for an efficient quantum simulation of the model, with a number of gate operations that scales only polynomially with system size. Employing numerical simulations, we illustrate that relevant phenomena can be observed in larger experimental systems, using as an example the production of particle--antiparticle pairs after a quantum quench. We investigate theoretically the robustness of the scheme towards generic error sources, and show that near-future experiments can reach regimes where finite-size effects are insignificant. We also discuss the challenges in quantum simulating the continuum limit of the theory. Using our scheme, fundamental phenomena of lattice gauge theories can be probed using a broad set of experimentally accessible observables, including the entanglement entropy and the vacuum persistence amplitude.


Introduction
In [1], we presented an efficient scheme that allows for the quantum simulation of realtime dynamics of lattice gauge theories and reported on its experimental realisation in a system of trapped ions. The purpose of the present paper is twofold: first, we give a detailed account of the theoretical basis behind the experimental demonstration in [1]. Second, we address questions relevant to future experimental work, including a careful analysis of the effect of imperfections and a discussion of the scalability of the approach to mesoscopic system sizes. We concentrate here on trapped-ion implementations of the type described in [2,3]. This platform allows one to realize a spin chain where (i) local operations can be performed with high fidelity and (ii) an all-to-all two-body interaction can be induced between the spin degrees of freedom using so-called Mølmer-Sørensen gates [4][5][6][7]. Our scheme uses these resources in an highly efficient way and is specifically tailored to realize Wilson's formulation of gauge theories on a discrete lattice, which provides an ideal starting point to investigate the dynamics of gauge theories within a non-perturbative framework [8][9][10].
At equilibrium, and in certain parameter regimes (e.g. at zero chemical potential) quantum Monte Carlo simulations of lattice gauge theories can be carried out very efficiently [9,10]. However, non-equilibrium properties, as relevant for a variety of highenergy physics phenomena including particle-antiparticle production at high-intensity laser facilities [11,12], are not accessible, due to the fundamental sign (or complex action) problem affecting simulations in real time [13]. In the last few years, several proposals for quantum simulations of real-time dynamics of lattice gauge theories have been put forward [14][15][16], based on a variety of platforms ranging from cold atoms in optical lattices [17][18][19][20][21][22][23], to superconducting circuits [24,25] and trapped ions [26,27]. The main difficulties in implementing gauge theories in quantum simulators stem from the fact that complex many-body interactions have to be realized, while at the same time local (gauge) symmetries have to be imposed on the system dynamics [14][15][16]. To address these challenges, we proposed [1] to use encoding techniques, which exploit an analytical elimination of gauge fields. This approach allowed us to realize a digital quantum simulation scheme [28] in a system of trapped ions, that simulates the Schwinger model [29,30], which describes quantum electrodynamics in (1+1) dimensions (1 spatial dimension + time).
The Schwinger model represents a simple, yet paradigmatic example of a U(1) gauge theory. It describes the coupling of fermions to a dynamical electromagnetic field in (1+1) dimensions, and exhibits a series of phenomena, such as chiral symmetry breaking and confinement [13,31], that play a key role in the current understanding of more complex theories such as quantum chromodynamics [10]. The real-time dynamics of the Schwinger model includes the spontaneous creation of particle-antiparticle pairs [32]. Despite the simplicity of the model, such phenomena are notoriously hard to address using numerical simulations, mostly due to the absence of unbiased meth-U(1) Wilson lattice gauge theories in digital quantum simulators 3 ods ‡. The Schwinger model provides a good starting point for studies of lattice gauge theories, since interesting physical insights can be gained using moderate experimental resources due to the reduced dimensionality. In quantum simulators, the real-time dynamics can be probed using a rich set of observables including entanglement entropies and vacuum persistence amplitudes. While the study of quantum information concepts such as entanglement in the context of high-energy physics is a rather recent development [39][40][41], vacuum persistence amplitudes play an important role for the theoretical understanding of spontaneous pair creation already in Schwinger's original work [29,42].
In this article, we explain in detail how an efficient implementation of the Schwinger model can be carried out by combining the toolkit of digital quantum simulation [2,3,7,43,44] with techniques used in numerical computations on classical computers (so-called encoding techniques [45]). In this model of (1+1)d quantum electrodynamics, electrons and positrons appear in the form of fermion fields that are defined on a lattice, and which interact via electric-field gauge bosons that are defined on the links between lattice sites (see Fig. 1). For implementing the model, the fermionic operators can be conveniently mapped to Pauli spin operators, which are the natural degrees of freedom in spin-based quantum simulators [2,3,7,43]. The gauge field operators, in contrast, are associated with an infinite-dimensional Hilbert space, which poses a challenge for their implementation on a quantum simulator with bounded Hilbert space [14]. To circumvent this difficulty, several proposals have considered quantum link models [46][47][48][49], where the gauge Hilbert space is truncated and gauge fields are represented by spin-operators of finite dimensionality. Here, we follow an alternative approach and realize Wilson's formulation of lattice gauge theories [8], which takes the full infinitedimensional Hilbert space of the gauge degrees of freedom into account.
As explained in detail below, our simulation scheme relies on analytically integrating out the gauge fields [45], which is reminiscent of the analytical elimination of the fermionic fields in Monte Carlo methods [10]. In this way, we obtain an effective description in terms of a pure spin model where the gauge fields no longer appear explicitly but rather enter in the form of long-range interactions with an exotic distance dependence. Previously, this approach has been employed to aid analytical or numerical calculations [45,[50][51][52][53]. In contrast, we use this idea for quantum simulation, i.e. the realisation of the Schwinger model in its encoded form in an actual physical system. A key difficulty in realising such a simulation is the anisotropic form of the long-range interactions, as illustrated in Fig. 2. Below, we will show how these interactions can be implemented efficiently in a digital quantum simulator that features single-qubit operations and entangling gates between arbitrary pairs of spins. These resources are naturally available in trapped-ion setups, where these gate operations can be performed with high accuracy [2,3,7,43]. This platform therefore provides an ideal match for the ‡ In parameter regimes characterized by high occupancies of the photon degrees of freedom, semiclassical simulations have recently been shown to provide accurate results [33], while tensor-network methods [14,[34][35][36] have been successfully applied for short-time dynamics [37,38]. U(1) Wilson lattice gauge theories in digital quantum simulators 4 realisation of the proposed scheme.
As mentioned above, one of the key challenges for quantum simulations of gauge theories [14,15,54] is the requirement that the dynamics must obey gauge invariance, i.e. they must take place in the subspace corresponding to states that are physically allowed in the simulated model. This requirement translates into local constraints that govern the interaction between matter and gauge fields. In the case of quantum electrodynamics, the constraints are imposed by the Gauss law, which in the continuum limit is given by ∇E = ρ, where E is the electric field and ρ is the charge density. In quantum simulation proposals where both matter and gauge fields are explicitly present, gauge invariance is typically imposed by suppressing processes that take the dynamics out of the allowed subspace [14], for example by enforcing energy penalties. Hence, the resulting dynamics is only gauge invariant up to some energy scale. Our approach has the advantage that, by construction, the dynamics takes place in the subspace where gauge invariance is automatically fulfilled. Instead of introducing 2N − 1 quantum systems to simulate N particles together with the accompanying gauge fields and restricting the dynamics to a much smaller Hilbert space, in this approach we simulate the dynamics of N fermions using only N spins. The combination of the mapping to a pure spin model and its realisation by means of a digital simulation scheme allows therefore for a very efficient use of resources, which renders the quantum simulation of the Schwinger model possible with present-day experimental means.
The remainder of this paper is organized as follows. Section 2 presents the model under consideration and discusses the details of our quantum-simulation scheme. We explain the encoding strategy, and describe how the resulting highly non-local Hamiltonian can be realized efficiently in a digital quantum simulator. In section 3, we illustrate the capabilities of our approach by showing how it can be used to simulate the spontaneous generation of particle-antiparticle pairs out of the vacuum of the bare fermionic particles. Here, we explain how the resulting decay of the vacuum can be monitored by studying the vacuum persistence amplitude and how the creation of entanglement during the pair creation process can be studied. In section 4, we discuss concrete implementations with ions confined in a linear Paul trap, as has recently been reported in [1]. In particular, we analyse the effects of imperfections and discuss the scalability of the approach. In section 5, we discuss the challenges in taking the continuum limit of the theory. Finally, in section 6, we present our conclusions and an outlook.

Digital quantum simulation of the Schwinger model
In this section, we introduce the model under consideration and explain how it can be mapped to a pure spin Hamiltonian with long-range interactions by eliminating the gauge fields exactly [45] (section 2.1). Afterwards, we describe how the resulting spin model can be realized efficiently by means of a digital quantum simulation scheme U(1) Wilson lattice gauge theories in digital quantum simulators  (section 2.2).

Mapping of the Schwinger model to a spin Hamiltonian with long-range interactions
We consider the Schwinger model, which describes the interaction between electrons and positrons via electromagnetic fields in one spatial dimension. In the following, we give an overview to this model in the continuum and on a lattice following the description in [45]. To this end, we introduce the vector potential at position x with temporal and spatial componentÂ 0 (x) andÂ 1 (x). In one spatial dimension, the electric field has only one componentÊ(x) = −∂ 0Â1 (x), where ∂ 0 is the partial derivative with respect to time.Ê(x) represents the canonical momentum conjugate toÂ 1 (x) with U(1) Wilson lattice gauge theories in digital quantum simulators . Matter fields are represented by two-component spinor In the temporal gauge, the Schwinger Hamiltonian in the continuum is given bŷ where ∂ 1 is the partial derivative with respect to x, m is the fermion mass andΨ ≡ Ψ † γ 0 . In one spatial dimension, the Dirac matrices γ 0 and γ 1 are given by the Pauli operators γ 0 =σ z and γ 1 = iσ y . Using natural units = c = 1, the coupling constant g = −e is given by the charge e of the elementary particles. This model can be formulated on a lattice where points in space are separated by a distance a, while time is continuous. In the following, we are using the so-called Kogut-Susskind Hamiltonian formulation [29,42,55] of the lattice Schwinger model. In this description, the continuous fieldsÂ 1 (x) andÊ(x) are replaced by conjugate variablesθ n = −agÂ 1 (x n ) andL n = 1 gÊ (x n ). The operatorsL n ,θ n are defined on the links connecting lattice sites n and n+1 as shown in Fig. 1(a) and commute canonically [θ n ,L m ] = iδ n,m . Particles are represented by Kogut-Susskind fermions, with one-component fermion field operators defined on each site n:Φ n = √ aΨ e − (x n ) for even n andΦ n = √ aΨ † e + (x n ) for odd n. The unit cell of this staggered lattice consists therefore of two sites and the presence of an electron (positron) is indicated by an occupied even (unoccupied odd) site, as sketched in Fig. 1(b). Accordingly, the interaction of matter-and gauge fields is described by the lattice Schwinger Hamiltonian where N is the number of lattice sites and m is the fermion mass; w = 1 2a and J = g 2 a 2 , where a is the lattice constant and g the fermion-light coupling constant. Using natural units = c = 1, the parameters w, J, m, and g have the dimension of inverse length, while a and t (time) have the dimension of length. The first term in Eq. (1) describes nearest-neighbour hopping and corresponds to the creation and annihilation of electronpositron pairs §. The second and the third term represent the rest mass and the electric field energy stored in the system, respectively. The resulting dynamics is constrained by the Gauss law . In the considered lattice formulation, it takes the form of a set of local constraints as illustrated in Fig. 1(b),(c). More precisely, physical states are eigenstates of the generators of the Gauss lawĜ n =L n −L n−1 −Φ † nΦn + 1 2 [1 − (−1) n ]. We will be interested in the zero-charge subspace, whereĜ n |ψ physical = 0. Note that the electric field operators for each link take integer eigenvalues L n = 0, ±1, ±2, ... . § This is illustrated in the upper half of Fig. 1(c). The fermion fields shown at lattice sites 5 and 6 represent empty space (vac,vac). A hopping process that swaps the fermion fields of these two adjacent sites leads to configuration representing an electron-positron pair, as shown at lattices sites 3 and 4.
In the continuum (a → 0), the Gauss law is given by ∂ 1 E(x) = gΨγ 0 Ψ. In three spatial dimensions, this takes the form ∇E = ρ, where ρ is the charge density. U(1) Wilson lattice gauge theories in digital quantum simulators 7 We aim at realising the Schwinger model in a spin system. This can be achieved by two transformations (see [1], Methods section). First, the fermionic field operatorŝ Φ n can be mapped to Pauli spin operators by a Jordan-Wigner transformation [56], Second, the operatorsθ n can be eliminated by a gauge transformation [45], where each spin operatorσ − n is multiplied by a phase that depends on all gauge field operatorsθ l to its left (l < n). In this way, the Schwinger model can be expressed in terms of spin operatorsσ n representing the matter fields and electric field operatorsL n , as shown in Fig. 1(c) (lower panel). In this formulation, the Gauss law takes the formL n −L n−1 = 1 2 [σ z n + (−1) n ]. For a given spin configuration that represents a certain fermion configuration and for a given choice of background field 0 , the electric fields are completely determined, as illustrated in Fig. 1(d). More specifically, the Gauss law allows one to express the electric field operators in the formL n = 0 + 1 2 n l=1 σ z l + (−1) l , such that the gauge fields do no longer appear explicitly in the description [45]. Instead, the electric field energy termĤ E lat = J N −1 n=1L 2 n = J N −1 n=1 0 + 1 2 n l=1 σ z l + (−1) l 2 gives rise to (i) a long-range spin-spin interaction H ZZ that corresponds to the Coulomb interaction between the charged particles and (ii) local energy offsets that lead to modified effective fermion masses. For simplicity, we will assume a zero background field (the following results can be straightforwardly generalized to arbitrary background fields). The resulting Hamiltonian can be cast in the formĤ S =Ĥ ZZ +Ĥ ± +Ĥ Z , witĥ As outlined in the introduction, the main challenge for realising the Schwinger model in its encoded form in an actual physical system is the implementation of the long-range interaction HamiltonianĤ ZZ , which features an exotic asymmetric distance dependence (see Fig. 2). Note that we could have encoded the electric field operators in the form L n = 0 − 1 2 N l=n+1 σ z l + (−1) l . This alternative encoding results in a long-range spinspin coupling of the same type, but with reverse directionality (i.e. every spin interacts U(1) Wilson lattice gauge theories in digital quantum simulators with a constant strength with all spins to its right and the coupling to spins on its left decreases linearly with distance). Both descriptions are equally valid.

Quantum simulation protocol
In the following, we describe a protocol that allows one to simulate the lattice Schwinger model in a 1D spin system with long-range interactions. The protocol consists of a digital quantum simulation scheme where we incorporated ideas put forward in [57]. Due to the complicated form of the HamiltonianĤ ZZ given in Eq.
(2), a standard digital simulation approach would require N 2 time steps, where N is the number of spins. In our protocol, the total number of time steps scales linearly with N and the realisation ofĤ ZZ costs only N − 2 time steps, which is optimal ¶. Our approach is quite general, and requires only single qubit operations that can be applied to individual spins and one type of two-body interactions,Ĥ 0 = J 0 n,lσ a nσ a l , where a can correspond to any direction on the Bloch sphere and J 0 is the coupling strength. In the following, we will explain the protocol forĤ 0 = J 0 n,lσ x nσ x l . Adapting the scheme to a = x requires only minor and straightforward modifications.
The simulation protocol is based on a time-coarse graining, where the effective interaction given by Eqs. (2)-(4) is obtained in a time-averaged description, while maintaining local gauge invariance at any stage. As illustrated in Fig. 3(a), the total simulation time ¶ The long-range spin-spin interaction given in Eq. (2) is fully determined by the N×N coupling matrix shown in Fig. 2, which has a rank of N − 2 and therefore N − 2 linearly independent components. This implies that the proposed simulation scheme, using only N − 2 independent gates, is optimal.  t sim is divided into several time windows of duration T , in the spirit of the Trotter decomposition [28,58]. During each of these time windows, a full cycle of the protocol that is described below is performed. Each cycle consists of three sections as shown in Fig. 3(b). In sections I and II,Ĥ ZZ andĤ ± are simulated employing the spin-spin couplingĤ 0 . In section III, only single particle rotations are performed realisingĤ Z .

Digital quantum simulation of the Schwinger model in Wilson's formulation
Section I is divided into N − 2 smaller time windows of length ∆t I as shown in Fig. 3(c). For each of the N − 2 time windows, a subgroup of spins is decoupled from the interaction, while the remaining spins interact according toĤ 0 . In the n th time window of section I, ions 1 to n + 1 participate in the interaction, such that the HamiltonianĤ x kσ x l is implemented. Since the realisation of Eq. (2) requiresσ z kσ z l -couplings, local single-particle rotations are added in the beginning and U(1) Wilson lattice gauge theories in digital quantum simulators 10 the end of section I, rotating the x spin component into the z direction. The resulting time-evolution operator for section I, realizes the desiredĤ ZZ for one time step T , with strength J = 2 ∆t I T J 0 . For a single time step, this is exact. Trotter errors will be discussed in section 4.2, where we address imperfections of the scheme.
In section II, the part of the Schwinger Hamiltonian involving nearest-neighbour interactionsĤ ± is realized. To this end, the underlying interactionĤ 0 , needs to be modified not only in range but also regarding the type of couplings. This is accomplished by dividing section II into N − 1 elementary time slots of length ∆t II (see Fig. 3(d)). Each of these time slots is used for inducing the required type of interaction between a specific pair of neighbouring atoms. This can be done by decoupling all but the selected pair of atoms from the evolution underĤ 0 . The selected pair of atoms undergoes a sequence of gate operations that transforms theσ x nσ x l -coupling into a σ + nσ − l +σ − nσ + linteraction and consists of four steps (cf. Fig. 3(e)): (i) a single qubit operation on the two selected spins n and n + 1, U = e i π 4 (σ z n +σ z n+1 ) (ii) an evolution of the qubits n and n + 1 under the HamiltonianĤ 0 during a time ∆t II /2, e −iĤ 0 ∆t II /2 (iii) another single qubit operation U † and finally (iv) another two qubit gate operation e −iĤ 0 ∆t II /2 . These four steps result in a time evolution with Note that this equation is exact. The time evolution operator associated with the described sequence of gate operations is given by e −iĤ (n,n+1) . Repeating these four steps for sites n = 1 . . . N − 1 yields, as long as The relative strength of the nearest-neighbour HamiltonianĤ ± and theσ z nσ z l -type couplingŝ H ZZ can be adjusted by tuning the ratio of the elementary time windows ∆t II /∆t I .
In section III, the single qubit terms of the typeĤ z given in Eq. (4) are implemented. This Hamiltonian is realized in a single time window of length ∆t III . All spins are acted upon simultaneously, but each one experiences a different coupling strength. Since we assume that single qubit operations can be performed on much faster time scales than gates operating on multiple spins, we use ∆t III ∆t II , ∆t I . Together, the operations in the time sections I,II, and III approximately realize the time-evolution operator of the Schwinger model for one time step, e −iĤ S T . Trotter errors due to the finite-time coarse-graining will be discussed below in section 4, where we address imperfections of the scheme.  where J and w quantify the electric field energy and the rate at which particleantiparticle pairs are produced, and m is the fermion mass, see Eq. (1). (b) After a fast transient pair creation regime, the increased particle density favours particleantiparticle recombination inducing a decrease of ν(t). This nonequilibrium interplay of regimes with either dominating production or recombination continues over time and leads to an oscillatory behaviour of ν(t) with a slowly decaying envelope. (c) The entanglement entropy S(t) quantifies the entanglement between the left and the right half of the system, generated by the creation of particle-antiparticle pairs that are distributed across the two halves. An increasing particle mass m suppresses the generation of entanglement.

Dynamics of particle production
The proposed quantum simulation scheme allows for the experimental study of a wide range of fundamental properties in U(1)-Wilson gauge theories that are of current interest. For example, strong efforts are under way at high-intensity laser facilities such as ELI and XCELS to observe a cascade of particle-antiparticle pairs generated out of the vacuum subject to extreme electric fields [11,12], and several theoretical proposals for the quantum simulation of particle production have been put forward in recent years [14,15,[59][60][61]. The preparation of the true vacuum (the eigenstate of the Schwinger HamiltonianĤ S for finite values of J, w, m) is challenging for current experiments, as discussed in section 5 below. Therefore, we demonstrate the capabilities of our scheme by studying the coherent real-time dynamics of the creation of particleantiparticle pairs out of the bare vacuum (the eigenstate of the Schwinger Hamiltonian for m → ∞) following a quantum quench, i.e., following a rapid change from m/w = ∞ U(1) Wilson lattice gauge theories in digital quantum simulators  Figure 5. Dynamics of particle production in comparison to the vacuum persistence amplitude. The evolution of the particle number density ν(t) and the rate function λ(t) of the vacuum persistence amplitude, for fermion mass m/w = 1 and two values of the electric field, (a) J/w = 0 and (b) J/w = 1. As J/w is increased, the cost for separating particle-antiparticle pairs rises, resulting in stronger recombination dynamics. This leads to smaller absolute values of ν(t) and larger oscillations. On the lattice, the oneto-one correspondence between ν (solid line) and λ (dashed), valid in the continuum, is qualitatively retained, even for comparatively small system sizes. The included results for ν(t) for smaller system sizes (squares: N = 12, circles: N = 18) visualize that the dynamics quickly converges with increasing N .
to a finite value. In the context of particle-antiparticle production, our approach provides the potential to study various interesting quantities in quantum simulation experiments. We consider here two key observables, the vacuum persistence amplitude of the unstable vacuum (section 3.1) and the entanglement generated during pair creation (section 3.2). While the former quantity requires only local addressability (which is already needed for implementing the simulation protocol), a measurement of the entanglement entropy is more ambitious due to the increased amount of resources required for reconstructing density matrices using, e.g., quantum state tomography. This section addresses the phenomenology for the perfect implementation of the proposed simulation scheme. A detailed analysis of the influence of errors and imperfections will be given in section 4.

Decay of the unstable vacuum
Since vacuum fluctuations promote the creation of particle-antiparticle pairs, the bare vacuum |vac (i.e. the state where particles are absent) is unstable. The particle number density ν(t) created out of the bare vacuum is measured by the observable Wilson lattice gauge theories in digital quantum simulators 13 with . . . = vac| . . . |vac denoting the average with respect to the initial vacuum state. After transforming the fermonic fields to spin operators by a Jordan-Wigner transformation (see section 2.1), the particle number density is given by and can be determined through local magnetization measurements (see Fig. 4(a)). In Fig. 4(b), the quantum real-time dynamics of ν(t) is shown, illustrating the instability of the vacuum. Initially, particles are produced quickly. After a sufficiently large particle density has been generated, particle-antiparticle recombination becomes favored, inducing a decrease of ν(t). This nonequilibrium interplay of regimes with either dominating production or recombination continues over time, leading to an oscillatory behavior of ν(t) with a slowly decaying envelope. Asymptotically, the system reaches a steady state with a balance between particle production and recombination. As shown in Fig. 4(a), increasing particle masses lead to a decrease in the particle production because of the increasing energy costs for pair creation. Similarly, with larger values of J/w, the higher cost of generating a field string between electrons and positrons reduces the density of generated pairs, see Fig. 5.
An important quantity in the context of dynamics in the Schwinger model is the vacuum persistence amplitude [30] which measures the deviation from the initial state during the simulated dynamics and quantifies therefore the aforementioned decay of the unstable vacuum. In the continuum limit it has been shown that the decay of the vacuum is directly related to the particle production ν(t) = λ(t) with λ(t) = −N −1 log (|G(t)| 2 ) [30]. In Fig. 5, we show a comparison between λ(t) and ν(t) for different parameters. While on the lattice their one-to-one relation ν(t) = λ(t) is broken, we find numerically that the similarity between λ(t) and ν(t) nevertheless remains clearly visible. Vacuum persistence amplitudes are not only important for spontaneous pair creation. They appear under different names in a variety of contexts in quantum many-body theory and are therefore also of interest in other types of quantum simulations (for example, in the theory of quantum chaos [62] the associated probability L(t) = |G(t)| 2 is also known as the Loschmidt echo and quantifies the stability of quantum motion [63]. It plays also a central role in dynamical quantum phase transitions far from equilibrium [64,65]. In quantum thermodynamics, the Fourier transform of G(t) is related to work distribution functions [66], which are the basic objects appearing in nonequilibrium fluctuation theorems [67] such as the Jarzynski equality [68]).

Entanglement dynamics in spontaneous particle production
A further important quantity for the theoretical characterisation of quantum manybody dynamics is given by the entanglement entropy [69]. In the following, we study the real-time entanglement production during pair creation. We focus on the entanglement between two contiguous blocks in the spin system, which is equivalent to the entanglement between the respective blocks in the original fermionic description including the gauge degrees of freedom (see Appendix A for details). Indeed, when the system is encoded, the Hilbert space takes a factorized form (while in the initial formulation, it is not factorized due to the presence of Gauss's law); this form recovers the construction proposed in [41]. Therefore, it is possible to measure the entanglement in the original model directly in a pure spin system, which has the advantage that the partition into subsystems is always gauge invariant thereby avoiding known difficulties with tracing out parts of the system [41]. In the following, we consider the real-time dynamics of the half chain entanglement entropy, with ρ A = tr B {ρ(t)} denoting the reduced density matrix of the first N/2 lattice sites, obtained by tracing out the remaining system B. S(t) quantifies the entanglement between the left and the right half of the system, which is generated by the creation of particle-antiparticle pairs that are distributed across the two parts. As already shown in Fig. 4(c), the generation of entanglement decreases with increasing mass m, since particle creation becomes energetically costly for a large fermion mass. The dependence of the entanglement production in the electric field energy J is particularly interesting. The effective coupling between particles and antiparticles, mediated by the gauge bosons, increases with their relative spacing such that it becomes energetically unfavourable to separate particle-antiparticle pairs over large distances. This constrains the dynamics by reducing the number of particle-antiparticle pairs that can be shared between the left and the right half of the system. Accordingly, this leads to a reduction of entanglement for increasing J. In Fig. 6, we show the time evolution of the half chain entanglement entropy S(t) for two different electric field strengths J/w = 0 and J/w = 0.2 for varying system sizes. For the case of free particles, i.e., without coupling to the gauge fields (J = 0), the entanglement entropy exhibits a linear growth in time, characteristic for free fermionic theories. In the thermodynamic limit N → ∞, the linear growth would continue for all times, but for finite system sizes N < ∞ it is cut off on a time scale wt ∝ N/2. If J = 0, the entanglement growth follows initially approximately the free case up to a time scale t J = J −1 , beyond which the entanglement production is substantially slowed down. As the effective potential experienced by a particle-antiparticle pair increases linearly in the distance for nonzero J, the electric field suppresses the separation of spontaneously generated pairs over large distances. This, in consequence, reduces the amount of entanglement that can be produced. As these examples show, the quantum simulation of the Schwinger model allows for the observation of an intricate interplay between different parameter regimes, which can be studied in quantities not accessible to conventional experiments.

Imperfections of the scheme and implementation in trapped ions
In the following, we describe how the proposed quantum simulation protocol can be implemented in a system of trapped ions (see section 4.1) and discuss the effects of imperfections. There are two types of errors: (i) those that are inherent to the scheme and (ii) those that are due to experimental imperfections. The former, discussed in section 4.2, arise since our digital quantum simulation protocol realizes the desired dynamics only in a time-averaged manner, which leads to a discretization error [28,58,70,71] also known as Trotter error. The latter depend on the concrete physical implementation. In section 4.3 and section 4.4, we discuss the sensitivity of the simulations to experimental imperfections by considering the trappedion implementation that has been realized in [1]. We identify dominant errors and show that the phenomena of interest are robust against the main sources of imperfections that generically occur in this type of setup.

Implementation of the simulation scheme using trapped ions
In trapped-ion setups, spin degrees of freedom are obtained by restricting the dynamics to two (meta-) stable Zeeman levels of the internal electronic level structure of the ions [2,3]. The use of focused laser beams acting on these levels allows one to realize single-qubit operations by inducing individually addressed AC-Stark shifts (Ĥ n AC ∼σ z n ) and Rabi flops (Ĥ n RF(θ) ∼ cos(θ)σ x n + sin(θ)σ y n ). Spin-spin interactions are realized U(1) Wilson lattice gauge theories in digital quantum simulators 16 using a global laser field coupling the internal levels to external vibrational degrees of freedom. Here, we assume an implementation based on the so-called Mølmer-Sørensen interaction [4][5][6][7] that is mediated by the motional center of mass mode [2] and provides an infinite range, all-to-all two-body couplingĤ 0 = J 0 n,lσ x nσ x l , as assumed in section 2.2. This effective description of the spin-spin interaction neglects the motional degrees of freedom of the chain of ions, which is valid if the spin-motion coupling after a gate operation is negligible and the duration of the interaction is much larger than the period of the harmonic motion of the ions in the trap. This condition is well fulfilled in the considered experimental setting [2,3], as the period of the harmonic motion is typically less than 1µs and the gate duration is longer than 50µs.
Our simulation protocol requires the decoupling of individual spins from the infiniterange couplingĤ 0 = J 0 n>lσ x nσ x l (see Fig. 3(c),(d)), which can be achieved in different ways. For example, one may strongly detune individual ions from the laser fields that induce the Mølmer-Sørensen interaction by applying strong addressed AC Stark shifts [72,73]. Alternatively, one can split the ion crystal into multiple chains during a single experiment, such that only the ions that take part in the long-range interaction form a connected crystal that interacts with the laser light. This flexible scheme requires the use of micrometer-scale ion traps which increases the experimental complexity considerably [74,75]. An alternative that can be implemented in a macroscopic trap, and which has been employed in [1], consists in the transfer of the population of the idling ions to additional electronic substates that are off-resonant with respect to the laser light inducing the gate operations [3]. This decoupling (recoupling) procedure will be referred to as hiding (unhiding) below.

Discretization errors
The digital quantum simulation scheme introduced in section 2.2, allows one to realize the time evolution under the HamiltonianĤ S by means of a stroboscopic sequence, which consists of the cyclic application of Hamiltonians that can be experimentally realized H 1 ,Ĥ 2 , ...Ĥ n (see Fig. 3). The Trotter error, i.e. the difference between the desired time evolution U S = e −iĤ S t and the evolution realized by the stroboscopic sequence U sim = e −iĤ 1 t/n ...e −iĤnt/n n is bounded by [28] where represents higher-order terms + . Hence, these errors are controllable and the accuracy of the Trotter decomposition can, in principle, be increased to any desired precision by increasing the number of time steps n. However, the implementation of the proposed scheme in trapped ions poses limits on the minimum length of the step size that can be used. In the presence of decoherence, this leads to practical limitations and ||O|| max is the maximum expectation value of the operator O (see [28]).  on the accuracy with which the dynamics can be realized. More specifically, in order to suppress undesired spin-motion coupling terms during the applied Mølmer-Sørensen gates (see section 4.1), we require a minimal length ∆t min of the basic time windows that are used, such that ω trap ∆t min 1. In the following, we consider typical experimental values, where ω trap takes values on the order of MHz [2,3] and evaluate the Trotter error numerically for ∆t I , ∆t II ω −1 trap (compare Fig. 3(c),(d)). Local operations can be performed about an order of magnitude faster than entangling operations [3], and thus, one-qubit rotations are assumed to be instantaneous in our model. As Fig. 7 and Fig. 8 demonstrate, the Trotter error can be made sufficiently small to obtain a good resolution of the relevant features (see also experimental results in [1]). Thus, this error intrinsic to digital quantum simulation is well controlled and not a limiting factor.

Main experimental errors in trapped-ion implementations
In this subsection, we discuss generic errors affecting the gate fidelity of the digital quantum simulation in trapped-ion systems. In the considered type of experiment [2,3], the two main error sources impairing the required gate operations are (i) fluctuating coupling strengths of the Mølmer-Sørensen interaction and (ii) collective dephasing (see below). The first imperfection, fluctuations in the coupling strength J = (J 0 + δJ), is mainly due to intensity and beam pointing fluctuations of the laser beam, resulting in slightly modified interactionsĤ 0 = (J 0 + δJ) n,lσ x nσ are slow on the time scale of a single simulation experiment, such that δJ can be considered to be constant within one experimental run. Consequently, the experiment effectively realizes an ensemble of coherent Hamiltonian evolutions with fluctuating interaction strength. The resulting expectation values are averages over trajectories corresponding to different values of δJ.
The second major source of imperfections, dephasing between the two Zeeman levels that encode the spin states, is induced by fluctuating magnetic fields. Since the ions are typically only micrometers apart, all spins experience approximately the same field fluctuations, and the effect can be described in terms of Hamiltonian perturbations of the typeĤ deph = δω nσ z n with randomly varying coefficients δω. Again, on the time scale of a single experiment, δω can be assumed to be time-independent such that we can take these experimental imperfections into account by averaging over an ensemble of trajectories for different values of δω. It is interesting to note that the considered ideal time evolutions (starting from the bare vacuum state shown in Fig. 4a) take place in the zero magnetisation subspace, which forms a decoherence-free-subspace with respect to collective dephasing * . * In the absence of spin flips or undesired state transfers, the quantum states after a full time step are therefore invariant under this type of noise. However, as explained in section 2.2, the active qubits undergo a rotation U y (see Fig. 3c) that changes the reference frame, such that the noise after this To estimate the influence of these generic imperfections, we perform a numerical simulation of the expected time evolution for the particle number density ν(t) and the rate function λ(t). To this end, we draw δJ and δω randomly from uniform distributions. Fig. 8 shows calculated data for a chain of ten ions with J/w = m/w = 1. The solid line corresponds to the ideal case, while the circles represent the expected results taking the discussed experimental imperfections into account. Here, we assumed fluctuations of the Mølmer-Sørensen coupling rate of 5% (δJ ∈ [−0.05, +0.05]J 0 ) and a collective dephasing rate of 2.5% (δω ∈ [−0.025, +0.025]J 0 ). Assuming a Mølmer-Sørensen coupling rate of J 0 = 4kHz, the circles correspond to Trotter steps where the elementary time window T (see Fig. 3) is 0.325ms long. These fluctuations lead to a damping of the amplitude of the oscillations, though the frequency is retained. As these results show, the qualitative agreement including these errors remains satisfactory through realistic evolution times.
Additionally, if part of the qubit register is in the hiding states during a many-body interaction (see section 4.1), different phase shifts should be taken into account for ions that participate in the Mølmer-Sørensen interaction and for those that are transferred to hiding levels. Fig. 8 contains a numerical simulation including this effect (crosses), where we consider the internal levels used in [1]. These results have been obtained by adding termsĤ deph = δω nσ z n + δω lσ z l for each time window, where the first (second) sum includes ions that occupy regular (hiding) states. As Fig. 8 shows, this leads only to minor corrections. The error model above captures the dominating sources of imperfections in the considered setting. The effect of imperfect local operations on non-hidden qubits is negligible compared to the errors discussed above, since they can be performed with a much higher accuracy than multi-qubit entangling gate operations. In particular, if only a finite set of local operations is required, fidelities larger than F local = 0.99 can be reached [3].

Error detection techniques
The use of hiding/unhiding techniques (see section 4.1) generates another source of imperfections in the simulated time evolution. If the applied hiding pulses fail to transfer a quantum state from the computational basis states to the corresponding hiding states |↑ → |h ↑ , |↓ → |h ↓ , or unhiding pulses fail to transfer the quantum states back, the many-body operation will not act on the ions as intended and thus induces manybody errors that are difficult to correct for. However, as we describe in the following, postselection techniques can be employed to detect and filter out this type of error without affecting the desired unitary evolution. Such a postselection scheme can for example be realized as follows. For each step in the protocol that involves a Mølmer-Sørensen gate acting on a subset of ions, suitable hiding operations are performed, rotation is effectively given by U yĤdeph U † y = δω nσ x n , an interaction that induces spin flips. If the pulse area of a pulse deviates from the target value, an imperfect state transfer is performed |↑ → cos( π 2 + δ) |↑ + sin( π 2 + δ) |h ↑ . If the qubit is measured after a series of imperfect pules, the resulting erroneous state can be associated with a failure probability p = sin(δ) 2 for each pulse. U(1) Wilson lattice gauge theories in digital quantum simulators 20 followed by the action of the quantum gate, and the corresponding unhiding operations. Following these steps, a measurement of the population residing in the hiding levels is performed. If population is detected in the hiding levels, the experimental run is discarded. This way, only events that involve a failure of both the hiding and the unhiding pulse on a single ion (in a single hiding-unhiding step) remain undetected, which, however is strongly suppressed † †. Apart from a controllable residual error that can be suppressed to any desired accuracy, decoupling errors lead therefore only to a reduction of the rate at which simulation data can be acquired, but do not impair the desired unitary evolution.
In [1], an alternative postselection scheme has been applied to ensure that the final state fulfills Gauss' law. This involves measurements of the total magnetisation M = nσ z n at the end of simulated time evolution. Due to charge conservation in the Schwinger model,M is a conserved quantity under the dynamics induced byĤ S . The simulations performed in [1] have been carried out in the zero-charge subspace, starting from the initial state |↑↓↑↓ ... . Using the employed decoupling and measurement techniques (see [1,3]), single hiding/unhiding errors can lead to measurement results that correspond to a nonzero magnetisation of the spin system. This subset of errors could therefore be detected and filtered out by postselection.
We remark that digital quantum simulation schemes allow one in principle to use quantum error correction schemes to ensure that the desired sequence of gates is carried out with high fidelity. The postselection techniques discussed here provide an alternative route that is easier to realize experimentally and allows one to filter out the dominant errors.

Continuum limit
In the previous sections, we have shown how to simulate the Schwinger model within its lattice formulation. As we will discuss in the following, the ultimate goal, the extrapolation from the lattice to the original continuum field theory, can in principle be reached, though it remains a challenging prospect. The continuum limit a → 0, where a is the lattice spacing, involves a rescaling [45] of the couplings in the lattice Schwinger HamiltonianĤ lat given in Eq. (1), while the rest mass m is independent of a. To correctly reproduce the original Schwinger model, the thermodynamic limit N → ∞ has to be taken before taking the limit a → 0. Practically, both for theoretical or experimental data, this specific order of limits can be implemented by first fixing a lattice constant a and then extrapolating the data to † † The probability for such an event is given by P 1 = p 2 N for each hiding-unhiding step, where p is the probability for a hiding or unhiding pulse to fail, assumed to be small. The probability for such an undetectable error to occur during a time window T scales therefore with N 3 p 2 . the thermodynamic limit. In this way, data can be obtained successively for decreasing lattice spacings a, which then finally may be extrapolated to a → 0.
Taking the continuum limit requires an initial state that correctly reproduces the long-wavelength properties of the continuum theory. While the bare vacuum state chosen in the previous sections is valuable to study the dynamics of the lattice Schwinger model, it exhibits a spatial modulation at the level of the lattice spacing far beyond the long-wavelength limit. In the following, we will therefore choose a different initial state that lies in the long-wavelength sector while still being accessible experimentally. More specifically, we start from the ground state for the parameters m > 0, w > 0 and g = 0, and then quench to g/m > 0. The system may be initialized in the desired state via adiabatic state preparation. For that purpose the system is initialized in a Néel-type state |↑↓↑↓ . . . , which can be prepared easily as well as with high fidelity [3] and which is also the ground state of the HamiltonianĤ m = (m/2) n (−1) iσz n . Afterwards, the Hamiltonian is adiabatically transformed fromĤ m toĤ =Ĥ ± + (m/2) n (−1) iσz n by engineering a time-dependent HamiltonianĤ(t) =Ĥ m + f (t)Ĥ ± with f (0) = 0 and f (t ) = 1 where t denotes the time of the end of the process. For a sufficiently slow transformation the state always remains in the ground state manifold of the instantaneous Hamiltonian according to the adiabatic theorem [76]. The error in following the adiabatic passage is mainly set by the minimal gap of the Hamiltonian H(t) during the sweep. Importantly, the considered transformation will not encounter a quantum phase transition and therefore no gap closing. As a consequence, the U(1) Wilson lattice gauge theories in digital quantum simulators considered adiabatic state preparation is very efficient. Instead of conventional adiabatic state preparation, which is based on a continuous deformation of the Hamiltonian over time, this procedure can also be realized in a digital quantum simulation device, as demonstrated recently with superconducting qubits [77]. When choosing the digital approach also for the preparation of the initial state, the overall runtime of the quantum simulation, i.e., the total number of gate operations, increases respectively. Fortunately, the number of gate operations per Trotter step for the preparation is N , which is smaller than for the digital simulation for the full Schwinger model which is 2N − 2.
In the following, we will consider the continuum limit exemplarily for the vacuum persistence probability L(t) = |G(t)| 2 with G(t) = ψ 0 |e −iĤst |ψ 0 . In the limit a → 0, the rate function has a well-defined limit, where L = aN is the total length of the system. In Fig. 9, we show the time evolution of κ(t) that is induced by switching-on of the electric-field coupling for different system sizes. The chosen parameters correspond to m/g = 1 and m/w = 1. As mentioned before, establishing the continuum limit requires a successive extrapolation to the thermodynamic limit for decreasing values of the lattice spacing a. In Fig. 9(b), we perform this extrapolation for a few specific time points. For the system sizes available numerically, we find that the leading order 1/N correction is not always sufficient for this purpose, and that also the next-to-leading order contribution 1/N 2 U(1) Wilson lattice gauge theories in digital quantum simulators 23 has to be considered to yield good fits for all data points. As next step towards the continuum limit, we decrease the lattice spacing by a factor of 2, yielding m/w = 0.5. The corresponding data is shown in Fig. 10. Compared to the previous case, finitesize effects become stronger, especially for larger times. For mt 3, the extrapolation to N → ∞ can be performed, as shown in the right panel of Fig. 10. Nevertheless, as becomes clear from these data, a fully reliable limit N → ∞ requires even larger systems, which is beyond the capabilities of the utilized numerics. These results show that it is in principle possible to obtain the continuum limit. The main long-term challenges in this context are the preparation of the initial state as well as the need for large system sizes, which becomes more severe for decreasing lattice constant a.

Conclusions and outlook
In this article we described and analysed a protocol for the digital quantum simulation of the Schwinger model, a U(1)-Wilson lattice gauge theory, that has recently been proposed and experimentally demonstrated in [1]. Our scheme is based on encoding the gauge degrees of freedom in a spin chain with long-range interactions, and thus exploits a strategy that has been used previously in numerical calculations. We have shown how the encoded Hamiltonian with its complex long-range interactions can be realized using resources available in state-of-the-art digital quantum simulators, requiring only addressable single-qubit manipulations and one type of two-qubit gate. By construction, the protocol retains exact gauge invariance of the dynamics. Furthermore, it realizes a quantum simulation of 2N − 1 degrees of freedom (N matter fields and N − 1 gauge fields) using only N physical qubits. As we have shown, the number of gate operations scales only linearly with system size, thus allowing for an efficient scaling of such experiments to larger chain lengths. We have performed a careful analysis of potential error sources. Through numerical finite-size scaling analyses, we found that already mesoscopic quantum simulators can reliably reproduce the real-time dynamics relevant for larger systems.
Such quantum simulations open exciting prospects, as they are able to extract observables that are not accessible to conventional experiments, such as the vacuum persistence amplitude and entanglement entropy. Thus, we can expect to obtain new insights, e.g., into the propagation of entanglement in out-of-equlibrium dynamics. In systems with short-range interactions, the Lieb-Robinson theorem restricts correlations to spread only within sound cones [78]. This restriction does not apply for long-range interactions [79], and different dynamical regimes can be observed, for example if interactions follow a power law [80][81][82]. It will be interesting to study how the possible dynamics is affected by the exotic type of long-range interactions that govern gauge theories. Vacuum persistence amplitudes appear in a wide range of contexts in quantum many-body theory ranging from quantum chaos [62] to dynamical quantum phase transitions far from equilibrium [64]. Therefore it is an interesting question to which extent further insights into the dynamics of lattice gauge theories can be gained by exploring the connections to these other concepts. Finally, in view of the long-term goal to simulate lattice gauge theories using controlled quantum systems, it will be very valuable to explore avenues to generalize the proposed scheme to two spatial dimensions or ladder geometries, and to non-Abelian gauge theories.

Appendix A. Entanglement in the encoded Schwinger model
We consider the Schwinger model in its encoded version, as described in section 2.1, and are interested in evaluating the entanglement between two adjacent regions in space. As explained below, the entanglement between two contiguous blocks of the reduced spin system (i.e. in the encoded model where the gauge fields have been analytically eliminated) is identical to the half-chain entropy in the original model (which includes both matter-and gauge fields). This is a consequence of the fact that the gauge degrees of freedom are fixed for a certain spin configuration and background field in combination with the requirement that only states within a fixed charge sector are considered.
As explained in detail in section 2.1, the encoding of the Schwinger model in a pure-spin Hamiltonian involves a gauge transformation that eliminates the gauge field operators which represent the vector potentials, as well as a Jordan-Wigner transformation that maps fermionic degrees of freedom to spin operators. Neither of these transformations has an effect on the correlations between the left and right half of the system (this can be understood by noting that either transformation acts on the left subsystem only with operators that are entirely contained within this part of the system). In the following, we discuss therefore the last step of the encoding, which entails the elimination of the gauge field operators that represent the electric fields in the continuum limit. We consider the situation illustrated in Fig. A1a, which involves U(1) Wilson lattice gauge theories in digital quantum simulators  Figure A1. Illustration of a quantum state that respects Gauss' law. As explained in the text, we introduce a bipartiton that is defined by a cut dividing the system into a left part and a right part. For analyzing the entanglement between these parts, we consider superpositions of states of the shown type, where the left (right) part of the system is described by |Ψ L/R = |Ψ spins L/R |Ψ E-fields L/R . spin operatorsσ n at lattice sites n and electric field operatorsL n that are defined on the links between two lattice sites and take integer eigenvalues L n = 0, ±1, ±2, .... We assume that the background electric field takes the value 0 = 0 (as in section 2.1) and consider a bipartition that is defined by a cut as shown in Fig. A1b. The Gauss lawL n −L n−1 = 1 2 (σ z n + (−1) n ) can be used to determine the values the electric field stepwise from left to right (compare Fig. 1, panels b and c), such that the state of the full quantum state of the system can be written in the form with coefficients λ ij . In this notation, the subscript L (R) refers to the left (right) part of the chain with respect to the chosen cut. For each side, the quantum state can be written as a tensor product of a spin state and a state that refers to the gauge fields. Note that the state referring to the gauge fields of the right part of the system depends on the spin-(and corresponding gauge field-) state of the left side, as indicated by the subscripts of the last ket. However, due to charge conservation, the Gauss law can also be used to determine the electric fields starting from the right end of the chain, where the background field also takes the value 0 = 0. Therefore, the full quantum state of the system can be expressed in the form This is a consequence of the fact that the electric fields of the right side can be inferred from the background field and spin configuration on the right side alone. The gauge field at the cut can be inferred by measuring the spins on either the right or the left side separately. In other words, the gauge field at the link bears no additional information once the spin state on either of the two sides is known. This applies also to the gauge field at the cut. Hence, the entanglement between the two parts of the system can be calculated using the reduced spin state

Entanglement in the lattice
It does not matter where the cut is placed with respect to the position to the nearest gauge field. The entanglement can be evaluated by including the gauge field nearest to the cut to either side. Similarly, methods that are based on doubling the link and the corresponding gauge field and assigning one copy to either part of the chain [41] lead to the same result. This argument can be straightforwardly generalized to mixed states.
We note that the entanglement with respect to the original model is not given by the entanglement in the reduced spin system if the requirement of charge conservation is not fulfilled. In this case, there exist for example separable reduced spin states which correspond to entangled states of the full system (involving both spin-and gauge fields). However, this is not a concern as long as charge conservation can be guaranteed to be respected by the considered time evolution. In the Schwinger model, changes in the fermion number are only possible through charge conserving particle-antiparticle creation or annihilation events. This is also the case for the discretized dynamics realized in our simulation scheme (Trotter errors do not violate charge conservation). Implementation errors such as spin flips (which correspond to the creation or annihilation of a single charge) can lead to states that do not respect charge conservation. In this case, the total magnetisation of the spin system is nonzero, and the corresponding states can be detected and filtered out by postselection (see Sec. 4.4).