Electron-Positron Annihilation Freeze-Out in the Early Universe

Electron-positron annihilation largely occurs in local thermal and chemical equilibrium after the neutrinos fall out of thermal equilibrium and during the Big Bang Nucleosynthesis (BBN) epoch. The effects of this process are evident in BBN yields as well as the relativistic degrees of freedom. We self-consistently calculate the collision integral for electron-positron creation and annihilation using the Klein-Nishina amplitude and appropriate statistical factors for Fermi-blocking and Boseenhancement. Our calculations suggest that this annihilation freezes out when the photon-electronpositron-baryon plasma temperature is approximately 16 keV, after which its rate drops below the Hubble rate. In the temperature regime near 16 keV, we break the assumption of chemical equilibrium between the electrons, positrons, and photons to independently calculate the evolution of the chemical potentials of the electrons and positrons while computing the associated collision integrals at every time step. We find that the electron and positron chemical potentials deviate from the case with chemical equilibrium. While our results do not affect the interpretation of precision cosmological measurements in elucidating the standard cosmological model, these out of equilibrium effects may be important for testing physics beyond the standard model.


I. INTRODUCTION
The chronology of the early universe is the story of a hot, dense, and entropy-dominated plasma that may have, at one time, included thermal populations of all particles in the Standard Model (the so-called "primordial soup"). While particle populations are well described by thermal and chemical equilibria, the standard cosmological model incorporates phase transitions that change the nature of the fundamental forces in the Standard Model as well as epochs where various processes freezeout of equilibrium. Freeze-out from equilibrium is not an instantaneous process, resulting in deviations of particle phase-space distributions as the approximations of local thermal and chemical equilibrium become invalid [1][2][3][4].
Upcoming high-precision cosmological observations [5,6] introduce the promise of the early universe as a testing ground for the Standard Model and the standard cosmological model. Out-of-equilibrium effects that occur when freeze-out is not instantaneous are an especially promising regime for revealing Beyond Standard Model (BSM) physics [7][8][9][10][11][12].
The well-studied epochs of weak decoupling and weak freeze-out occur when the temperature of the photonelectron-positron-baryon plasma, T , is approximately an MeV [13]. Weak decoupling (T ∼ 3 MeV) occurs when neutrinos no longer efficiently exchange energy with the plasma so they are no longer in thermal equilibrium. Weak freeze-out (T ∼ 0.7 MeV) occurs when the charged current lepton-nucleon processes inefficiently inter-convert protons to neutrons, freezing-out the neutron-to-proton ratio as they fall out of chemical equilibrium with the leptons. These processes occur over many Hubble times (expansion timescales), requiring a self-consistent treatment of the Boltzmann evolution of neutrino distribution functions and a nuclear network to describe with the consequent Big Bang Nucleosynthesis (BBN) [1]. The results of these computationallyintensive calculations may be observable in the energy density of neutrinos (as parameterized by N eff ) and in primordial abundances as the next generation of cosmic microwave background (CMB) observatories and thirtymeter telescopes come on-line [5,6].
Careful analysis of this epoch reveals a variety of Standard Model processes that affect the value of N eff and BBN yields because weak decoupling and weak freezeout are not instantaneous. The better known of these effects is the increase of N eff from its standard, instantaneous decoupling value of 3 to the accepted value of 3.046 due to a variety of effects including electron-positron annihilation, out-of-equilibrium neutrino up-scattering off the photon-electron-positron-baryon plasma, and finite temperature QED radiative corrections to the in-medium electron and photon masses [1,[14][15][16]. Other work has looked into the effects of neutrino quantum kinetic behavior [2,3] and non-standard cosmological models including the introduction of non-zero lepton numbers [17]. BBN yields are also affected through their influence on neutrino distribution functions and the time-temperature relationship throughout the BBN epoch.
These in-depth studies are necessary to disentangle out-of-equilibrium Standard Model effects from BSM physics, both of which may have an observable signal in N eff and the BBN yields. During the BBN epoch, electron-positron creation becomes kinematically suppressed compared to annihilation, driving the equilibrium process -e − e + γγ -to the right, resulting in the number densities of electrons and positrons plummeting. In the standard cosmological model, BBN yields and N eff are calculated as if this process always occurs in local thermal and chemical equilibrium [18].
The goal of this work is to closely examine the pro-cess of electron-positron annihilation as it freezes out from its local chemical equilibrium. We will determine if and when the electron and positron number densities freeze-out, investigate the consequent effects on cosmological observables, and speculate on whether the effects of freeze-out should affect the standard BBN calculations. To do so, we calculate the collision integral for electron-positron annihilation using the Klein-Nishina amplitude with the appropriate Fermi-blocking/Boseenhancing statistical factors. Appendix B details our methodology to self-consistently and efficiently calculate these integrals using methods similar to those in Ref. [1]. Throughout this work we use the convention that h = c = k B = 1. In Section II, we use physical principles, including chemical equilibrium between the electrons and positrons, to introduce the standard cosmological computation to evolve the electron and positron distributions. Next, in Section III, we calculate the collision integral associated with electron-positron annihilation to infer when freeze-out occurs. We subsequently use the collision integral to separately evolve the electron and positron distributions without the assumption that the charged leptons maintain chemical equilibrium in Section IV. In Section V we will discuss the results and present conclusions. Appendix A details the equations of motion derived and solved in this work.

II. CHEMICAL EQUILIBRIUM
As evidenced by the cosmic microwave background, the primordial plasma is well described by local thermodynamic equilibrium. Scattering between photons, electrons, positrons, and baryons remains fast compared to the Hubble rate, maintaining thermal distributions well into the BBN epoch. This means the photons, electrons, positrons, and baryons share a single plasma temperature T , and the electron/positron distributions are Fermi-Dirac. When further assuming the process e − e + γγ to be in local chemical equilibrium, the chemical potentials of the involved particles maintain µ e − + µ e − − 2µ γ = 0, stipulating µ e + = −µ e − . For chemical potentials µ e ± , we define the electron and positron degeneracy parameters η − ≡ µ e − /T and η + ≡ µ e + /T . Chemical equilibrium is implemented by setting η − = η and η + = −η.
Given local thermodynamic and chemical equilibrium, the electron and positron momentum/energy distributions are where E p = p 2 + m 2 e . The evolution of the electron and positron distributions is encapsulated in the evolution of T and η. To evolve these quantities, we introduce two physical principles -conservation of comoving entropy, and the effects of weak interactions on the relative number of electrons and positrons.
The early universe is homogeneous and isotropic, precluding spatial heat flows. In the absence of out-ofequilibrium decays, the total entropy in a comoving volume is conserved, Here a is the scale factor in the Friedmann-Lemaître-Robertson-Walker metric and entropy density is The quantities ρ and P are total energy density and pressure of all particles in the plasma, while electron and positron number densities are denoted n e ± . The scale factor evolves according to the Friedmann equation, where m pl is the Planck mass. The number densities of electrons and positrons are affected by the expansion of the universe, creation and annihilation, and charged-current weak interactions. The product na 3 is proportional to the number of particles in a comoving volume, which remains unchanged by the expansion of the universe. Creation and annihilation rates are fast compared to the weak rates, but are computationally expensive to calculate. In order to avoid the need to self-consistently calculate these rates with the evolution, we introduce the quantity (n e − − n e + )a 3 , which is insensitive to both the expansion of the universe and electromagnetic processes. Electromagnetic processes preserve the relative number of electrons and positrons, so the evolution of this quantity does not depend on the collision integral, even when the assumption of chemical equilibrium is broken.
Weak interactions that inter-convert neutrons and protons are the only processes that cause (n e − − n e + )a 3 to change. During the BBN epoch, these processes are e + + n p +ν e n p +ν e + e − .
When the above reactions proceed to the right, (n e − − n e + )a 3 increases, while the reverse reactions cause it to decrease. The rate of change of this quantity is d dt (n e − − n e + )a 3 = n n a 3 (λ νen + λ e + n + λ n ) − n p a 3 (λ e − p + λν ep + λν e e − p ).
The rates λ i , with i indicating the reactants, are computed in Ref. [19]. Rather than directly computing the number densities of neutrons (n n ) and protons (n p ), we can recast Eq. (9) in terms of the electron fraction resulting in where n b is the number density of baryons, and the quantity n b a 3 is constant. The right-hand-side comes from the rate of change of Y e [19], which accounts for the production and destruction of electrons and positrons by the weak interactions.
Together, the conservation of comoving entropy (Eq. 3), the evolution of the difference in comoving numbers of electrons and positrons due to weak interactions (Eq. 11), the evolution of Y e (Eq. 12), and the Friedmann equation (Eq. 5) form a system of coupled differential equations for dependent variables T , η, Y e , and a (see Section 1 of Appendix A). It should be noted that the definition of Y e , Eq. (10), is an integral relationship between all four dependent variables, so while the four aforementioned differential equations are not independent, we find it more efficient to simultaneously use them to evolve the four variables.
We elect an initial plasma temperature T i = 10 MeV, sufficiently high that we are assured that the electronpositron annihilation rate is much greater than the Hubble rate. At this temperature, the weak interactions (Eqs. [6][7][8] are in chemical equilibrium so that the initial electron fraction is [13] The initial degeneracy parameter, η i , is chosen to be consistent with the observed baryon-to-photon ratio from CMB observations [20], where ζ(3) ≈ 1.20206 and the subscript f indicates the values at the end of the calculation. (Here, we assume there is no BSM physics that changes the baryon-tophoton rates between the keV scale and recombination [7].) Since n b a 3 is conserved, the initial baryon number density is In the standard cosmological model [13,21], which accounts for the heating of the plasma when the rest mass energy in electrons and positrons is liberated by annihilation. The initial baryon number can be related to the initial degeneracy parameter through integration of the Fermi-Dirac distribution, assuming that T m e and η i 1. Eqs. (13)(14)(15)(16)(17) produce the initial conditions to evolve T , η, Y e , and a through the epoch where annihilation is important. Figure 1 presents the results of this chemical equilibrium solution. The number of particles in a comoving volume for each species in the photon-electron-positronbaryon plasma is plotted as a function of decreasing plasma temperature. When the plasma temperature drops below about an MeV, pair production becomes energetically suppressed, causing the numbers of electrons and positrons to fall precipitously, with positrons tracking electrons until electron numbers abruptly level off and positron numbers continue to plummet. As seen in Fig.  1, this difference is driven by rapid growth in η.
Local chemical equilibrium should hold so long as the creation/annihilation rates for e − e + γγ remain fast relative to the Hubble rate. The rapid decrease in electron and positron numbers causes the equilibrium creation/annihilation rate to similarly fall off. The following section utilizes the equilibrium solution to calculate the annihilation rate to assess the possibility of annihilation freeze-out.
One minor issue with this calculation is the lack of a BBN nuclear network calculation. BBN would "protect" neutrons from undergoing free neutron decay by locking up the free neutrons primarily in 4 He nuclei. In the standard cosmological model, this occurs at T ∼ 100 keV and results in Y e ≈ 0.88 at late times (while in our calculation Y e → 1 as free neutron decay eventually converts every neutron into a proton). To ensure that this approximation has negligible effects on our results, we performed the same calculation with a constant Y e equal to its initial value throughout the evolution, and found no significant changes in the conclusions. Having bracketed the solution that includes BBN, we expect that self-consistently following the nuclei through BBN will not affect the conclusions of this work.

III. ANNIHILATION RATES
To determine when the annihilation process e − e + → γγ freezes out, we compute the fractional rate of change of the number density of electrons and positrons, due solely to annihilation, The time derivative of the number densities can be calculated directly from the collision integral, Here, C a is the collision integral for the annihilation pathway, which is one of many terms in the Boltzmann equation for the rate of change of f e ± . Appendix B outlines the computation of the collision integral. Figure 2 presents R e ± calculated using T , η, Y e and a from the chemical equilibrium solution. Due to the paucity of positrons, R e − drops below the Hubble rate at 16 keV, so at this time the universe expands much faster than the rate at which electrons annihilate with positrons. However, R e + remains many orders of magnitude faster than the Hubble rate as it is much more likely for a positron to annihilate on the significantly more abundant electrons, so essentially every positron will annihilate with an electron.

IV. OUT OF EQUILIBRIUM
With annihilation freeze-out estimated to occur at 16 keV, we lift the assumption of local chemical equilibrium and reexamine this period more closely. Electrons and positrons continue to efficiently scatter with the plasma, so their distributions are still Fermi-Dirac, but allowing for deviations from chemical equilibrium implies that the electron and positron degeneracy parameters evolve separately. With independent electron and positron degeneracy parameters, we need to re-evaluate the four coupled equations of motion used in the chemical equilibrium solution, Eqs. (3,5,11,12), to account for the separate degeneracy parameters (see Section 2 in Appendix A). In addition, we introduce the time evolution of the sum of electron and positron numbers in a comoving volume, (n e − + n e + )a 3 , which is modified by weak interactions as well as creation and annihilation processes: The first two terms on the right-hand-side account for the creation and absorption of electrons and positrons in the weak interactions, and N accounts for out-of-equilibrium effects from creation/annihilation, To determine N , we self-consistently calculate the collision integral for creation/annihilation, C(p), where C a (p) is the collision integral for the annihilation pathway as in the previous section, and C c (p) is for the creation pathway. It follows that Further details on the collision integral calculations can be found in Appendix B. Together, Eqs. (3,5,11,12,22) yield a system of coupled differential equations for dependent variables T , η − , η + , Y e , and a. As with the chemical equilibrium solution, the definition of Y e , Eq. (10), is an integral relationship between all the dependent variables. While this statement again means that the differential equations are not independent, we simultaneously solve the five differential equations for the five variables, as we did in the chemical equilibrium solution. Figure 2 shows that the creation and annihilation rates are many orders of magnitude greater than the Hubble rate for T > 16 keV, which would require far too many time steps to integrate these equations from T = 10 MeV as we did with the chemical equilibrium solution. However, because we expect chemical equilibrium to hold for T > 16 keV, we can use the chemical equilibrium solution to determine initial conditions for this non-equilibrium solution. Figure 3 presents the solution to these fully nonequilibrium conditions between 24 and 12 keV, lasting Deviation from η + = −η − indicates a shift away from local chemical equilibrium, so the quantity ∆η = η − + η + parameterizes the magnitude of any such out-of-equilibrium effect. The initial conditions set the initial value of ∆η = 0, but the non-equilibrium solution immediately deviates from this initial condition. Throughout the evolution, ∆η η − , and the electron degeneracy parameter in the non-equilibrium solution (η − ) differs from the equilibrium solution (η in Figure 1) by an amount on the order of ∆η. As a result, the positron-to-electron ratio and the electron chemical potential as shown in Figure 3 slightly differ from the corresponding values in the chemical equilibrium solution with a fractional difference on the order of ∆η.
To test the validity of taking the initial conditions directly from the chemical equilibrium solution, we repeated the process of determining initial conditions at 25 keV, then evolved the non-equilibrium equations. Once again, the initial conditions are chosen such that the initial ∆η = 0 and the non-equilibrium solution immediately deviates from this value. We found that this solution was consistent with the results presented in Fig. 3 for 24 keV < T < 12 keV.
There are two possible observable consequences of this out-of-equilibrium process: BBN yields and N eff . The light element abundances formed in BBN are sensitive to the time-temperature relationship and the neutron-toproton ratio during the BBN epoch (∼ 100 keV), both of which are influenced by the electron and positron distribution functions. The observationally inferred value of N eff is affected by the electron/positron number densities and out-of-equilibrium entropy transfer between the neutrinos and the plasma. We find that electron-positron annihilation freeze-out occurs at a sufficiently late time that there should be negligible effects on BBN yields and N eff .
Out-of-equilibrium effects are small (∆η ∼ 10 −6 ) at T ≈ 24 keV, yet the BBN yields are set before this time, at higher plasma temperatures when ∆η is many orders of magnitude smaller. (see, e.g., Ref. [1] for a discussion of BBN with out-of-equilibrium Boltzmann neutrino transport.) As a result, we expect the BBN yields calculated with the standard assumption of thermal and chemical equilibrium in the electron/positron distributions to be insensitive to the loss of local chemical equilibrium observed in this calculation.
The value of N eff = 3 is predicated on a number of simplifying assumptions [21] that neglect a number of effects. Out-of-equilibrium scattering between neutrinos and electron/positrons distorts neutrino spectra, transfers entropy from the plasma to the neutrinos, and increases the total entropy of the universe [1]. Non-zero electron/positron mass and finite-temperature QED affect the equation of state for electrons, positrons, and photons in the plasma. Together, these effects increase N eff by ∆N eff ≡ N eff − 3 ≈ 0.046 [1][2][3]14]. The initial non-zero electron/positron degeneracy parameters are responsible for a negligible increase, ∆N eff ≈ 10 −19 , due to a higher initial plasma entropy. On the other hand, the non-zero final number density of electrons (and negligible number density of positrons) results in a negligible decrease, ∆N eff ≈ −10 −7 , because some of the entropy remains in the electron seas at recombination.

V. DISCUSSION AND CONCLUSIONS
In this work, we found that the electrons and positrons remain in local chemical equilibrium throughout the BBN epoch and that this chemical equilibrium persists until the plasma temperature is approximately 16 keV, when the electron annihilation rate drops below the Hubble rate. This annihilation freeze-out occurs after the BBN yields have been set and when the number densities of electrons and positrons are sufficiently low such that the value of N eff is unaffected. To reach this conclusion, we self-consistently calculated the collision integral for electron-positron annihilation using the Klein-Nishina amplitude and the appropriate Fermi-blocking and Boseenhancing statistical factors to take into account out-ofequilibrium effects.
The currently accepted value of N eff is predicated on the complete annihilation of the electrons and positrons from thermal number densities that are nearly equal to the photon number density to essentially zero [21]. The slight changes in the degeneracy parameters of electrons and positrons when the chemical equilibrium approximation is abandoned does not change this conclusion -so as long as the standard cosmological model is assumed, annihilation freeze-out should have a negligible effect on the value of N eff .
BBN calculations assume chemical equilibrium between the electrons and positrons with the photons to model the creation/annihilation process (e − e + γγ) [18]. We found that this approximation holds well throughout the BBN epoch, so BBN calculations in the standard cosmological model are not affected by the specifics of the annihilation freeze-out process.
BSM physics or non-standard cosmological models may push the annihilation freeze-out earlier (or nuclear freeze-out earlier) so that they may influence these cosmological observables. Two possibilities that would require BSM physics to alter the standard cosmological model are a higher baryon-to-photon ratio during the BBN epoch and a higher effective electron mass.
While the baryon-to-photon ratio is directly measured from the CMB, this value is indicative of the baryon-tophoton ratio in the epochs immediately preceding photon decoupling, not necessarily at earlier epochs [7]. Although the deuterium and helium BBN yields and observations appear to form a concordance between the baryon-to-photon ratio at recombination and the BBN epochs, the lithium problem may suggest a more complicated evolution of this quantity. A higher photon-tobaryon ratio during the BBN epoch would increase the annihilation freeze-out temperature, albeit by a slight amount unless it is many orders of magnitude larger. This would drastically alter the BBN yields, but also require a late generation of entropy to reduce the baryon-to-photon ratio to the CMB measured value [22]. Such an event would be a signal of BSM physics, and such a significant entropy generation would be strongly constrained by N eff measurements [22] and possibly cascade nucleosynthesis that would determine the BBN yields [23].
Non-standard interactions between electrons and positrons and the plasma could increase the effective mass of the electrons and positrons [21]. A higher effective in-medium electron mass would cause the precipitous decline in electron and positron numbers seen in Figure  1 to occur at a higher temperature, in turn, resulting in the strong decrease of R e ± seen in Figure 2 occuring at a higher temperature as well. The result is a higher annihilation freeze-out temperature. The increase in the effective electron/positron mass due to finite temperature QED corrections is too small to significantly change the annihilation freeze-out temperature, so BSM interactions, and likely particles, would be required to produce such an effect, which brings up the same interesting ques-tions about the effects on BBN yields and N eff .
In conclusion, we have computed the collision integral for electron-positron annihilation throughout the weak decoupling and BBN epochs of the early universe. We find that electrons and positrons fall out of local chemical equilibrium with the photons at a sufficiently late time such that there will be no signature of this freeze-out on N eff or BBN yields. However, one should be cautioned that non-standard cosmologies may result in non-trivial consequences of this annihilation freeze-out.

ACKNOWLEDGMENTS
We'd like to thank G. Fuller for useful discussions. CK and LT acknowledge support from NSF Grant PHY-1812383. EG acknowledges support from NSF grant PHY-1630782 and the Heising-Simons Foundation, grant 2017-228.
with E p = p 2 + m 2 e , and T is the common temperature of all particles in the plasma. This should be an excellent assumption since electromagnetic scattering between the elements of the photon-electron-positron-baryon plasma occurs at a rate much greater than the Hubble rate throughout the calculation.
To self-consistently determine the equations of motion for T (t) and η(t), we need to calculate the energy density, ρ, pressure, P , and number density, n, of electrons and positrons from the distribution functions: These quantities can be expressed in terms of dimensionless functions by defining = p/T and x = m e /T , where E = √ 2 + x 2 . Further, define the partial derivatives of these functions with i = 1, 2. Each of these functions are, by design, dimensionless and positive definite.

Chemical equilibrium
In the chemical equilibrium solution, η e − = η and η e + = −η. The total entropy density is In the preceding equation, the '−' superscript corresponds to using the electron degeneracy parameter as the appropriate argument for each function, e.g., I − 1 = I 1 (x, η), while the '+' superscript corresponds to using the positron degeneracy parameter, e.g., I + 1 = I 1 (x, −η). While the entropy in the neutrino seas change through weak interactions beyond weak decoupling, it is a small effect that does not impact our results. We therefore assume that neutrinos are decoupled from the plasma, so that the total entropy in the neutrino seas in a comoving volume, s ν a 3 , is separately conserved.
The conservation of entropy in a comoving volume produces the differential equation: The evolution of the quantity (n e − − n e + )a 3 , Eq. (11), yields the differential equation: The quantity, n b a 3 , is the number of baryons in a comoving volume is constant and equal to its initial value. The evolution of Y e due to weak interactions, Eqs. (6)(7)(8) is Each of the weak rates are integrals over the electron distributions [19], so they depend on T and η. The evolution of the scale factor is given by the Friedmann equation, where T cm = T i (a i /a), and we have neglected the baryon energy density. These four coupled differential equations, Eqs. (A11-A14) define the evolution of the four dependent variables, T , η, Y e and a using the initial conditions discussed in Section II.
It should be noted, however, that Y e can be calculated directly (without a differential equation) from We chose to treat Y e as a dynamical variable with its own differential equation, but confirm that the above relationship between Y e , T , a, and η hold.

Out of chemical equilibrium
When we no longer assume chemical equilibrium, we are left with two independent degeneracy parameters, η e − = η − and η e + = η + . The total entropy density is We use the same notation as previously used in this Appendix, where '−' corresponds to electrons, e.g., I − 1 = I 1 (x, η − ), and '+' corresponds to positrons, e.g., I + 1 = I 1 (x, η + ). The conservation of entropy in a comoving volume yields The most easily distinguished difference between this equation and the corresponding one in the equilibrium scenario, Eq. (A11), are the separate degeneracy parameters, η − and η + , and their derivatives.
Likewise, in this scenario, the evolution of (n e − − n e + )a 3 yields the differential equation: The evolution of Y e and a, Eqs. (A13,A14), remain unchanged. However, it should be noted that the weak rates depend on both η − and η + . Now that there are two independent degeneracy parameters, we require an extra independent differential equation. We introduced the evolution of (n e + + n e − )a 3 , Eq. (22), which results in the differential equation: As before, the n b a 3 is left explicitly uncancelled in the expressions above because it is constant. The net electronpositron annihilation rate, N , is defined in Eq. (25). The five coupled differential equations, Eqs. (A13,A14,A17-A19) define the evolution of five dependent variables, T , η − , η + , Y e and a, whose initial conditions and solutions are discussed in Section IV.

Appendix B: Simplification of the Collision Integral
This appendix details the reduction of the collision integral for electron-positron annihilation from a ninedimensional integral to three-dimensional in a manner similar to the appendix in Ref. [1]. For the annihilation process, e − + e + γ + γ, the Klein-Nishina amplitude is [24] |M| 2 = 2e 4 P · K 2 P · K 1 + P · K 1 P · K 2 + 2m 2 where e is the elementary charge, m e is the electron mass, P and Q are the four-momenta of the electron and positron, respectively, and K 1 and K 2 are the photon four-momenta. We define the four-momenta of the electron as P = (E, p) with particle energy E and three-momentum p such that P ·P = E 2 −p 2 = m 2 e , with p = |p|, and the same convention for positrons. Likewise, for the photon, K = (k, k) with photon three-momentum k and energy k = |k|.
The collision integral is and is the statistical factor for the creation (F c ) and annihilation (F a ) of electrons and positrons with the appropriate Fermi blocking and Bose enhancement terms. For the case F = F a , we define C a as the collision integral for solely the annihilation pathway (e − e + → γγ) and F = F c defines C c for the creation pathway (γγ → e − e + ).
Integrating over k 2 reduces the delta function to a one-dimensional energy-conserving delta function: where K 2 = (k 2 , p + q − k 1 ), and k 2 is no longer an independent variable of integration, but is related to the other integration variables through momentum conservation, and defining the new integration angle θ 1 angle between the vectors k 1 and p + q. We need to evaluate the four-dot products which requires us to parameterize the dot products p · q and p · k 1 in terms of the variables of integration. First, define θ q as the angle between p and q so that p · q = pq cos θ q (B8) |p + q| = p 2 + q 2 + 2pq cos θ q 1/2 .
The vector p + q is in the plane generally defined by p and q, but k 1 is not necessarily in this plane. This means k 1 is defined by both θ 1 , the angle between it and p + q, and an azimuthal angle φ 1 measured out of this plane. To determine p · k 1 , we need to introduce another angle, ψ, between p and p + q. It follows that p · k 1 = pk 1 (cos ψ cos θ 1 − sin ψ sin θ 1 cos φ 1 ) , with cos ψ = p · (p + q) p|p + q| = p + q cos θ q p 2 + q 2 + 2pq cos θ q .
The delta function is non-zero on the interval [a, b] only if which will constrain the range of the other integrals once the delta function is resolved. To determine these ranges we need to consider two cases: Case 1: k 1 < |p + q|. The inequality has two parts, first: This is always true because |p+q| < p+q (triangle inequality), and momenta are always less than their corresponding energy. Further, the second inequality: The right-hand side is always positive, as stated previously. This reduces the range of the k 1 integral such that its integrand is non-zero. Case 2: k 1 > |p + q|. This case differs from Case 1 only in the first inequality, and again, the second inequality yields k 1 > k min . Putting the cases together, we find that the delta function is non-zero when k min < k 1 < k max , so upon evaluating the du integral the delta function fixes cos(θ 1 ) in terms of other variables: The collision integral then becomes C = 1 (2π) 4 16E p q 2 dq d cos q E q |p + q| kmax kmin dk 1 2π 0 dφ 1 M (p, q, θ q , k 1 , θ 1 , φ 1 ) F (E p , E q , k 1 , k 2 ), where, once again, k 2 = E p + E q − k 1 . Note that only M depends on φ 1 , so we can perform the integral in the square brackets analytically. To simplify our notation, we define the functioñ M (p, q, θ q , k 1 ) = 1 (2π)(2e 4 ) 2π 0 dφ 1 M (p, q, θ q , k 1 , θ 1 , φ 1 ).
The amplitude, Eq. (B1), depends only on the four-dot products P · K 1 and P · K 2 . These dot products can be written in the form A + B cos φ 1 . The functionM , above, can be analytically integrated as long as A > B. These four-dot products can be shown to be positive definite, which proves that A > B, by considering the annihilation