Soft photon radiation and entanglement

We study the entanglement between soft and hard particles produced in generic scattering processes in QED. The reduced density matrix for the hard particles, obtained via tracing over the entire spectrum of soft photons, is shown to have a large eigenvalue, which governs the behavior of the Renyi entropies and of the non-analytic part of the entanglement entropy at low orders in perturbation theory. The leading perturbative entanglement entropy is logarithmically IR divergent. The coefficient of the IR divergence exhibits certain universality properties, irrespectively of the dressing of the asymptotic charged particles and the detailed properties of the initial state. In a certain kinematical limit, the coefficient is proportional to the cusp anomalous dimension in QED. For Fock basis computations associated with two-electron scattering, we derive an exact expression for the large eigenvalue of the density matrix in terms of hard scattering amplitudes, which is valid at any finite order in perturbation theory. As a result, the IR logarithmic divergences appearing in the expressions for the Renyi and entanglement entropies persist at any finite order of the perturbative expansion. To all orders, however, the IR logarithmic divergences exponentiate, rendering the large eigenvalue of the density matrix IR finite. The all-orders Renyi entropies (per unit time, per particle flux), which are shown to be proportional to the total inclusive cross-section in the initial state, are also free of IR divergences. The entanglement entropy, on the other hand, retains non-analytic, logarithmic behavior with respect to the size of the box (which provides the IR cutoff) even to all orders in perturbation theory.


Introduction
The emission of soft radiation is ubiquitous in scattering processes in gauge theories and gravity. The hard asymptotic particles are accompanied by an infinite number of soft photons and gravitons. Much of the structure of the IR dynamics is dictated by symmetry, independently of the details of the UV completion of the theory [1,2,3,4,5,6,7,8,9,10]. It is important to provide measures of the entanglement between the hard and soft degrees of freedom [11,12,13,14,15,16,17], and to uncover implications for the structure of the S-matrix (in both gauge theories and gravity) [18,19,20,21,22], black hole physics [23,24,25,26,27] and holography [28,29,30,31].
In [17] we computed the leading perturbative entanglement entropy between the soft and the hard particles produced in two-electron scattering processes in QED. The initial and final states were dressed with clouds of soft photons, in accordance with the Faddeev-Kulish construction [32,33,34,35]. We found that tracing over the entire spectrum of soft photons, including those in the clouds, leads to (logarithmic) IR divergences in the perturbative expansion of the Renyi and the entanglement entropies. The coefficient of the logarithmic divergence was found to be independent of the dressing and to contain physical information. In a certain kinematical limit, it is related to the cusp anomalous dimension in QED [36].
A question that arises concerns the order of limits. We place the system in a large box of size L, imposing, therefore, an infrared cutoff λ of order 1/L. For the perturbative calculations, we first expand to a given order in perturbation theory, keeping λ finite, and take the continuum, λ → 0 limit at the end. The logarithmically divergent quantities in this limit are the entanglement and Renyi entropies per unit time per particle flux.
Thus, because of the IR divergences, the perturbative expansion of the entanglement and Renyi entropies breaks down in the strict λ → 0 limit. Rather, its validity is guaranteed at small, but fixed λ, taking the electron charge e to be sufficiently small. Equivalently, we take the size of the box to be large but finite. 1 These results suggest that the continuum, λ → 0 limit does not commute with the perturbative expansion. An analogy to keep in mind concerns the perturbative expansion of the conventional Fock basis scattering amplitudes. At any finite order in perturbation theory, these amplitudes are nonvanishing and exhibit IR divergences due to virtual soft photons. To all orders, however, these IR divergences exponentiate leading to the vanishing of these amplitudes [37,10]. This result can also be understood as a consequence of symmetries associated with large gauge transformations [1,2].
The structure of the reduced density matrix is indeed very suggestive. At any finite order in perturbation theory, the matrix elements that are off-diagonal with respect to momentum indices, are non-zero, containing IR logarithmic divergences in λ (for both the dressed and undressed cases). To all orders in the electron charge however, the logarithmic IR divergences exponentiate, leading to the vanishing of these off-diagonal elements in the continuum limit [12,13,14]. The diagonal elements are given in terms of inclusive Bloch-Nordsieck rates associated with dressed box states. These are free of any IR divergences in λ, order by order in perturbation theory [38,37,10]. However they scale inversely proportional with powers of the volume of the box. They also tend to zero in the continuum limit, albeit less fast than the off-diagonal elements in the momentum.
So, to all orders in perturbation theory, the reduced density matrix assumes an almost diagonal form in the continuum limit -there are nonvanishing off-diagonal elements with respect to particle polarization indices. The number of nonvanishing eigenvalues grows large, revealing strong entanglement between the soft and hard degrees of freedom and decoherence [11,12,13,14,15,16]. See also [39,40] for earlier work.
As we will see, there is an exponentiation pattern for the IR divergences appearing in the perturbative expansion of the Renyi entropies. As a result, to all orders in the electron charge, the Renyi entropies are 1 It was suggested in [15,16,17] to distinguish between the soft cloud photons, "dressing" the asymptotic charged particles, and additional radiative photons, comprising the soft part of the emitted radiation. The extra soft radiated photons have energies greater than the energy characterizing the soft cloud photons. Indeed, the dressed amplitudes describing the emission of photons with energy less than or equal to the energy of the photons present in the clouds are suppressed. Restricting the partial trace over the additional soft radiated photons leads to a well-behaved perturbative expansion for the Renyi and entanglement entropies. In particular, these are finite, free of IR divergences, order by order in perturbation theory [17]. We will not study this second type of partial tracing in this work. free of logarithmic IR divergences, rendering the continuum limit well defined. The entanglement entropy, though, has non-trivial volume dependence. It is non-analytic in terms of dimensionless parameters given by products of the size of the box L and various energy scales characterizing the scattering process.
Below we outline the plan and summarize the main results of this paper. We first generalize the perturbative computation of [17] to generic scattering processes, and study the structure of the leading entanglement entropy. We consider an initial Faddeev-Kulish state with arbitrary (but relatively small) numbers of electrons and positrons, and trace over the entire spectrum of soft photons in the final state. The reduced density matrix has one large eigenvalue, governing the behavior of the Renyi entropies for integer m ≥ 2 and the entanglement entropy at leading order. The coefficient of the IR logarithmic divergence of the leading entanglement entropy exhibits certain universality properties. In particular, it is independent of the dressing of the incoming charged particles with soft photons. The dominant contributions arise from amplitudes associated with two-particle interactions, at small scattering angle or scattering angle close to π. In the relativistic, high energy limit, these dominant contributions are proportional to the cusp anomalous dimension in QED, extending the result of [17] to more general scattering processes. We also find that the entanglement entropy per unit time per particle number density scales proportionally with the number of charged particles in the initial state. These results are described in great detail in Section 2.
We proceed in Section 3 to compute the next to leading order corrections to the Renyi and entanglement entropies. Since the leading order analysis reveals that the singular part of the entanglement entropy is independent of the dressing function, we focus in a Fock basis computation concerning electron-electron scattering for simplicity. At next to leading order, two types of IR divergent terms appear in the expressions for Renyi entropies. The first is proportional to the sum of the squares of the amplitudes for the emission of two soft photons, and it is of order (ln λ) 2 . The second is proportional to the sum of the squares of the amplitudes for the emission of one soft and one hard photon, and it is of order ln λ. A similar term to the latter one appears at leading order, but the coefficient is controlled by the sum of the squares of the amplitudes for single soft photon emission [17]. We comment on the pattern of the exponentiation of these logarithmic IR divergences. In addition, we obtain the next to leading order corrections to the coefficient of the non-analytic part of the entanglement entropy.
The analysis of the two-electron scattering case to all orders in perturbation theory is carried out in Section 4. We derive an exact expression for the large eigenvalue of the density matrix in terms of a series of hard scattering amplitudes. At any finite order in perturbation theory, this eigenvalue is plagued with logarithmic IR divergences due to virtual soft photons. To all orders, however, the IR divergences exponentiate leading to a finite, nonvanishing result. The shift of the large eigenvalue from unity becomes arbitrarily small in the continuum, large volume limit, determining completely the behavior of the Renyi entropies to all orders. The contributions from the small nonzero eigenvalues (to the Renyi entropies) remain subleading, despite the fact that their number grows with the volume of the box. Thus the Renyi entropies to all orders are free from any IR divergences in the continuum limit. The Renyi entropies per unit time per particle flux are proportional to the total, inclusive cross-section in the initial sate. We comment on the dependence of the total cross-section on the reference energy scale separating the soft from the hard parts of the Hilbert space, as well as other characteristic energy scales associated with the scattering process. The entanglement entropy, on the other hand, retains the non-analytic, logarithmic behavior with respect to the size of the box L, even to all orders in the continuum limit.
Throughout, we follow the notations and conventions of [17] -see also Appendix A. The following leading soft theorems [10,37,38,41,42,43,44,45] concerning soft photon emission are applied lim | q|→0 and lim | k|→0 A multielectron/positron dressed state α = {e i , p i , s i } is given as a product of a conventional Fock basis state |α and a coherent state |f α describing the photons in the cloud The dressing function is given by [32,33] In particular, the dressing function f µ α ( k, p) is singular as the photon momentum k vanishes. Also, t 0 is a time reference scale and c µ is a null vector, c 2 = 0, satisfying ck = 1. Only soft photons, with energies below the infrared reference scale E d , are present in the cloud. Taking E d < 1/t 0 , we may approximate the phase e −ipk t0/p 0 in Eq. (1.5) with unity -see [16] for comprehensive discussions. S-matrix elements between asymptotic dressed states are free of IR logarithmic divergences, order by order in perturbation theory [32,33,34,35,46,47,48]. The dressing energy scale can be taken to be sufficiently small so as the leading soft theorems can be readily applied in various computations of amplidutes.
In [17] we computed the overlap between coherent photon states associated with generic charged states α = {e i , p i , s i } and β = {e i , p i , s i }. See also [16]. Calling the β particles outgoing and the α particles incoming, we define η i to be +1 for all outgoing particles and −1 for all incoming particles. Then when the net charges of the incoming and outgoing states are equal, Q α = Q β , and to all orders in the electron charge e, we find is the magnitude of the relative velocity of particle j with respect to i. The overlap admits the following expansion Notice that when the momenta of the states β and α differ, B βα is nonzero and positive [10]. Therefore, to all orders in the electron charge, the overlap vanishes in the λ → 0 limit. As we already remarked, the infrared cutoff λ is taken naturally to be the inverse of the size of the box: λ = 2π/L. The particle momenta take discrete values and the corresponding single particle states are unit normalized.
Using a reference energy scale E, which we take to be smaller than the mass of the electron, we separate the asymptotic Hilbert spaces into soft and hard factors where H H includes particle states with energy greater than E and H S includes photon states with total energy less than E. We take the energy scale E d , characterizing the cloud photons to be arbitrarily small, so that we can apply soft theorems to simplify various expressions: λ < E d < E. We take the scale Λ characterizing soft virtual photons to be of the order of E d (or E for Fock basis calculations

Scattering and entanglement
In this section we extend the perturbative analysis of [17] to generic scattering setups. We will work at leading order in perturbation theory. Higher order corrections are discussed in the following sections.
Consider an incoming Faddeev-Kulish state α with m electrons and n positrons. The corresponding bare state |α is taken to be a momentum eigenstate. Then the incoming state factorizes in H H × H S as follows Initially, there is no entanglement between the soft and hard degrees of freedom since |α d is a product state. Entanglement arises as a result of scattering. In order to ensure the applicability of the perturbative analysis, we keep the total number of charged particles in the initial state to be relatively small. The outgoing state is given by where we denote by iT the non-trivial part of the S-matrix. Since the S-matrix is unitary, the following relation holds Inserting a complete basis of dressed states, the outgoing state can be written as where T βα = d β|iT |α d and T βγ,α = d βγ|iT |α d are S-matrix elements between dressed states. Summation over the final state indices β and γ is implied. The indices stand for both the momenta and the polarizations of the particles. The set {|β } includes all allowable final states with no photons. In addition, it includes final states with the smallest number of photons that can be produced as a result of a given number of electron/positron annihilations. These photons are hard. 2 For example, if l annihilations take place, the corresponding state contains m − l electrons, n − l positrons and 2l hard photons. Notice that 0 ≤ l ≤ min(m, n). The state |βγ contains an additional photon, which can be either hard or soft, emitted from a charged particle.
The undressed amplitudes T βα and T βγ,α 3 are given in terms of the conventional connected and disconnected Feynman diagrams. The leading contributions to T βα arise from disconnected diagrams (or connected in some special cases) of order e 2 if no photons are present in |β , and of order e 2l if l pairs of photons are present. Accordingly, the leading contribution to T βγ,α is of order e 3 or e 2l+1 . We denote the order of T βα with e N (where N depends on the details of β).
We will also use the following relations between dressed and undressed amplitudes [33,16,17] and As shown in [47], the function F βα ( q γ , r ( q γ )) is smooth and nonsingular in the limits λ, ω γ = | q γ | → 0 (and of order the dressing scale E d ). The volume factor appears due to the relative normalization between box and continuum states (see Section 2.3). The ellipses in Eq. (2.4) stand for higher order contributions, which do not contribute to the leading Renyi and entanglement entropies. The outgoing density matrix (which is pure) is |Ψ out Ψ| out .
To obtain the hard density matrix ρ H , we trace over the entire spectrum of soft degrees of freedom. These include cloud photons with energies less than the dressing scale E d , and additional radiative photons with energy less than E. 4 Therefore, the reduced density matrix is given by Using Eq. (2.4) and implementing the partial trace, we can expand ρ H in terms of ket-bra operators associated with conventional, Fock-basis hard states. The terms relevant for the leading order computation of the Renyi and the entanglement entropies are derived explicitly in Appendix B. A useful list of formulae to calculate the various partial traces is included in Appendix A. Let us stress that for our perturbative calculations, we expand to a given order in perturbation theory, keeping the volume of the box finite, and take the continuum λ → 0 limit at the end. The matrix elements of ρ H are given in terms of dressed amplitudes, which are free of IR divergences in λ (at least perturbatively), as well as overlaps of coherent photon states describing the clouds -see Eq. (B.1). The diagonal elements are proportional to inclusive Bloch -Nordsieck rates associated with dressed box states. They are free of any IR divergences in λ, order by order in perturbation theory [17]. The off diagonal elements are proportional to the overlaps f β |f β , which at any finite order in perturbation theory, induce logarithmic divergences in λ -see Eq. (1.9). Finally, the unitarity relation Eq. (2.3) ensures that Trρ H = Trρ = 1.

The large eigenvalue of ρ H and the leading order Renyi entropies
In Appendix B, we compute the Renyi entropies for integer m ≥ 2 perturbatively, extending the analysis of [17] to the generic case. The leading order Renyi entropies are of order e 6 , irrespectively of the number of charged particles in the initial state. This is because disconnected diagrams contribute to leading order. Either a pair of electrons or positrons scatter and the rest propagate intact, or an electron/positron pair gets annihilated to two photons (and the rest of the charged particles propagate intact). Indeed to leading order, we find for m ≥ 2 is an order e 6 quantity, depending on the amplitude for single real photon emission in the energy range λ < ω γ < E. Therefore, the Renyi entropies for integer m ≥ 2 take the form These results are in accordance with the fact that ρ H has one large eigenvalue, which is equal to 1 − ∆ to order e 6 . All other nonvanishing eigenvalues are of order e 6 (or higher). Their sum must be equal to ∆. To prove this, let us set as dictated by the first two lines of Eq. (B.1), and write the reduced density matrix in the form where G is an order e 6 matrix. 5 The coefficients C β and C βγ are given in terms of dressed amplitudes and overlaps of coherent cloud photon states in Eq. (B.2) and Eq. (B.3), respectively. Specifically we obtain where the ellipses stand for higher order terms. Using this expression (as well as energy conservation), we find that G annihilates |Φ to order e 6 . 6 Thus |Φ is an eigenstate of ρ H with eigenvalue The last equality is shown explicitly in Appendix B. The other nonvanishing eigenvalues coincide with the nonvanishing eigenvalues of the matrix G. 7 Each is at least of order e 6 . Their sum is equal to ∆, as can be verified by computing explicitly the trace of G to order e 6 , Eq. (B.30) and Eq. (B.31).

Entanglement entropy to leading order
The entanglement entropy is given by the expression where the sum runs over the eigenvalues of the reduced density matrix. As we have seen, to leading perturbative order (e 6 ), there is a "large" eigenvalue, λ Φ = 1 − ∆, and the rest non-zero eigenvalues add up to ∆. We set λ i = e 6 a i when i = Φ, where a i are order one quantities, to obtain The leading entanglement entropy is of order e 6 ln e 6 . The coefficient of the non-analytic part is obtained to be equal to ∆. We will consider the dominant singular part in the limit λ → 0. We get Using soft theorems, we obtain for the singular part of the leading entanglement entropy e s e s η s η s p s p s (p s q γ )(p s q γ ) (2.19) where the bare amplitude T βα is calculated at tree level. The sum over β is now restricted to final states for which the leading non-trivial (disconnected) diagrams are of order e 2 . This formula is universal, applicable for generic initial states. It is also independent of the dressing.

Continuum limit
To take the continuum limit, we need the relative normalization between continuum states and box states Up to negligible terms, of order E 2 d (at most), which we drop. 7 The corresponding eigenstates must be orthogonal to |Φ .
The relevant final states {|β } comprise configurations of m electrons and n positrons, arising from the scattering of any two charged particles in the initial state, as well as configations of m − 1 electrons, n − 1 positrons and 2 photons, resulting from the annihilation of an electron/positron pair. In either case the number of particles in the final state is N = m + n. We denote the corresponding invariant amplitudes by iM(α → me − + ne + ) and iM(α → (m − 1)e − + (n − 1)e + + 2γ), respectively. We shall consider generic cases where the momenta of the initial particles differ from each other. For simplicity we average over the polarizations of the initial particles.
The singular part of the leading entanglement entropy takes the following form in the continuum limit The integral over the soft photon momentum can be carried out using the relation yielding a logarithmically divergent factor in the IR cutoff λ. The kinematical factor B βα is positive, given in terms of the momenta of the initial and final charged particles by Eq. (1.7). Therefore, the leading entanglement entropy takes the form This expression generalizes the two-electron scattering case studied in [17]. Notice also that it is Lorentz invariant. We would like to understand the scaling with the number of charged particles N in the initial state.
In the generic case, non-trivial disconnected diagrams in which only two particles interact with each other contribute at leading order. In Appendix C we explain in detail how to obtain the dominant leading contributions to the squares of the amplitudes iM(α → me − +ne + ) and iM(α → (m−1)e − +(n−1)e + +2γ) in the continuum, large volume limit. We also implement the integrations over the momenta of the final particles in Eq. (2.23).
Disconnected Feynman diagrams for a three particle initial state that contribute to the entanglement entropy at leading order are shown in figure 1. Detailed diagrams and explicit expressions are included in Appendix C. We end up with the following expression for the leading entanglement entropy Here, T is the time scale of the experiment (arising from an overall energy conserving δ-function). The sums run over all possible pairs of interacting incoming particles. We may always choose to work in the center of mass frame of the two interacting particles, where p i = − p j . The scattering angle is denoted by θ and the center of mass energy by E ij = 2E i . The magnitude of the relative velocity is given by is the two electron result of [17], averaged over all possible interacting pair configurations.
The entanglement entropy per unit time per particle density is equal to and diverges logarithmically in the continuum, λ → 0 limit. Notice the scaling with the number of particles in the initial state, implying that the entanglement entropy is further enhanced as compared to the two-electron scattering case.
Of course the appearance of logarithmic divergences signals the breakdown of perturbation theory. Rather the perturbative analysis is valid when the size of the box is kept large, but fixed, taking the coupling e to be sufficiently small. As shown in [17], upon restricting the partial trace to be over the additional radiated soft photons, which energies in the range E d < ω γ < E, and thus excluding the soft photons in the clouds, regulates the IR logarithmic divergences, and leads to well-behaved results in the continuum limit. This type of tracing alleviates the decoherence of the reduced density matrix [15,16].

Relation to the cusp anomalous dimension
The integrand in Eq. (2.24) diverges for forward (θ = 0) and backward (θ = π) scattering. Since forward or backward scattering cannot be distinguished from no scattering, we introduce an effective lower and upper cutoff on the scattering angle, θ 0 ≤ θ ≤ π − θ 0 where θ 0 is a small angle, to regularize the integral. As a result, the dominant contributions to the leading entanglement entropy arise from the regions θ θ 0 and θ π − θ 0 , close to the cutoffs.
These dominant contributions exhibit certain universal behavior, irrespectively of the details of the initial state. To see this, recall that the squares of the amplitudes for the processes e − e − → e − e − , e + e + → e + e + and e + e − → e + e − in the center of mass frame are given, respectively, by where E and p = | p| are the energy and momentum of each one of the incoming particles, while θ is the scattering angle. As before, we have averaged over the polarizations of the initial particles for simplicity. The above three amplitudes diverge for forward and backward scattering. For the process e + e − → γγ, we obtain which remains finite and well-behaved at θ = 0 and θ = π. As usual, it is useful to define the Mandelstam variables We will examine the behavior of the dominant terms at θ 0 and π − θ 0 in a particular double scaling limit. We take the momentum to be very large and the cutoff angle very small, p → ∞ and θ 0 → 0, keeping the quantity ξ = 4p 2 sin 2 (θ 0 /2) fixed and large. When θ = θ 0 , u −s → ∞ in this limit, with t = −ξ remaining fixed and large. Similarly, when θ = π − θ 0 , t −s → ∞, with u = −ξ fixed. For the first three invariant amplitudes squared, we obtain in this limit Evidently, the amplitude squared for the process e + e − → γγ is suppressed in this limit, as compared to the other three. So this contribution can be neglected. 9 The kinematical factors B 1 and B 2 appearing in Eq. (2.24) are given by while B 3 is given by In the double scaling limit, we obtain at both θ 0 and π − θ 0 This factor is precisely equal to the cusp anomalous dimension in QED, e 2 Γ(ϕ)/4π 2 , via the relation ξ = 2m 2 (cosh ϕ − 1), which controls the vacuum expectation value of a Wilson loop with a cusp of angle ϕ [36,54]. Collecting all factors together, we arrive at a universal formula for the leading perturbative differential entanglement entropy, per unit time, per particle density, in the double scaling limit where the sums are over all charged particles in the initial state;ê i is +1 for positrons and −1 for electrons, and the cusp angles ϕ ij are defined via the relation ξ ij = 4p 2 i sin 2 (θ 0 /2) = 2m 2 (cosh ϕ ij − 1). Thus, the coefficient of the logarithmic singularity is given in terms of the cusp anomalous dimension of QED, averaged over all pairs of interacting particles via the above equation.

Next to leading order analysis for electron/electron scattering
In this section we extend the perturbative analysis to the next to leading order and obtain O(e 8 ) corrections to the Renyi and entanglement entropies. For simplicity, we focus on a Fock basis computation concerning the case of two-electron scattering. Thus we set the dressing function f µ equal to zero. As we have seen from the leading order analysis, the singular, logarithmically divergent part of the Renyi and entanglement entropies is independent of the dressing.

3.1
The large eigenvalue of ρ H and the Renyi entropies to order e 8 Two-photon outgoing states, |βγ 1 γ 2 , contribute to next to leading order and must be properly taken into account. For such a state to be in the soft part of the Hilbert space H S , the sum of the energies E γ1 + E γ2 must be below the reference energy scale E. Therefore, The outgoing state is given by where T βγ1γ2,α = βγ 1 γ 2 | iT |α is the amplitude for the two incoming electrons (described by the state |α ) to scatter and produce the electron/positron state |β , emitting at the same time two photons (|γ 1 γ 2 ). The set {|β } comprises of electron/positron states with zero number of photons. To order e 8 , we must consider states with two electrons and states with two electrons and an electron/positron pair. 10 The amplitudes T βα are of order e 2 for the two-electron case and e 4 for the case of three electrons and a positron. Likewise, the amplitudes T βγ,α , describing the emission of an additional photon, are of order e 3 and e 5 , respectively. Finally, T βγ1γ2,α is of order e 4 if |β is a two-electron state and of order e 6 if |β is a state of three electrons and a positron. States with more than two photons will not affect the entanglement and Renyi entropies to order e 8 . So we do not write them explicitly in the equation above. The reduced density matrix is obtained upon taking a partial trace over H S : ρ H = Tr H S |Ψ out Ψ| out . We define the state and write the reduced density matrix as where to order e 8 G is given by 10 As it turns out, the latter do not contribute to the Renyi and entanglement entropies to this order.
The first sum in the first line includes cross-terms between two-electron states and states with three electrons and a positron. Using the fact that at finite λ, T * αγ,α = T * αγ1γ2,α = 0 by energy conservation and β T βα T * βγ,α = 0 by incompatibility of the corresponding energy-momentum δ-functions, we see that G annihilates |Φ to order e 8 . So |Φ is an eigenstate of ρ H with eigenvalue Φ|Φ . Thus to this order there is a "large" eigenvalue given by Using unitarity Eq. (2.3), the eigenvalue can be written as follows States with three electrons and a positron do not contribute to ∆ to order e 8 , and so the sums over β in the expression above can be restricted to two-electron states. The rest of the non-zero eigenvalues of ρ H coincide with the non-zero eigenvalues of G. These split into two groups, the eigenvalues of order e 6 (which may receive corrections of order e 8 ) and the eigenvalues of order e 8 . The sum of all of these must be equal to ∆, as can be verified by computing TrG.
The traces Tr(ρ H ) m for integer m ≥ 2 are given by This equation can be verified by direct computation of the traces to order e 8 -see Appendix D. The corresponding Renyi entropies take the form Notice that to all orders, the traces satisfy the inequalities 0 < Tr(ρ H ) m ≤ 1. So the perturbative result Eq. (3.8) breaks down at large m ∼ 1/∆. To regulate the perturbative expansion of the entanglement entropy, we consider instead the generating functional [55] The entanglement entropy is recovered in the limit w → 1 The Taylor expansion of G(w, ρ H ) is convergent for |w| < 1 and gives In fact, setting w = e −κ , for κ ∼ ∆, suppresses the contributions of the traces for large integer m ∼ 1/∆ and regulates the perturbative expansion. Using Eq. (3.8), we obtain to order e 8 at fixed w The leading perturbative entanglement entropy computed in the previous section is recovered to be −∆ ln e 6 . Because of the logarithmic cutoff term, there are ambiguities in the analytic, cutoff independent part. The perturbative expansion does not commute with the w → 1 limit.

Entanglement entropy to order e 8
We proceed now to discuss next to leading order corrections to the entanglement entropy. In particular, we will study the structure of the non-analytic part (in the electron charge e) of the entanglement entropy, which depends on the non-vanishing, "small" eigenvalues of ρ H . As we remarked in the previous section, these eigenvalues coincide with the (non-zero) eigenvalues of G. They split into two groups, of order e 6 and e 8 , respectively, and their sum is given by i =Φ λ i = 1 − λ Φ = ∆. Let ∆ 6 and ∆ 8 be the order e 6 and the order e 8 parts of ∆, respectively.
In Appendix E, we describe how to obtain the shifts of the order e 6 eigenvalues of G using second order perturbation theory. Let us denote the sum of these shifts with∆ 8 . It is then clear that the sum of the order e 6 eigenvalues of G is given by ∆ 6 +∆ 8 , while the sum of all the order e 8 eigenvalues by ∆ − (∆ 6 +∆ 8 ) = ∆ 8 −∆ 8 . These sums determine the non-analytic part of the entanglement entropy Thus the non-analytic part is given by 3 ln e 6 (3.15)

IR divergences to order e 8
The IR divergent terms in the Renyi and entanglement entropies can be obtained by examining the singular part of ∆ (which gives the shift of the "large" eigenvalue from unity). When the energy scale E is sufficiently small, we can apply soft theorems for real photon emission to obtain (in the continuum limit) Now this expression must be computed to order e 8 . So the amplitude squared in the first term should be calculated to order e 6 , and, therefore one-loop diagrams contribute. As a result, logarithmic IR divergences due to virtual photons appear. Soft theorems for virtual photons give [10] where scattering amplitude squared (T βα T * βα ) Λ does not include contributions from virtual photons with energy below the cutoff scale Λ ∼ E. The factor ln(Λ/E) is negligible when the reference scales E and Λ are of the same order. Substituting into the expression for ∆, yields The first term is in accordance with an exponentiation pattern. Notice however the appearance of the last term which originates from two-photon final states in which one photon is soft and the other is hard. The logarithmic divergence is now proportional to the amplitude squared for real (hard) photon emission, and a different kinematical factor. So the exponentiation pattern of the IR logarithmic divergences is an intricate one.

To all orders entanglement in electron/electron scattering
In this section we attempt to generalize some of the results of the previous sections to all orders in perturbation theory so as to gain insight into the exponentiation pattern of the IR logarithmic divergences, and to link our perturbative results with the results of [11,12,13,14] (see also [15,16]) on the decoherence of the density matrix ρ H . As in Section 3, we focus on a Fock basis computation associated with twoelectron scattering for simplicity.
We denote by |γ = |γ 1 γ 2 , . . . , γ n a generic multiphoton box state. The number of photons n is nonzero and finite, but arbitrary. It will be convenient to order the photons in ascending order, according to their energy: E 1 ≤ E 2 ≤ · · · ≤ E n . If the total energy of the photons is less that the reference energy scale E, If no such m exists, then the state |γ is in H H . Photons with energy scaling with the IR cutoff λ will appear in the soft part of the Hilbert space in the continuum, λ → 0 limit. 11 Adopting these notations, the outgoing state can be written as at any order in perturbation theory. The set {|β } comprises of states with 2 + l electrons, l positrons states and zero number of photons, with l = 0, 1, . . . . The reduced density matrix is given by ρ H = Tr H S |Ψ out Ψ| out . As before, we define in terms of which the reduced density matrix is given by The expression for G is valid to all orders in perturbation theory.
Explicitly λ Φ is given by The sum of the rest of the eigenvalues of ρ H is given by We remark that the amplitudes T βα , T βγ,α appearing in the above expressions for ∆ and λ Φ are hard, in that the final state does not include soft photons with energies scaling with λ in the continuum limit. At any finite order in perturbation theory, these amplitudes are plagued by IR logarithmic divergences, due to virtual soft photons running in the loops [10], when the momenta of the two-electron states α and β differ. Since the corresponding transition processes do not involve the emission of soft radiation, there are no similar divergent contributions from real soft photons to λ Φ (and ∆) to cancel the virtual IR divergences.
However, the term T αα + T * αα is free of any IR divergences, order by order in perturbation theory, since this is related to the total inclusive cross-section in the state α by unitarity, Eq. (2.3): The corresponding transition processes include the emission of real soft photons. As is well known, such inclusive transition rates are IR finite, with the IR divergences due to virtual soft photons cancelling against the ones due to real photon emission [10]. We emphasise that T αα is a box transition amplitude, scaling inversely proportional with the volume of the box in the λ → 0, continuum limit: where T is the time-scale of the experiment, iM αα is the invariant amplitude at infinite volume 12 and E 12 is the total energy in the center of mass frame of the two electrons of the state |α . In particular, where v ij is the relative velocity of the two electrons in the state α and Σ α is the total, inclusive crosssection in this state. As a result, at any finite order in perturbation theory IR logarithmic divergences appear in the expressions for the Renyi and the entanglement entropies, as we have seen explicitly at leading and next-to-leading order in Section 2 and Section 3.
To all orders, however, the virtual IR logarithmic divergences exponentiate [10] leading to the vanishing of the hard amplitudes T βα , T βγ,α when the momenta of β and α differ. The vanishing of these hard amplitudes can be also understood as a consequence of symmetry [1,2]. Taking into account the scaling of the box amplitude T αα in the large volume limit, we conclude that to all orders, λ Φ becomes equal to a diagonal element of ρ H free of any IR divergences in λ. Indeed, to all orders, the reduced density matrix ρ H assumes a diagonal form in terms of the momentum indices exhibiting decoherence [11,12,13,14]. The diagonal elements (in momentum) are given in terms of inclusive Bloch-Nordsieck rates associated with box states, and are free of IR divergences. Assuming that Σ α is finite, we see that to all orders, the reduced density matrix is dominated by a large eigenvalue, λ Φ , in the continuum, λ → 0 limit. The rest of the eigenvalues must be small, scaling inversely proportional with powers of the volume of the box. Despite the fact that the number of the small eigenvalues also grows with the volume, the traces Tr(ρ H ) m for integer m ≥ 2 are given by and the corresponding Renyi entropies take the form The Renyi entropies per particle flux, per unit time remain finite in the limit. They are proportional to the total cross-section Σ α , free of IR divergences. Let us discuss the dependence of Σ α on the reference scale E, taking into account the scaling of the box amplitudes with the volume and the vanishing of the purely hard amplitudes to all orders in the continuum limit. Eq. (4.8) gives The total cross-section Σ α is proportional to a sum of inclusive Bloch-Nordsieck rates, of the form, for the initial two electrons α to scatter and produce the hard state βγ, emitting at the same time any number of soft photons with total energy less than E. Notice that in the continuum limit, the sums over β and γ give rise to multiple integrals over the momenta of the hard final particles with appropriate volume factors (as well as discrete sums over the polarization indices) -see Section 2.3. Each such contribution scales proportionally with T /V in the continuum limit, irrespectively of the number of particles in the final state. The inclusive rates Eq. (4.18) have been computed to all orders in [10], taking into account the cancellation of the IR divergences due to virtual soft photons and real soft photon emission where b βα 1 − π 2 B βα 2 /3 [10]. Therefore, we obtain The ratio E/Λ remains fixed in the λ → 0 limit. To all orders, the small parameter controlling the expansion of the Renyi entropies in the continuum limit is set by the inverse of the product of the size of the box L and the center of mass energy of the initial state. The entanglement entropy is not analytic in the volume of the box. (More precisely it is not analytic in the small dimensionless parameter discussed above). The non-analytic part is governed by the small eigenvalues, via the expressions −λ i ln λ i , for i = Φ. Ignoring the spin polarization structure, we may approximate the small eigenvalues with the diagonal elements of ρ H , which are given in terms of inclusive rates where |γ is a hard photon state. This diagonal element scales with T 2 /V N f in the large volume limit, where N f is the number of particles in the final state, taking into account the energy-momentum conserving δ-functions. It gives a contribution −D βγ,βγ ln(D βγ,βγ ) to the entanglement entropy. So terms logarithmic divergences in the size of the box (or the IR cutoff λ), of the form ln V N f /T 2 , in the expression for the entanglement entropy per particle flux per unit time persist to all orders in the continuum, large volume limit.
One way to recover the entanglement entropy is to take the w → 1 limit of the generating functional G(w, ρ H ), Eq. (3.10), whose Taylor expansion is given in terms of the traces Tr(ρ H ) m (and hence the exponentials of the Renyi entropies) -see Eq. (3.12). How is it then that the all orders entanglement entropy (per particle flux, per unit time) retains non-analytic, logarithmic behavior with respect to the size of the box or the IR cutoff? Notice that for any given large m, the scaling of Tr(ρ H ) m with the volume, Eq. (4.14) (consequently the scaling of the corresponding Renyi entropy, Eq. (4.15)), is valid for sufficiently large volume. At large but fixed volume, however, this scaling breaks down for large enough integers m > 1/∆, since the traces satisfy the strict inequalities 0 < Tr(ρ H ) m ≤ 1. A similar behavior was observed in the case of the perturbative results of the previous sections (see the discussion around Eq. (3.10)), leading to the non-analyticity of the entanglement entropy in the electron charge e. In the case at hand, the expansion parameter scales inversely proportional with the size of the box L. As a result, the w → 1 limit does not commute with the large L expansion of the traces. Indeed when |w| < 1, the contribution of the large m traces in Eq. (3.12) is suppressed. Using Eq. (4.14) and the Taylor expansion of the generating functional, we obtain at large but fixed volume However, this yields a diverging result in the limit w → 1, indicating that we cannot reverse the order of limits. Similarly, another way to recover the entanglement entropy is to analytically continue the expression for the Renyi entropies for integer m ≥ 2 to general real values, and take the limit S ent = lim m→1 S m . However, this limit does not commute with the large volume expansion of the Renyi entropies, as can be seen from Eq. (4.15). Essentially, the entanglement entropy is non-analytic in the eigenvalues of the reduced density matrix. Since to all orders the small eigenvalues scale inversely proportional with the volume of the box, the entanglement entropy retains non-analytic behavior in the inverse of the volume or the IR cutoff.

Conclusions
In this work we study generic scattering processes in QED in order to establish measures of the entanglement between the soft and hard particles in the final state. The initial state consists of an arbitrary number of Faddeev-Kulish electrons and positrons. To regulate the computation of the entanglement and the Renyi entropies, we place the system in a large box of size L, which provides an infrared momentum cutoff λ ∼ 1/L. The perturbative computations are carried out at fixed λ, for sufficiently small coupling e. The density matrix for the hard particles is constructed via tracing over the entire spectrum of soft photons, including those in the clouds dressing the asymptotic charged particles. The leading perturbative behavior of the Renyi and entanglement entropies is governed by a large eigenvalue, which we compute at leading order, and next to leading order for two-electron scattering, in perturbation theory. The leading entanglement entropy and Renyi entropies for integer m ≥ 2 are logarithmically divergent with respect to the IR cutoff λ for all scattering cases we consider. We find that the coefficient of the logarithmic divergence is independent of the particle dressing and exhibits certain universality properties, irrespectively of the detailed properties of the initial state. The dominant contributions arise from two-particle interactions at forward and backward scattering. In the relativistic limit, these dominant contributions are proportional to the cusp anomalous dimension in QED. We see that the entanglement entropy in momentum space can be IR divergent, with the coefficient of the divergence containing physical information. We also study next to leading order corrections to the non-analytic part (in the electron charge e) of the entanglement entropy, and comment on the pattern of the exponentiation of the IR logarithmic divergences.
Since the perturbative expansion does not commute with the continuum, λ → 0 limit due to the appearance of IR divergences, we proceed to study the behavior of the Renyi and entanglement entropies to all orders in perturbation theory. We focus on Fock basis computation associated with two-electron scattering for simplicity. We derive an exact expression for the large eigenvalue of the density matrix in terms of hard scattering amplitudes, which is valid at any finite order in perturbation theory. At any finite order in the electron charge e, this eigenvalue is plagued with IR logarithmic divergences due to soft virtual photons. To all orders in perturbation theory, the IR divergences exponentiate, leading to a finite result. In the large volume limit, the shift of this eigenvalue from unity governs the behavior of the Renyi entropies for integer m ≥ 2. The Renyi entropies, per unit time, per particle flux, turn out to be proportional to the total, inclusive cross-section in the initial state, which is free from any IR divergences. The contributions of the small eigenvalues remain subleading, despite the fact that their number grows with the volume in the continuum limit. The entanglement entropy though retains non-analytic, logarithmic behavior with respect to the size of the box, even to all orders in perturbation theory, in the continuum limit. This reveals strong entanglement between the soft and hard particles produced in the scattering process.
It would be interesting to extend our analysis to gravitational and black hole scattering processes, and study the correlation between the soft photons and gravitons with the hard particles produced in the process. Interesting correlations could be uncovered if we examine four-dimensional celestial scattering amplitudes involving boost eigenstates, which are sensitive to both UV and IR physics. The factorization structure of the S-matrix found in [18] could be important in order to quantify the entanglement between the soft and the hard degrees of freedom, and the information carried the unobserved soft particles.

A Notations and conventions
We work in the Lorenz gauge ∂ µ A µ = 0, in which the electromagnetic gauge field satisfies the wave equation ∂ ν ∂ ν A µ = 0. We adopt a mostly plus signature metric. In terms of creation and annihilation operators, the gauge field is given by with the Greek indices taking values from 0 to 3. We adopt the conventions of [56]. The polarisation vector µ r ( k) satisfies the following orthonormality relations where ζ 0 = −1, ζ 1 = ζ 2 = ζ 3 = 1. The commutators of the photon creation and annihilation operators are given by Notice also the Gupta-Bleuler condition on physical states Finally the non-trivial anticommutators associated with the electron/positron creation and annihilation operators are, respectively The indices r, s, . . . of the spinor creation and annihilation operators take the values 1, 2. We hope they will not be confused with the same symbols used to label the polarisations of the electromagnetic field. In the continuum, single particle states are normalized so that they satisfy the Lorentz invariant norm On the other hand, box single particle states have discrete momenta and they are unit normalized. We also list some partial traces over the soft part of the Hilbert space, H S , associated with dressed ket-bra operators, which are useful in obtaining the reduced density matrix ρ H -see [17] for explicit derivations. We first have where f β |f β is given in Eq. (1.6), in the large volume, continuum limit. Taking | q γ | < E d gives Finally, when both | q γ |, | q γ | < E d , we obtain Partial traces for the cases in which two or more soft photons are present in the initially undressed states are higher order in the dressing function f µ β ( q γ ).

B Leading order analysis
The terms in the expansion of the reduced density matrix Eq. (2.7) relevant for the leading order computation of the Renyi and the entanglement entropies are where Notice that the leading contributions in C β are of order e N , in C βγ of order e N +1 , in D β,β of order e 2N , in D βγ,β of order e 2N +1 and in D βγ,β γ of order e 2N +2 . Also in the expression for C β , the second term and terms in the ellipses vanish when β = α [17]. So C α = T αα to all orders.

B.1 Traces to leading order
To facilitate the perturbative analysis, we set Since ε is of order e 2 , we need to expand (ρ H ) m to cubic order in order to obtain the leading perturbative contributions to the Renyi entropies. The fact that (ρ 0 ) 2 = ρ 0 and the cyclic property of the trace limit the number of structures we need to consider. As in the case of two-electron scattering, only linear and quadratic structures contribute to leading order (e 6 ). The cubic structures in ε turn out to vanish (to order e 6 ).
To linear order in ε, it suffices to consider the term ερ 0 whose trace is The ellipses include terms of higher order than e 6 and do not contribute to the entanglement entropy at leading order. The quadratic terms include ε 2 , ε 2 ρ 0 , ερ 0 ερ 0 . The trace of the first term is The third line is absent for the Fock basis computation. For the second quadratic trace we obtain Notice the appearance of off diagonal elements in the perturbative expansions for both Trε 2 and Trε 2 ρ 0 . At any finite order in perturbation theory, the off diagonal elements are non-zero and contain IR logarithmic divergences in λ. Both of these traces contribute to the Renyi and the entanglement entropies to leading order (e 6 ).
Finally we have which does not contribute to order e 6 by unitarity. The cubic terms are ε 3 , ε 3 ρ 0 , ε 2 ρ 0 ερ 0 , ερ 0 ερ 0 ερ 0 . However these terms do not contribute to leading order. Firstly which does not contribute to order e 6 by unitarity.
For the next terms we have which do not contribute to the leading entanglement entropy.
To order e 6 , only Tr(ερ 0 ), Tr(ε 2 ) and Tr(ε 2 ρ 0 ) contribute. To proceed, we use the perturbative expansion for the overlap of the coherent states Eq. (1.6). Recall that we first expand to a given order in perturbation theory, keeping the volume of the box and the cutoff λ finite. The continuum limit (where the volume of the box is taken to infinity and the cutoff λ to zero) is taken in the end. If the logarithmic IR divergences do not cancel, the perturbative expansion breaks down in the continuum limit. We will also need the relations between dressed and undressed amplitudes Eq. (2.5) and Eq. (2.6). We will use these relations to express the traces in terms of undressed amplitudes, which we can calculate more easily via the Feynman diagrams.
For the first trace we have The second term in parentheses is zero due to energy conservation. In the continuum limit, the first term in parentheses yields up to λ independent multiplicative factors. Since F αα ( q, r ( q)) is finite and nonsingular in the limits | q|, λ → 0, the contribution of this term is at most of order E 2 d , and can be neglected relative to other contributions.
For the quadratic traces we obtain Since F βα is of order e 3 and f α is of order e, T βα must be calculated at tree level. So it does not lead to any IR divergences in the limit λ → 0. In the continuum limit, the last lines of expressions Eq. (B.19) and Eq. (B.20) are given by the following integral (up to smooth factors in the limit λ → 0 and volume factors) The integral is not divergent in the limit | q| → 0. In fact it is of order E d making a negligible contribution to the Renyi and the entanglement entropies. Let us calculate the quantity Tr(ρ H ) 2 to order e 6 .

B.2 The large eigenvalue to leading order
The large eigenvalue of ρ H is given by As shown in [17] (and also discussions in the previous subsection), the terms in the last line give suppressed contributions, of order the characteristic energy of the photons in the clouds E d , and can be dropped. Using unitarity and the relations between dressed and undressed amplitudes (Eq. (2.5) and Eq. (2.6)), λ Φ can be written as The other nonvanishing eigenvalues coincide with the nonvanishing eigenvalues of the matrix G, and each is of order e 6 , at least. Their sum is equal to ∆. Let us verify this by computing explicitly the trace of G to order e 6 : Using unitarity, this can be written as C The leading entanglement entropy in the continuum limit In the generic case, non-trivial disconnected diagrams in which only two particles interact with each other contribute to the singular part of the leading entanglement entropy, Eq. (2.23), in the continuum, large volume limit. Each disconnected line associated with a noninteracting incoming (outgoing) particle i (j), with momentum p i,in ( p j,out ), contributes to the Feynman diagram with a factor (2π) 3 2E i,in δ 3 ( p i,in − p j,out ). We consider the case where the momenta of the initial particles differ from each other. The amplitude iM(α → me − + ne + ) can be written as a sum of the following three terms. The first term arises from Feynman diagrams in which only two electrons interact, while the rest of the initial particles do not: We have denoted the momenta of the electrons by p l and of the positrons by k l . For brevity, we have omitted the inclusion of polarizations indices, since we average over these in the end. In figure 2 we show a disconnected diagram contributing to I 1 to leading order, for an initial state of 3 electrons and 2 positrons. Also, iM 1 is the scattering amplitude for the process e − e − → e − e − , calculated at tree level, and iM 1t is the t-channel amplitude. The permutation of m (n) objects is denoted by σ (σ ). The second term arises from Feynman diagrams in which only two positrons interact: where iM 2 is the ampltude for the process e + e + → e + e + and iM 2t is the t-channel amplitude. See figure  3. Finally, when one electron and one positron interact, we get the following contribution where iM 3 is the total amplitude for the process e − e + → e − e + , calculated at tree level. See figure 4 for the t-channel and s-channel diagrams. To leading order, the amplitude iM(α → (m − 1)e − + (n − 1)e + + 2γ) takes the form where iM 4 is the tree-level amplitude for the process e − e + → 2γ. The momenta of the photons are denoted by q 1 and q 2 . Relevant diagrams are shown in figure 5.
To calculate the entanglement entropy, we need to determine the square of the absolute value of the amplitudes iM(α → me − +ne + ) and iM(α → (m−1)e − +(n−1)e + +2γ), and integrate over the momenta of the final particles (Eq. (2.23)). Each amplitude is obtained as a sum of Feynman diagrams in which only a pair of particles interact, as we explained before. Taking the absolute value squared, a number of terms appear, each containing a product of 2(m + n − 2) δ-functions (in three-dimensional momentum space). There is also an overall 4-momentum δ-function squared, imposing energy-momentum conservation p 1,out p 2,out q 1,out q 2,out k 2,out p 1,out p 2,out p 3,in k 1,in k 2,in 1 Figure 5: A disconnected Feynman diagram contributing to iM(3e − + 2e + → 2e − + e + + 2γ). The u-channel diagram is obtained by exchanging the outgoing photon momenta between the outgoing photon lines.
-see Eq. (2.23). In most cross-terms that are not perfect squares, the δ-functions are incompatible with each other. So, upon integrating over the momenta of the noninteracting particles, such cross terms give vanishing contributions.
The dominant contributions come from perfect squares (of the form αα * ). For such terms, the sign of the diagram does not contribute. For such terms, the products of δ-functions from each factor are identical. If we integrate over the momenta of the n + m − 2 = N − 2 noninterracting particles, we get a factor [(2π) 3 δ 3 (0)] N −2 = V N −2 . There are also cross terms with nonvanishing contributions, but these are suppressed in the continuum, large volume limit, since they are proportional to a smaller power of the volume V . We will drop such subleading terms. Notice however that the precise volume dependence is non-trivial. A detailed example with three charged particles in the initial state is presented in the following subsection.
For a given pair of interacting incoming electrons, i and j, there are m!n!/2 dominant terms, each proportional to |iM 1 ( p i , p j ; p i , p j )| 2 , where we denote the corresponding, outgoing interacting particles by i , j . Likewise, we get m!n!/2 dominant terms, proportional to |iM 2 ( p i , p j ; p i , p j )| 2 , for a pair of interacting positrons. When the i-th electron scatters against the j-th positron, we get m!n! dominant terms, each proportional to |iM 3 ( p i , p j ; p i , p j )| 2 . Finally, there are (m − 1)!(n − 1)! dominant terms produced, proportional to |iM 4 ( p i , p j ; p i , p j )| 2 , when the i-th electron and the j-th positron annihilate to produce two outgoing photons. In this latter case i , j denote the two outgoing photons.
There are two more integrations involved in Eq. (2.23), over the momenta of the interacting outgoing particles. Integrating over p j , imposes momentum conservation and gives rise to an additional volume factor. We can always choose to work in the center of mass frame for the two interacting incoming particles i and j ( p i = − p j ). The center of mass energy of these particles is denoted by E ij = 2E i and the relative velocity by v ij = 4| p i |/E ij . We also denote the scattering angle by θ. Integrating over the magnitude of p i results in a factor (2π)δ(2E i − 2E i ) = T , equal to the time scale of the experiment, and a factor of the relative velocity v ij .

C.1 Three-particle scattering
In this appendix we illustrate the continuum, large volume analysis of Section 2.3 for a particular case concerning a three-electron initial sate. We take the momenta of the incoming particles to be different. The leading perturbative entanglement entropy takes the form The leading contributions to the amplitude iM(α → 3e − ) (arising from disconnected tree-level Feynman diagrams in which only a pair of two electrons are interacting) are given by Upon taking the absolute value squared of this amplitude, we get terms that are perfect squares, like for instance, and cross-terms. Let us calculate the contribution from the perfect square above. We work in the center of mass frame of the interacting particles 1, 2. The center of mass energy is denoted by E 12 = 2E 1 , the relative velocity of the incoming particles by v 12 = 4| p 1 |/E 12 and the scattering angle by θ. Integrating over the final momenta p 3 , p 2 , and | p 1 |, produces a factor V 2 T v 12 /4, giving a net contribution There are cross-terms with incompatible δ-functions, for example, Upon integration over p 3 , we get a factor (2π) 3 δ 3 ( p 3 − p 2 ), which vanishes for generic states in which the incoming three particles have different momenta. Another example that gives a vanishing contribution is the following We work in the center of mass frame of the two interacting particles (1 and 2). Integrating over p 3 , p 2 imposes p 3 = p 2 = p 3 . This is incompatible with momentum conservation unless p 3 = 0. But then generically the energy will not be conserved. A cross-term with a non-vanishing but suppressed contribution in the large volume, continuum limit is Integration over p 3 , p 2 imposes p 3 = p 3 and p 2 = p 2 . Integration over p 1 , using one of the overall momentum δ-functions imposes p 1 = p 1 and yields a factor of V T 2 . The contribution of this term is of order T 2 /(V 2 E 2 1 E 2 E 3 ) times dimensionless amplitudes and kinematical factors. Letting T L/c, where L is the size of the box and c the speed of light, the ratio of this contribution to Eq. (C.8) is of order 1/(L 2 E 2 E 3 ), which becomes vanishingly small when the IR cutoff 1/L is taken to zero.
There are also linear combinations of |β states that are eigenstates of G (0) with eigenvalues of order e 8 (or higher). For example, up to order e 6 , momentum conservation requires that two-electron |β states with vanishing net momentum are annihilated by G (0) . Presumably, certain linear combinations of such zero-momentum two-electron states can be eigenstates of G (0) with eigenvalues of order e 8 . Certain linear combinations of |β states involving three electrons and a positron can also be eigenstates of G (0) with eigenvalues of order e 8 . Notice however that the perturbation G (1) does not mix these eigenstates with order e 6 eigenstates, or the order e 6 eigenstates among themselves. To find the e 8 corrections to the leading eigenstates of G (0) , we must apply second order perturbation theory.
The other eigenstates of G (0) , with eigenvalues of order e 8 , comprise linear combinations of single photon states (and energies E α > E βγ > E α − E) only. We shall denote these eigenstates by βγ and the corresponding eigenvalues by λ (0) βγ . The perturbation G (1) mixes these eigenstates with the order e 6 eigenstates of G (0) .
Assuming that the order e 6 eigenstates of G (0) have been found to be non-degenerate, we can apply second order perturbation theory to determine the order e 8 corrections to the eigenvalues