Computational forms for binary particle interactions at different levels of anisotropy

Particle interactions are key elements of many dynamical systems. In the context of systems described by a Boltzmann equation, such interactions may be described by a collision integral, a multidimensional integral over the momentum-phase space of the interaction. This integral is often simplified by assuming isotropic particle distributions; however, such an assumption places constraints on the dynamics of the system. This paper presents computational forms of the collision integral for relativistic, binary interactions at three levels of anisotropy, including a novel form in the isotropic case. All these forms are split into two parts, an absorption and an emission spectrum, which may be pre-calculated via numerical integration for simulation purposes. We demonstrate the use of these forms by comparison with the analytically integrated, isotropic emission spectrum of electron-positron annihilation, which are shown to agree to numerical precision. The emission spectrum is then further extended to axisymmetric particle distributions, where two-dimensional spectral maps can be generated to provide new insight.


I N T RO D U C T I O N
In the context of high-energy astrophysical sources, accurate modelling of processes, such as electron-positron pair annihilation/creation, Compton scattering, synchrotron, etc., is critical.Prior works on the emission and absorption spectra of such processes have focused mainly on isotropic interactions (Weaver 1976 ;Dermer 1984Dermer , 1985Dermer , 1986 ; ;Baring 1987 ;Coppi & Blandford 1990 ), with only small ventures into the world of anisotropy (Stepney & Guilbert 1983 ).The assumption of isotropy in these works is often a compromise between simplicity and complexity as the dimensionality of the integrals rapidly increases with anisotropy.
An example application is the emission modelling of astrophysical jets, where the current state-of-the-art models (Potter 2018 ;Cerruti et al. 2021 ;Lucchini et al. 2022 ;Zacharias et al. 2022 ;Klinger et al. 2023 ) evolve shells of isotropic non-thermal plasma travelling at a bulk Lorentz factor out from the central black hole.The constraint of isotropy prevents self-consistent 'hydrodynamical' evolution of the jets as plasma shells cannot interact with each other.In addition, isotropy is only a good approximation if the dynamical timescales of the system are much greater than isotropizing time-scales, which is not guaranteed at all energies.Work in the field of pairplasma winds (Aksenov, Milgrom & Usov 2004 ;Aksenov, Ruffini & Vereshchagin 2007, 2010 ) has demonstrated self-consistent evolution of anisotropic particle distributions using a kinetic approach 1 show-E-mail: christopher.everett@physics.ox.ac.uk 1 Treatment of integration bounds and domain conditions of the emission and absorption terms tends to be rather ambiguous; ho we ver, their approach appears to be consistent with this work's equations ( 12) and ( 13) if integrated o v er p 1 momenta.
ing that the additional dimensionality associated with anisotropy is within the grasp of modern computing.Ho we ver, a clear frame work in which to treat general anisotropic interactions is missing in the literature.
Although inspired by astrophysical sources, the work presented here is applicable to a wide range of physical systems where reactions, collisions, and interactions play an important role.To that aim, effort has been made to keep the forms of the emission and absorption spectra given in this paper as general as possible, while presenting them in a way that is easy for the reader to apply to any desired interaction.To maintain a reasonable scope for this paper, the forms developed here pertain to a specific subset of interactions, though may be easily extended upon in future works.The subset consists of all binary interactions that are reversible, i.e. take the form 12 34.No assumptions are made about the species of particles involved in the interactions, i.e. they can be massive or massless, nor about their energy/momenta; therefore, the forms given in this papers are relativistic and may be evaluated in an arbitrary frame of reference.
The format of the rest of this paper is as follows.In Section 2 , a general o v erview, conte xt, and definitions of reaction rates, emission, and absorption spectra are provided.In Section 3 , a framework in which to treat systems that exhibit axisymmetry in their particle distributions is developed, and then Section 4 examines the isotropic cases.Section 5 identifies and examines some singular cases of the inte grals dev eloped in the previous two sections that are of particular note and develops new forms dedicated to these cases.In Section 6 , an example application of this work to electron-positron annihilation is given, by comparison between the computational forms developed in this work and the known analytical form of the isotropic emissivity.The comparison demonstrates the replicability of the isotropic emissivity without need for difficult and time-RASTAI 00, 1 (2024) consuming algebraic integration, and hence may be applied to a wider range of problems.The comparison is then further expanded by examining the axisymmetric emission spectrum for two test cases, an axisymmetric case and a cosmic ray case -generating a two-dimensional (2D) spectral map of the emissivity for each that pro vides no v el insight.

OV E RV I E W
The dynamics of a system that involves interaction between particles may generally (ne glecting e xternal forcing) be described in kinetic theory by a set of relativistic Boltzmann equations (one for each species in the system): where ) is the collision integral, and f ( x μ , p ) is the distribution function defined by f ( x μ , p )d 3 x d 3 p = N , the number of particles.The on-shell condition relating energy to momentum m 2 = ( p 0 ) 2 − p • p constrains f such that it is only a function of time, position, and momentum, not energy.Henceforth, the spacetime positional dependence x μ of terms will be dropped for brevity.
The collision integral C( p 1 ) 3 is itself an inte gral o v er momentum states (hence its name) involving the product of distribution functions f ( p i ) and transition rates W (Groot, Leeuwen & Weert 1980 ;Cercignani 2002 ): 4 where δ (4) is the 4D Dirac delta function, σ 12 | 34 ( s, t) ≡ is the Lorentz-inv ariant dif ferential cross-section, which for a binary interaction is only dependent on the Mandelstam variables s, t defined in Appendix A , and p * in is the magnitude of the incoming particle momentum in the centre of momentum frame.In the case of the forward interaction 12 → 34, where λ( s, m 2 1 , m 2 2 ) is K äll én's function (K äll én 1964 ) and F 12 ( s) is the invariant flux of incoming particles.The additional factor of 1 1 + δ 34 has been included in equation ( 2 ) so as not to o v ercount ( p 3 , p 4 ) momentum states when the particles are indistinguishable (see Cercignani 2002 , section 2.1 for a discussion).A similar factor accounting for the indistinguishability of ( p 1 , p 2 ) states is not needed as this cancels with the 'true' collision integral for that particle being C( p 1 ) + C( p 2 ).Integrating C( p 1 ) o v er the Lorentz invariant momentum-space volume element d 3 p 1 p 0 1 is equal to the reaction rate 2 The quantity p 0 = E/c, i.e. has dimensions of momentum, with c = 1 and a metric signature ( + − −−) being used throughout this text. 3Numbered subscripts refer to particles in the binary interaction 12 34, which may be of the same species. 4For the interactions rele v ant to this work, all distributions are assumed to have Boltzmann statistics.
R 1 , which is the net rate of change of particles of type 1 per volume of space: Performing the same integral over equation ( 1 ) yields where n = f ( p )d 3 p is the number density and v is the bulk velocity.The aforementioned form may be more familiar than equation ( 1) particularly when the nature of the distribution of particles o v er momentum space is not needed or rele v ant to the system.
For practical purposes and to align with previous works, the collision spectrum is defined here as Following the notation of Baring ( 1987 ), the collision spectrum splits into two portions: The absorption spectrum accounts for the loss of type-1 particles with particular momentum p 1 due to the forward reaction 12 → 34 (also referred to as the 'loss term' in literature), while the emission spectrum accounts for the gain of type-1 particles into a particular state p 1 from the reverse reaction 34 → 12 (also referred to as the 'gain term').
The absorption integral does not depend on p 3 , p 4 apart from in the transition rate equation ( 3 ), which can integrated over to give where σ 12 | 34 ( s) is the cross-section for the interaction 12 → 34. 5ubstituting this in equation ( 10) gives yields the simple form this is a form that, in theory, can be integrated depending on the nature of the distribution functions.
The emission integral is less simple as now there is dependence on p 3 , p 4 via their distribution functions and hence the simplification of equation ( 11) cannot be used.Substituting equation ( 3 ) into equation ( 9), noting p * in = F 34 ( s) √ s as this is the reverse reaction, and where is the Heaviside step function.Equations ( 12) and ( 13) are the simplest general forms of the emission and absorption spectrum when the distribution functions are assumed to be fully anisotropic in momentum space. 6

A X I S Y M M E T R I C M O M E N T U M S PAC E
Systems which exhibit a global preferred direction, e.g.ˆ r for a spherically expanding shock or B in a magnetized flow, permit their distribution functions to be described as locally 7 axisymmetric in momentum space.Local axisymmetry occurs when any momentumorthogonal forcing term, e.g.v × B , averages to zero on rele v ant time-and/or length-scales.
Working in spherical coordinates for momentum space ( p, θ, φ), where p = | p | , momentum-space axisymmetry about the polar axis yields that all properties derived from the distribution function are independent of φ.With this property, it is convenient to replace the distribution functions f ( p ), absorption spectrum T spe ( p ), and emission spectrum S spe ( p ), with their axisymmetric counterpart e.g.f ( p, θ ) = 2 πp 2 f ( p ), which equates to integrating over the momentum-space element p 2 d φ.In other words, f ( p ) equal to a constant describes a distribution of particles whose density per unit energy is constant, whereas f ( p ) being a constant describes a distribution that increases in density with energy.

Axisymmetric absorption spectrum
Under the assumption of axisymmetry, the absorption spectrum (equation 12 ) has no dependence on φ 1 , φ 2 apart from the terms containing s.A change of variables can be instated from d φ 2 → d s, to decouple d φ 1 , then integrating equation ( 12) over p 2 1 d φ 1 yields the axisymmetric absorption spectrum T spe ( p 1 , θ 1 ): with The s integration bounds in equation ( 14) must be further constrained by the physical region of the Mandelstam variables (see 6 In principle, the delta function can be used to integrate over one of the six integrals involved.This approach has not been chosen, as the choice of integral to eliminate is situation dependent.Ho we ver, Appendix B provides a useful form, whereby the emission spectrum is inte grated o v er p 1 in order to remo v e the delta function. 7Locally is taken to refer to a local spatial axis; e.g. the distribution of particles from a spherically symmetric implosion may be described as being locally axisymmetric about the local radial direction but requires a connection coefficient to translate the momentum-space coordinates from one spatial point to another. Appendix A ), in this case there is only one rele v ant condition which is equation ( A5 ).
When the particle species are restricted to be massless, this form is equi v alent to Stepney & Guilbert ( 1983 , equation B5), when inte grated o v er d p 1 d cos θ 1 .

Axisymmetric emission spectrum
Under the assumption of axisymmetry, the emission spectrum (equation 13 ) has no dependence on φ 1 , φ 3 , φ 4 apart from the terms containing s, t and u .Integrating over p 2 1 d φ 1 and changing variables from d φ 1 d φ 3 d φ 4 → d sd td u gives The delta function can be integrated over using the s, t, or u integrals.
Choosing the u integral with and the u is given by equation ( A4 ).When integrating the delta function, the domain of u integral dictates that the emission spectrum is non-zero only for Additional constraints on the Mandelstam variables, described in Appendix A , also apply, in particular equations ( A5 -A8 ).

I S OT RO P I C M O M E N T U M S PAC E
Systems may be described as locally isotropic in momentum space if all forcing terms may be averaged to zero and sufficient collisions have occurred on rele v ant time-and/or length-scales.
Again with working spherical momentum-space coordinates, isotropy is described by distribution functions that are independent of θ and φ.With this property it is convenient to replace the distribution functions f ( p ), absorption spectrum T spe ( p ), and emission spectrum S spe ( p ) with their isotropic counterparts, e.g.f ( p) = 4 πp 2 f ( p ), which equates to integrating over the momentum-space element p 2 sin θd θd φ= p 2 d Ω.

Isotropic absorption spectrum
Under the assumption of isotropy, equation ( 12) then has no dependence on θ 1 , θ 2 , φ 1 , φ 2 apart from the terms containing s. Integrating o v er p 2 1 d Ω 1 yields the isotropic absorption spectrum T spe ( p 1 ): The angular variables may be freely rotated to a primed frame aligned to p 1 with d Inte grating o v er these decoupled components gives the simplest form: with and the condition equation ( A5 ) still applying.This result can be made identical to Baring ( 1987 , equation 27) with the further substitutions ε * 2 = s/ 4 and f ( p

Isotropic emission spectrum
Under isotropy, the emission spectrum equation ( 13) has no dependence on θ 3 , θ 4 , θ 1 , φ 3 , φ 4 , φ 1 apart from the terms containing s, t and u .Inte grating o v er momentum-space solid angle p 2 1 d Ω 1 yields the isotropic emission spectrum: ) Akin to the absorption spectrum (Section 4.1 ), the angular variables may be rotated to a primed frame where they are aligned to where = φ 4 − φ 1 .Making the change of variables d θ 1 d θ 4 d φ 4 → d td sd u , decouples Ω 3 and φ 1 from the system and the delta function can be integrated over to remove u dependence and yield: where Due to the delta function integration, u is bound by equation ( 16e ) and is given in terms of s and t by equation ( A4 ).The additional constraints given by equations ( A5 -A8 ) also apply to equation ( 23a ).
The result of equation ( 23) is identical to Baring ( 1987 , equation 30), without the need to consider conversions of collision angles to the centre of mass frame, making this form easier to apply in general.

S I N G U L A R C A S E S O F N OT E
An observant reader may have noticed that the forms presented in equations ( 14) and ( 16) and equations ( 18) and ( 23) are not well defined at the points where their integration bounds are equi v alent, i.e. s + = s − , t + = t − , or u + = u − , which occur when particles are aligned with the axis of symmetry or have zero momenta.These are coordinate singularities arising from the definition of the polar axis in momentum-space spherical coordinates.Here, we will describe how to deal with such points under the formalism detailed in previous sections by giving specific examples of the emission spectrum for three cases of particular use: the cosmic ray, lab, and collider.

Cosmic ray case
Cosmic rays may be considered to be plain parallel rays incident on to the Earth's atmosphere.As such, their momentum vector p 3 is best described as being aligned to the polar axis travelling downwards, i.e. θ 3 = π .Given the cosmic rays are incident from a single direction, their emission spectra will in general be axisymmetric.Ho we ver, if treated using equation ( 16), the s, t, u integrals in this case are not well defined as s + = s − and t + = t − , resulting in an unphysical value of zero for the emission spectrum for all output cases if numerical integration is used.To treat this point accurately, begin with the general case equation ( 13), replace components with their axisymmetric counterparts (as described at the beginning of Section 3 ) with the inclusion of a delta function as f ( p 3 , θ 3 ) = f ( p 3 ) δ( cos θ 3 + 1).With the integration over the delta function taking place first, the usual change of variables can be instated to yield the result: and u ± given by equation ( 16d ).Note that the factor of 1 + δ 34 has been remo v ed as the distributions of type-3 and type-4 particles have been distinguished even though they may be of the same species.For the cosmic ray case, it is also convenient to consider the population of target particles, i.e. type 4, to be isotropic.For such a case, the RASTAI 00, 1 ( 2024) emission spectrum is given by with t given by equation ( 24c), u by equation ( 24d), s ± given by equation ( 23b), and u ± by equation ( 23d ).

Lab case
The lab case is when one of the incident particles is at rest, e.g.p 4 = 0 and f ( p 4 ) = n 4 δ( p 4 ) , where n 4 is the number density of particles of type 4. In this case, both axisymmetric and isotropic p 3 cases can be considered.For the axisymmetric case, the emission spectrum is given by with and t ± given by equation ( 16c ).For the isotropic case, with s given by equation ( 26b ) and with The value of t in equation ( 26a) is subject to the additional constraint of with t ± given by equation ( 23c ) due to one of the delta function integrations.
The form of the emission spectrum is given by where and The remaining delta function indicates that p 1 and θ 1 are no longer independent variables.If p 3 = p 4 and it is assumed that incident particles are travelling at a known momentum, then f ( p 3 , θ 3 ) = n 3 δ( p 3 − p 3 ) δ( cos θ 3 + 1) and f ( p 4 , θ 4 ) = n 4 δ( p 4 − p 3 ) δ( cos θ 4 − 1) and the emission spectrum is then given by and t given by equation ( 28c ).

E X A M P L E : E L E C T R O N -P O S I T R O N A N N I H I L AT I O N E M I S S I O N S P E C T RU M
The angle-averaged (isotropic) emissivity E iso is a measure of the rate of emission of particles of particular momentum from an interaction as a function of the momentum/energy of the initial state.This was analytically generated by Svensson ( 1982 , equations 41 and 55-58) 8 for the particular case of photon emissions from electronpositron annihilation, i.e. e + e − → γ γ .The process of deriving the emissivity in this case involved significant algebraic manipulation and careful analysis of boundary conditions.As an example use case of the framework presented in Section 4 , Svensson's work can be replicated numerically with ease and may further be extended to the axisymmetric case using the methods presented in Section 3 .The angle-averaged emissivity is related in general to equation ( 23) by The Lorentz-invariant differential cross-section for the electronpositron annihilation to two photons is given by (Berestetskii, Lifshits & Pitaevskii 1982 , equation 88.4) 26) for the first panel and numerical integration of equation ( 31) for the latter two panels including the fractional error compared with Svensson ( 1982 , equations 55-58).The angular dependence of the line with a ‹ is explored in Fig. 2 .
where σ T is the Thompson scattering cross-section and m e is the electron mass.By numerically integrating 10 equation ( 31 ) using equation ( 32), Svensson's work is reproduced in Fig. 1 .The notation of Svensson's work has been converted to be consistent with this work.The incident electron/positron energy is denoted by p 0 e ± or p 0 e ∓ , with the use of ± there to represent the exchange symmetry of the emissivity.The fraction errors shown in Fig. 1 are within the precision goal of the numerical scheme, demonstrating the reliability and use of equation ( 23) o v er tedious analytical integration.Ho we ver, the fractional errors do tend to be positive, suggesting that numerical integration is underestimating the actual value, and can fluctuate around maximal points of the emissivity; hence, a small amount of care must be taken around these points.The emissivity in the first panel, where p 0 e ± = 1, is a singular case of equation ( 23) and has instead been calculated using equation ( 26).This form is analytical, and when applied to the case here of electron-positron annihilation is equi v alent to Svensson ( 1982 , equation 41).
Using the equations derived as part of this work, Svensson's isotropic emissivity may be extended to the axisymmetric case with an angle-dependent emissivity related to equation ( 16): The spectral lines of Fig. 1 may be expanded using equation ( 34 ) to 10 Numerical integration is performed using a GlobalAdaptive scheme in Mathematica version 13.2, with a PrecisionGoal of 3. reveal spectral 'islands', which display dependence on momentum and angle from the axis of symmetry.A single spectral line from Fig. 1 is selected (labelled with a ‹) and expanded in Fig. 2 for a specific set of incident angles.The incident electron-positron pair is travelling within two cones in momentum space, one with p 0 e ± = 2 at an angle of θ ± = 3 π/ 4 while the other at angle of θ ∓ = π/ 4 with a higher energy of p 0 e ∓ − 1 = 10 2 .Fig. 2 shows the angledependent axisymmetric emissivity of emitted photons as a function of p γ and θ γ .The emissivity peaks in regions closely located around what would be considered a 'glancing' collision, where outgoing photons match the trajectory of the incoming pair.This effect lines up with the peaks of the isotropic emissivity observed in Fig. 1 .The breadth of the emissivity in emitted photon energy is narrower than the isotropic case, being truncated at low momenta, with the smallest RASTAI 00, 1 (2024) momentum being around 10 0 , compared with ≈ 10 −0 . 8, indicating that given the incident angles, these photon momenta are no longer physically obtainable.Two sweeping ridges can be observed in Fig. 2 , a feature not present in the isotropic case, which relate to likely pairs of emitted photon states.Consider two photons emitted at the peaks of the emissivity: As one, e.g. the higher energy photon, deviates away from the peak in phase space, the other must travel in order to conserve momentum and energy.If it travels down the peak to the north-west of Fig. 2 , its cone of emission (remembering that these photons are emitted symmetrically about the axis) has converged closer to the axis of symmetry and it has lost energy.To conserve momentum, the lower energy photon must begin to travel north-east in order to broaden its cone of emission and gain energy.
For an alternate example of photon emissivity, Fig. 3 presents the particular case of high-energy positrons ( p 0 e + − 1 = 10 2 ) incident along the axis of symmetry, on an isotropic population of lowenergy electrons with energy p 0 e − = 2.This emissivity is generated from equation ( 25), as in a cosmic ray case, and peaks around the location of incident positron; ho we ver, unlike the case in Fig. 2 , the electrons have no preferred orientation in momentum space, leading to a broadened out secondary peak of photon emission around the electron energy ( p γ ≈ 2).It may also be noted that even though the incident electrons are isotropic, the emergent photons are not, even at the electron's energy where emission angles of greater than θ γ ≈ 150 • are unphysical.This may have utility in calculating the angular distribution of particles in air shower simulations.

C O N C L U S I O N
New computational forms for emission and absorption spectra have been generated for relativistic, binary interactions with varying levels of anisotrop y.These tak e the form of multidimensional integrals, of minimal dimensionality.Their use is analysed by comparison with previously disco v ered inte gral forms and is demonstrated to agree to numerical accuracy to the specific case of isotropic electronpositron annihilation, for which an analytically integrated form of the emission spectrum is known.
The specific case of electron-positron annihilation is briefly expanded upon to demonstrate the potential use of the new computational forms for axisymmetric emission spectrum presented in this w ork.Tw o cases are examined: photon emissions from an axisymmetric state of electron-positron pairs, and a cosmic ray case where high-energy positrons are incident on an isotropic distribution of low-energy electrons.A spectral island is generated for each, denoting the emissivity of photons as a function of momenta and angle from the axis of symmetry.
The computational forms presented here are expected to be useful in work requiring accurate description of collisional systems with anisotropic processes.Further papers will explore such scenarios, with an application to dynamical modelling of particle production and high-energy emissions in astrophysical jets. 9

Figure
Figure Electron-positron annihilation angle-averaged emissivity as a function of emergent photon and pair energy , p γ and p 0 e ± , respectively .This figure is a reproduction of Svensson ( 1982 , fig. 5) using the analytical form of equation (26) for the first panel and numerical integration of equation (31) for the latter two panels including the fractional error compared withSvensson ( 1982 , equations 55-58).The angular dependence of the line with a ‹ is explored in Fig.2.

Figure 3 .
Figure 3. Angle-dependent emissivity of photons from high-energy cosmic ray positrons at an angle of θ e + = π at an energy of p 0 e + = 1 + 10 2 incident on an isotopically distributed, low-energy population of electrons with p 0 e − = 2.